Predicted boron-carbide compounds: A first-principles study De Yu Wang, Qian Yan, Bing Wang, Yuan Xu Wang, Jueming Yang, and Gui Yang Citation: The Journal of Chemical Physics 140, 224704 (2014); doi: 10.1063/1.4882071 View online: http://dx.doi.org/10.1063/1.4882071 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/140/22?ver=pdfcov Published by the AIP Publishing Articles you may be interested in First principle study of elastic and thermodynamic properties of ZrZn2 and HfZn2 under high pressure J. Appl. Phys. 115, 083514 (2014); 10.1063/1.4867221 First principle study of elastic and thermodynamic properties of FeB4 under high pressure J. Appl. Phys. 114, 183517 (2013); 10.1063/1.4829926 First-principles structural design of superhard materials J. Chem. Phys. 138, 114101 (2013); 10.1063/1.4794424 Structure and mechanical properties of tantalum mononitride under high pressure: A first-principles study J. Appl. Phys. 112, 083519 (2012); 10.1063/1.4759279 Pressure-induced phase transition and mechanical properties of molybdenum diboride: First principles calculations J. Appl. Phys. 112, 013522 (2012); 10.1063/1.4733954

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Sat, 23 Aug 2014 04:50:45

THE JOURNAL OF CHEMICAL PHYSICS 140, 224704 (2014)

Predicted boron-carbide compounds: A first-principles study De Yu Wang, Qian Yan, Bing Wang, Yuan Xu Wang,a) Jueming Yang, and Gui Yang Institute for Computational Materials Science, School of Physics and Electronics, Henan University, Kaifeng 475004, People’s Republic of China

(Received 21 November 2013; accepted 27 May 2014; published online 12 June 2014) By using developed particle swarm optimization algorithm on crystal structural prediction, we have explored the possible crystal structures of B-C system. Their structures, stability, elastic properties, electronic structure, and chemical bonding have been investigated by first-principles calculations with density functional theory. The results show that all the predicted structures are mechanically and dynamically stable. An analysis of calculated enthalpy with pressure indicates that increasing of boron content will increase the stability of boron carbides under low pressure. Moreover, the boron carbides with rich carbon content become more stable under high pressure. The negative formation energy of predicted B5 C indicates its high stability. The density of states of B5 C show that it is p-type semiconducting. The calculated theoretical Vickers hardnesses of B-C exceed 40 GPa except B4 C, BC, and BC4 , indicating they are potential superhard materials. An analysis of Debye temperature and electronic localization function provides further understanding chemical and physical properties of boron carbide. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4882071] I. INTRODUCTION

The search for new superhard materials is an important field of research in materials science and technology due to their importance in fundamental science and technological applications. As well known, diamond is the hardest material, which has been applied in a variety of industries.1 However, applications of diamond due to its poor resistance to oxidation as well as reaction with ferrous metals. The demand for superhard materials in industrial usages has stimulated research in the synthesis of alternate hard materials. It is known that superhard materials are usually made of light elements such as B, C, N, and O. Therefore, great efforts have been devoted to synthesizing various compounds constituted by these light elements. Novikov et al. reported that they synthesized cubic boron nitride(c-BN), and it has a hardness of 62 GPa.2 Diamondlike BC5 phase with the measured hardness of 71 GPa has been synthesized at 24 GPa and about 2200 K,3 which has a similar x-ray diffraction pattern as diamond. Zinin et al. synthesized a cubic BC3 phase by direct transformation from graphitic phases at a pressure of 39 GPa and temperature of 2200 K in a laser-heated diamond anvil cell.4 Ab initio calculation is a powerful tool to explore the atomic structures and physical properties of materials. Several theoretical studies have been performed to investigate the structures and physical properties of B-C system. A density functional theory (DFT) study has predicted that BC5 is metallic and superconducting with a critical temperature of 45 K, comparable to the critical temperature of MgB2 .5 Yao et al.6 noted that a tetragonal I-4m2 structure is the most stable structure of BC5 . A tetragonal BC2 phase originating from the cubic diamond structure was predicted by first-principles calculations and its theoretical Vickers hardness is 56 GPa.13 Nakae et al.14 investigated the stability of BC5 relative to a a) E-mail: [email protected]

0021-9606/2014/140(22)/224704/8/$30.00

mixture of BC3 and graphite under pressures using density functional theory, and found that high pressure synthesize is very useful for BC5 . Li el at.15 explored the crystal structures of the synthesized superhard diamond-like BC5 through ab initio evolutionary algorithm, and their simulated XRD patterns and Raman modes suggest that Pmma-1 and Pmma-2 structures are the best candidate structures of superhard BC5 . Liu et al.16 have explored the crystal structures of synthesized diamond-like BC3 (d-BC3 ) with particle swarm optimization (PSO) algorithm, and found that the simulated Raman modes of Pmma-b phase are in agreement with the experimental data. Recently, Mikhaylushkin et al. predicted two stable boron carbides (BC3 and BC5 ) and found they are more stable than those previously reported.17 However, information about the atomic structure of B-C system with rich boron content has not been well studied up to now. It is known that crystal structures of materials are the basis for understanding of any related properties. Therefore, it is necessary to predict their possible configuration and the corresponding electronic and mechanical properties. In this work, we have extensively investigated the possible structures of B-C system by using ab initio particle swarm optimization (PSO) algorithm on crystal structural prediction,18–22 which requires only the chemical compositions for a given compound at specified external conditions, unbiased by any known structural information. We also studied their corresponding electronic and mechanical properties.

II. COMPUTATIONAL DETAIL

To search for potential crystal structures, the PSO technique implemented in the Crystal Structure Analysis by Particle Swarm Optimization (CALYPSO) package was employed at 0 GPa and 80 GPa with 1–6 formula units (f.u.) each simulation cell.18 For the PSO calculations, the maximum

140, 224704-1

© 2014 AIP Publishing LLC

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Sat, 23 Aug 2014 04:50:45

224704-2

Wang et al.

distance between any two k-points in the k-point sampling grid is less than 0.04 Å−1 . The underlying density functional theory calculations were performed by using the Vienna ab initio simulation package (VASP).24 The projectoraugmented wave (PAW) method with a plane-wave basis set was used to describe the interaction between ion cores and valence electrons.23 The valence-electron configurations considered are 2s2 2p1 for B and 2s2 2p2 for C. Geometry optimization was performed using the conjugate gradient algorithm method with a plane-wave cutoff energy of 500 eV. The generalized gradient approximation (GGA) as parameterized by Perdew, Burke, and Ernzerhof25 was used to describe the exchange-correlation function. For the hexagonal structures,  centered k mesh was used, and for other structures, Monkhorst-Pack k-points were used, and the k-point sampling grid was set to be 0.02 Å−1 . The cutoff energy and k-point sampling grid were tested to ensure that all structures are well-converged to better than 5 meV/atom. The structures were relaxed with respect to both lattice parameters and atomic positions. For the relaxed structures, the force on each ion and the target pressure are less than 0.005 eV/Å and 0.1 kbar, respectively. The smearing methods for relax structure calculations and density of states calculations are the Methfessel-Paxton scheme and the tetrahedron method with ¨ corrections, respectively. The strain-stress method was Blochl used to obtain the elastic constants. From the calculated elastic constants Cij , the polycrystalline bulk modulus B and shear modulus G were further estimated using the Voigt-ReussHill approximation.26 In addition, the Young’s modulus Y and Poisson’s ratio ν were obtained by the equations Y = (9GB)/(3B+G), and ν = (3B-2G)/(6B+2G), respectively. The phonon calculations were performed within the framework of force constants method. Density functional perturbation theory (DFPT) with VASP was employed to calculate real-space force constants. The phonon frequencies were calculated from the force constants using Phonopy code.27 The lattice constants of supercells are all larger than 10 Å to obtain accurate force constants. For calculating force constants, the Gaussian smearing method was used.

J. Chem. Phys. 140, 224704 (2014)

FIG. 1. Crystal structures of boron carbide in B-C systems. The black and white spheres represent B and C atoms, respectively. (a) BC4 : No. 131, the B atom is at 2e (0, 0, 1/4), and the C atom is at 4i (0, 1/2, 0.3545), 4h (1/2, 1/2, 0.9308); (b) B3 C: No. 107, the B atom is at 4b (1/2, 0, 0.7298), 2a (1/2, 1/2, 0.0748), and the C atom is at 2a (0, 0, 0.3581); (c) B2 C: No. 146, the B atom is at 9b (0.4553, 0.3287, 0.4250), 9b (0.0114, 0.3373, 0.6256), and the C atom is at 9b (0.3284, 0.3213, 0.5922); (d) BC: No. 136, the B atom is at 4f (0.7375, 0.7375, 0), and the C atom is at 4g (0.8667, 0.8667, 1/2).

graphite), where H for elemental boron was calculated based on the structure of the α-B12 at 0 GPa, From Table I, it can be seen that the formation energies of B-C system are positive except of B5 C, indicating their metastable stability. The negative enthalpy of B5 C means that it is possible to

III. RESULTS AND DISCUSSIONS A. Structure and stability

We performed variable-cell structure prediction simulations using PSO methodology for B-C system containing one, two, three, four, five, and six f.u. in the simulation cell at 0 and 80 GPa. The earlier proposed structures of I41 /amd BC2 ,13 P4m2 BC3 ,16 and I-4m2 BC5 6 were successfully reproduced, validating our methodology in application to B-C system. Besides the known candidate structures, some new stable structures were found as shown in Figs. 1 and 2. The relaxed structure parameters, density, and atom volume of boron carbides with different boron concentration are tabulated in Table I. The calculated results show that the equilibrium lattice parameters of BC2 , BC3 , and BC5 are in excellent agreement with the previous data,6, 13, 16 indicating the reliability of present calculations. The formation energies (H) of B-C system were calculated by H = H(Bm Cn ) - mH(solid B) - nH(solid

FIG. 2. Crystal structures of R3c-B4 C (1 × 2 × 1 supercell) (a) and P3m1B5 C (1 × 1 × 2 supercell) (b). The black and white spheres represent B and C atoms, respectively. For B4 C, the B atom has two sites: 36f (0.2149, −0.4657, −0.9679), 36f (0.1212, −0.3668, −0.8716), and the C atom is at 18e (0.2587, 0, 3/4). For B5 C, the B atoms has four sites: 6i (0.4823, 0.5181, 0.25074), 6i (0.30477, 0.15238, 0.24874), 6g (0.33482, 0.33482, 0), 2d (1/3, 2/3, 0.62355), and C atoms have two sites: 2d (1/3, 2/3, 0.3669), 2c (0, 0, 0.6295).

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Sat, 23 Aug 2014 04:50:45

224704-3

Wang et al.

J. Chem. Phys. 140, 224704 (2014)

TABLE I. Calculated formation energy per unit (H), optimized equilibrium lattice parameters a and c (Å), density (ρ in g/cm3 ), and atom volume (V in Å3 /atom) of B-C system.

B4 C

B3 C B2 C BC BC2 BC2 a BC3 BC3 b BC4 BC5 BC5 c

H

a

P3m1 P63 /mmc P6m2 R3/c R3 P4/nmm R3m R3m I4mm R3 P42 /mnm I41 /amd

− 0.697 1.079 0.560 0.322 0.778 0.859 1.521 1.814 0.558 0.571 0.397 0.745 0.562 1.106

5.5460 4.2443 4.2418 5.9471 5.3697 2.7948 5.9255 6.4570 2.7746 4.8546 4.3761 2.5234 2.5200 2.5104 2.5015 2.5621 2.5174 2.5250

P4m2 P42 /mmc I4m2

1.477 1.291

c

13.7413

8.7913

7.2229 9.1113 2.7591 11.9262 11.9190 3.9186 3.915 10.0945 11.3923 11.323

ρ

V

2.731 2.076 2.074 2.709 2.380 2.709 2.335 2.528 2.688 2.733 2.890 3.060 3.060 3.162

6.79 8.93 8.94 7.00 7.82 6.87 7.97 7.22 6.95 6.89 6.61 6.33

2.958 3.267

6.63 6.02

6.17

a

Reference 13, VASP. Reference 16, VASP. c Reference 6, VASP. b

synthesize it under ambient pressure. The metastable character is also found in the recent study of the synthesized c-BC5 .3 We hope that the predicted B-C phases could be synthesized using the high-pressure and high-temperature method under appropriate temperature and pressure ranges. B4 C and B5 C are well-known boron rich carbides and are widely used for polishing, grinding, and surface hardening. They have been synthesized by some experimental workers and were found to be p-type semiconductors.7–9 For example, Pasquale and Kelber synthesized semiconducting boron carbide (B10 C2 Hx ) films formed by the crosslinking of B10 C2 H12 icosahedral.8 They found that the chemical environment of the icosahedral cage B atoms changed during crosslinking of the icosahedral. Moreover, the B–H bond was destroyed by electron bombardment. Therefore, it is desirable to study the structure of B5 C and B4 C. The structure of B4 C are attracting great interest and have been explored by various theoretical methods.10–12 Many previous works suggested that B4 C crystal were derived from the B15 structure which was formed by linking of B12 icosahedral with BBB atom chain. B4 C can be formed by replacing BBB chain by C atoms with different atom arrangements. Ivashchenko and Shevchenko investigated three possible structures of B4 C crystal and found that the 15-atom rhombohedral cell with the polar arrangements atoms is most stable among various structures.10 However, this rhombohedral structure is not stable and will transform to a monoclinic structure after optimizing. Such monoclinic lattice distortion may decrease the total energy and open a wide bandgap in its band structure. Hence, it is interesting to explore the possible structure of B4 C with the rhombohedral symmetry. By using the PSO method, we explored possible structure of B4 C and found that the structure with the lowest energy has a rhombohedral symmetry (R3c)

B4C-B4C B2C+2B-B4C

0 Δ H (eV)

B5 C

Symmetry

4 B3C+B-B4C B5C-B-B4C

-4

BC+3B-B4C 4BC2-7C-B4C 4BC3-11C-B4C

-8

4BC4-15C-B4C

-12

4BC5-19C-B4C

20

40 60 Pressure (GPa)

80

100

FIG. 3. Calculated enthalpy under pressure of boron carbides relative to the R3c phase of B4 C.

(shown in Fig. 3(a)). This rhombohedral structure also contains boron icosahedral, and these adjacent icosahedral are connected by carbon atoms. Moreover, the carbon atom forms a B-C-B chain with the boron atoms in polar sites. Each carbon atom is coordinated by four boron atoms with bond lengths of 1.59 and 1.62 Å. The boron atom is coordinated by one carbon atom and five boron atoms. The B–B bond lengths are in the range of 1.72−1.86 Å. Our predicted stable structure of B4 C is different from that proposed by previous work,10 and the present predicted R3c-B4 C keeps rhombohedral symmetry after optimizing. In addition, in this R3c structure, the carbon atoms do not replace the polar-site boron atoms. In Ref. 10, Ivashenko suggested that, in the lowest energy structure, the polar-site boron atoms were replaced by carbon atoms. However, their proposed structure (labeled as POL-15 ) cannot keep rhombohedral symmetry after optimizing. The present PSO method also predicted other three rhombohedral structures (R3, R3m, and R3m). The R3 structure has B11 C icosahedral, and the icosahedral is far away from next one. Previous X-Ray diffraction work suggested that B4 C may have a rhombohedral structure with R3m space group. Our predicted R3m structure is higher in energy than other three rhombohedral structures. For B5 C, our predicted stable structure has a trigonal symmetry (P3m1). The boron atoms form 16-boron-atoms cages which share faces, and these boron cages are linked by carbon atoms along the c axis. Each carbon atom is coordinated by four boron atoms with band lengths of about 1.60 Å. There are two types boron atoms in B5 C. The first-type boron atom is coordinated by one carbon atom and five boron atoms with bond lengths of 1.60–1.91 Å. The second type one is coordinated by seven boron atoms with bond lengths of 1.57– 1.93 Å. The shortest B–B bond in B5 C has a much smaller length than that in R3c-B4 C. Thus, for boron carbides with boron rich content, P3m1-B5 C is possibly more stable than R3c-B4 C, which is consistent with the calculated formation energies in Table I. This trigonal B5 C structure may be different from the structure of the experimentally synthesized B10 C2 Hx films in Ref. 8. In such films, the B10 C2 Hx cages are cross-linked after electron bombardment. For B3 C, five carbon atoms and 13 boron atoms form boron-carbon cages which are linked along the c direction. The boron (carbon) atoms in polar site are coordinated by one carbon (boron) atom and four boron (carbon) atoms with B–C bond lengths of 1.57 and 1.67 Å. The B–B bond lengths in these cages are 1.78 and 1.96 Å. For B2 C, six boron atoms

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Sat, 23 Aug 2014 04:50:45

224704-4

Wang et al.

J. Chem. Phys. 140, 224704 (2014)

form a trigonal prism. These trigonal prisms are linked by carbon atoms. There are two types boron atoms in B2 C. The firsttype boron atoms are coordinated by two carbon atoms and four boron atoms with B–C bond lengths of 1.64 and 1.65 Å. The second-type boron atoms are coordinated by three carbon atoms and tow boron atoms with B–C bond lengths of 1.61, 1.64, 1.70 Å. Each carbon atom of B2 C is coordinated by five boron atoms. The B–B bond lengths in B2 C are 1.76 and 1.84 Å which are slight smaller than those in B3 C. For BC, six boron and six carbon atoms form face-shared decahedrons linking along the c direction. In BC, each boron atom is coordinated by four carbon atoms with bond lengths of 1.59 and 1.69 Å, while each carbon atom is done by five boron atoms with bond lengths of 1.59, 1.65, and 1.69 Å. In BC4 , four boron atoms and six carbon atoms form boron-carbon cages which are linked by carbon atoms along the c direction. Each boron atom is coordinated by four carbon atoms with bond lengths of 1.66 Å. There are two-types carbon atoms in BC4 , The first-type carbon atoms are coordinated by two boron atoms and two carbon atoms with C–C and B–C bond lengths of 1.49 and 1.66 Å, respectively. The second-type carbon atoms are located at polar sites and are coordinated by three carbon atoms with bond lengths of 1.40 and 1.49 Å, and the shortest C–C bonds are formed by two adjacent carbon atoms at the polar sites of boron-carbon cages. As shown in Fig. 1(a), the two adjacent boron-carbon cages along the c direction are linked by only one C–C bond, which may induce a weak ability to against shear elastic deformation and decreases the hardness of BC4 . Because boron carbides were often synthesized under high pressure, we calculated the enthalpy of these predicted boron carbides under pressure range of 0–100 GPa. In order to compare B4 C with other predicted boron carbides, we plotted the enthalpy of each predicted boron carbides relative to that of the R3c phase of B4 C in Fig. 3. As seen in this figure, B5 C

is always more stable than R3c-B4 C. R3c-B4 C has a stronger stability than B3 C, B2 C, and BC with the studied pressure range. With more increasing of carbon content, the predicted boron carbides become less stable than B4 C under low pressure. However, BC2 , BC3 , BC4 , and BC5 become more stable than R3c-B4 C as the pressure increasing to 29 GPa, 27 GPa, 53 GPa, and 16 GPa, respectively. This means that boron carbides with rich carbon content can be synthesized under high pressure. For example, BC5 and BC3 has been synthesized at 24 GPa and 39 GPa, respectively.3, 4 B. Lattice vibration and elastic properties

The mechanical stability is a necessary condition for a crystal to exist. Accurate elastic constants can help us to understand the mechanical properties and also provide very useful information to estimate the hardness of a material. The calculated elastic constants of the B-C system are listed in Table II. It is clearly seen that all studied compounds satisfy the mechanical stability criteria,28 indicating that they are elastically stable. On the other hand, the positive eigenvalues of the elastic constant matrix for each considered compound further prove that they are elastically stable. The large values of C11 , C22 , and C33 for these boron carbides suggest that they are difficult to be compressed along a-axis, b-axis, and c-axis, respectively. Dynamic stability is important for a new structure, as the appearance of soft phonon modes will lead to the distortion of crystal. It is necessary to check dynamical stabilities of the predicted structures of boron carbides. Therefore, we calculated the phonon dispersion curves at high symmetry points of predicted boron carbides at 0 GPa as shown in Figs. 4– 6. As seen in these figures, no imaginary phonon frequency were detected in the whole Brillouin zone of the boron carbide structures except R3-B4 C, indicating that they are all

TABLE II. Calculated elastic constants (in GPa), bulk modulus (B, B0 in GPa), shear modulus (G in GPa), Young’s modulus (Y in GPa), Poisson’s ratio ν, Debye temperature (D in K), and hardness (Hv in GPa) of B-C system.

B5 C(P3m1) B5 C(P63 /mmc) B5 C(P6m2) B4 C(R3c) B4 C(R3) B4 C(P4/nmm) B4 C(R3m) B3 C(I4mm) B2 C BC BC2 BC3 BC4 BC5 Diamond

C11

C12

C13

C33

C44

C66

B(B0 )

G

B/G

Y

ν

D

Hv

565 365 354 478 477 595 367 586 603 502 701 654 730 870 1047

136 124 125 111 73 35 235 50 66 263 53 158 40 76 128

74 99 102 193 76 76 110 87 131 74 217 180 147 169

684 352 358 328 319 815 530 787 564 919 629 778 858 840

249 49 37 185 113 212 157 213 245 244 388 443 62 477 560

215 121 115 188 204 125 66 107 269 297 230 261 140 318

265(265) 192(192) 192(191) 252(251) 188(186) 261(251) 241(239) 264(261) 270(259) 302(300) 333(333) 353(355) 329(325) 378(386) 434

268 91 78 176 161 251 133 242 273 254 338 371 155 444 564

0.99 2.11 2.45 1.43 1.17 1.04 1.81 1.09 0.99 1.19 0.99 0.95 2.12 0.85 0.77

601 236 206 428 376 570 337 556 613 595 758 824 402 957 1181

0.12 0.30 0.32 0.22 0.17 0.14 0.27 0.15 0.13 0.17 0.12 0.11 0.30 0.08 0.05

1276 999 737 1016 1047 1655 949 1753 2033 2224 2199 2004 1282 1961

50 9 6 24 30 45 14 42 51 39 58(56a ) 65(66b ) 13 82(80c ,71d ) 107(115e )

a

Reference 13, theory. Reference 16, theory. c Reference 6, theory. d Reference 3, experiment. e Reference 37, experiment. b

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Sat, 23 Aug 2014 04:50:45

224704-5

Wang et al.

J. Chem. Phys. 140, 224704 (2014)

40

Frequency/THz

50 40

(b)

30

(c)

30

30 20

20

10

10

20 10 0

ΓA

HK

Γ

ML H (d)

Frequency/THz

40

(a)

0

Z

AM

ΓZ

Γ

RX

40

(e)

0

30

30

30

20

20

20

10

10

10

0

Z

AM

Γ Z RX Γ

0

Z

A M

Γ Z R X Γ

Z

AM

ΓZ

RX

40

0

Γ (f)

Z

AM

ΓZ

RX

Γ

FIG. 4. Calculated phonon dispersion curves at 0 GPa: (a) P63 /mmc-B5 C, (b) P4/nmm-B4 C, (c) B3 C, (d) B2 C, (e) BC, and (f) BC2 .

dynamically stable except R3-B4 C. Because boron carbides was often synthesized under high pressure, we also calculated the phonon dispersion curves of these boron carbides at 50 GPa and found that they are all stable at 50 GPa except R3-B4 C. The phonon dispersion is also a powerful tool to analyze the ability to against elastic deformation of a material. The higher maximum frequency of acoustic branch means the higher ability to against elastic deformation. As seen from Figs. 4 and 5, the maximum frequency of acoustic branch for BC5 , BC3 , and BC2 are all above 20 THz, indicating their high hardness. For BC4 , two ultrasoft degenerate acoustic modes appear along the -Z direction, meaning a relative weak ability to against shear deformation in BC4 . Hence, the hardness of BC4 should be relative low among these studied boron carbides, which is consistent with the calculated hardness shown in Table II. The ultrasoft modes along the -Z direction of BC4 may be related to its bonding feature along the c direction. As seen in Fig. 1(a), the two adjacent boron-carbon cages along the c direction are connected by only one C–C bond. Such bonding behavior is possibly related to the appearing of ultrasoft acoustic modes along the -Z line. For long wavelength limitation, acoustic wave is elastic wave. As seen from Fig. 6, for long wave-length limitation, the transverse acoustic wave velocity of B5 C is much larger than that of B4 C. In

Frequency/THz

40

(a)

other words, B5 C has a higher shear wave velocity than B4 C. Thus, B5 C should have a stronger ability to against the shear deformation than B4 C, which suggests that B5 C may have a higher hardness than B4 C. The results of hardness shown in Table II confirms this conclusion. The acoustic vibration branches along the L-H line in B4 C are also softer than those in B5 C, indicating a possible higher hardness of B5 C than B4 C. The difference in longitudinal acoustic wave velocity of B4 C and B5 C is not so significant. For B5 C, two separated optic modes appears in high frequency zone, which may be related to its relative high hardness. Bulk and shear moduli are still the most important parameters for identifying the hardness of a material. High hardness requires that bulk modulus and shear modulus must be as large as possible. The bulk modulus, shear modulus, Young’s modulus, and Poisson’s ratio of the boron carbides have been estimated from the calculated elastic constants. In order to confirm the reliability of the present calculation method, the third order Birch-Murnaghan equation of state was employed to calculate bulk moduli (B0 ). The calculated results show that the bulk moduli (B) are in good agreement with the measured B0 , demonstrating the reliability of the present theoretical method. Being well known, superhard materials always have a high bulk modulus to resist the volume change caused by an applied load. From Table II, it can be seen that

50

(b)

40

30

40

(c)

30

30 20

20 20

10 0

10

10 Z

AM

Γ Z RX Γ

0

Z

AM

ΓZ

RX

Γ

0

Z

AM

ΓZ

RX

Γ

FIG. 5. Calculated phonon dispersion curves at 0 GPa: (a) BC3 , (b) BC4 , and (c) BC5 .

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Sat, 23 Aug 2014 04:50:45

Wang et al.

224704-6

J. Chem. Phys. 140, 224704 (2014)

40

B4C (R-3c)

30 20 10 0 Frequency / THz

40

Γ

A

H

Γ

K

M

L

H B4C (R3)

30 20 10 0 40

Γ

A

H

Γ

K

M

L

H

B5C (P-3m1)

30 20 10 0

Γ A

H K

Γ

M L

H

FIG. 6. Calculated phonon dispersion curves of R3c-B4 C, R3-B4 C, and P3m1-B5 C at 0 GPa.

all ground-state boron carbides have high bulk moduli (above 250 GPa), which indicates that these materials are difficult to be compressed. It also indicates that these boron carbides can be grouped into low compressible materials. Moreover, with the increasing of carbon composition, the bulk modulus increases. Among the predicted boron carbides, BC5 has the largest bulk modulus. Compared with bulk modulus, the shear modulus of a material quantifies the resistance to the shear deformation and is a better indicator of potential hardness. Except B4 C and BC4 , other studied boron carbides have high shear moduli (above 240 GPa), indicating their strong resistance to shape change at a constant volume. BC5 has the largest shear modulus (444 GPa), suggesting that it can withstand the shear strain to the largest extent. In addition, the shear moduli of B5 C, B2 C, BC2 , BC3 , and BC5 exceed their bulk moduli. The ratio value of B/G is commonly used to describe the ductility or brittleness of materials, with 1.75 as the critical value.28 A B/G value higher (or lower) than the criteria is considered to be ductile (or brittle). From Table II, it can be seen that the B/G values of all the predicted phases are under the critical value except BC4 , implying their ductile nature. In addition, the calculated B/G ratio of the predicted boron carbides is larger than that of diamond, suggesting that they are more ductile than diamond. Materials with small Poisson’s ratio are more easily compressed than sheared (small B/G), whereas those with high Poisson’s ratio resist compression in favor of shear (large B/G).29 ν is also expected to decrease in cross-linked structures, which opposes transverse contraction upon tensile loading. The relatively low values of B/G and hence ν for hard materials in general is an indication of the high degree of directional covalent bonding. From Table II, the calculated diamond has the highest directional covalent bonding. BC4 has

the lowest directional covalent bonding. Young’s modulus (Y) is used to provide a measure of the stiffness of a solid, which is defined as the ratio of stress and strain. A large value of Y manifests a stiff material. BC5 has the largest Young’s modulus, which is comparable to that of diamond, indicating its high hardness. All of these excellent mechanical properties strongly suggest that the predicted phases are potential candidate of ultra-incompressible and hard materials except BC4 . Mechanical properties (in particular, hardness of solids) can be related to thermodynamical parameters such as Debye temperature, specific heat, thermal expansion, and melting point.30 In this concept, Debye temperature (D ) is one of the fundamental parameters for the characterization of materials and the microhardness of a solid.31, 32 The calculated D of boron carbides from their elastic moduli and density of materials are listed in Table II. We can see that they have high Debye temperature, in which BC has the highest D . Hardness is a highly complex property, which depends on the loading force and on the quality of samples (i.e., the presence of defects such as vacancies and dislocations). Because Vickers hardness is experimentally measured as a function of the applied loading forces, the saturated hardness value (or experimental load-invariant indentation hardness) is usually considered to be the hardness value of a given material. In order to figure out the theoretical hardness of B-C system, we employ the first-principle model of Chen et al.33, 34 to get the hardness of boron carbide. The results were shown in Table II, and the previous results are also listed for comparison. This is the empirical formula based on the Pugh modulus ratio k = G/B.35 The formula for calculation is as follows: Hv (GPa) = 2(k 2 G)0.585 − 3.

(1)

From Table II, the calculated hardness of BC2 , BC3 , BC5 , and diamond agree with the previous theoretical and experimental values,3, 6, 13, 16, 36, 37 demonstrating the reliability of the present theoretical method. Except B4 C, BC, and BC4 , the hardnesses of other boron carbide exceed 40 GPa, indicating that they are potential superhard materials.

C. Electronic structure and chemical bonding

In order to gain better understanding of the electronic structure and chemical bonding of boron carbide, the total and partial density of states (DOS) are calculated (shown in Fig. 7). Previous experimental works suggested that B5 C and B4 C are p-type semiconductor.7, 8 As shown in Fig. 7, our predicted P3m1-B5 C is a p-type semiconductor with a narrow bandgap of 0.4 eV, which agrees with the experimental result.8 The DOS around the Fermi energy is mainly attributed to the B p states. For B4 C, there are two DOS peaks from −1 to 1 eV, which mainly come from the B p states and the C p states. The DOS shape of B p and C p from −1 eV to 1 eV are very similar, indicating a strong hybridization between them. Hence, the B and C atoms may form strong covalent bonding in R3c-B4 C. Above 0.9 eV, there is a small gap of 0.5 eV, indicating the p-type semiconducting behavior of R3c-B4 C. The DOS peaks below −9 eV are mainly attributed to B s states and C s states. The states above −9 eV mainly originate from

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Sat, 23 Aug 2014 04:50:45

224704-7

Wang et al.

J. Chem. Phys. 140, 224704 (2014)

4

15 BC

B4C

10 2

Density (states/eV)

5 0 -9 15

-6

-3

0

3

6

0 -9 15

-6

-3

0

3

B2C

10

10

5

5

0 -9 4

-6

-3

0

3

6

0 6

0

B3C

BC4

4 2 2 0 -9

-6

-3 0 E-EF (eV)

6 B5C

3

6

0 -9

-6

-3 0 E-EF (eV)

Total B-p B-s C-p C-s

3

6

FIG. 7. Calculated total and partial DOS of boron carbides. The Fermi level is at zero.

B p and C p states with slight contributions of B s and C s states. As seen in Fig. 7, B2 C is semiconducting with a very narrow bandgap. B3 C, BC, and BC4 are metallic due to their finite values on the Fermi level. The appearing of pseudogap around the Fermi energy indicates that BC is very stable, which is confirmed by its low formation energy as shown in Table I. For B3 C, the major states near the Fermi level stems from B p electrons, which is the principal cause for the metallicity. In BC4 , the major states near the Fermi level stems from B p electrons and C p electrons, and its C s and B s electrons are located in low energy region. Moreover, the partial DOS profiles for both B p and C p in BC are very similar, indicating the significant hybridization between the two orbitals. This fact shows a strong covalent interaction between the B and C atoms. With the increasing of carbon composition, C p electrons play an increasingly important role in the DOS of B-C system. To gain more detailed insight into the bonding character of boron carbide, we have calculated their electronic localization function (ELF), which offers a reliable measure of electron pairing and localization. ELF is based on the HartreeFock pair probability of parallel spin electrons and is widely used to describe and visualize chemical bonding in molecules and solids.38 The ELF is employed to understand the localized bonding in real space, which is scaled between 0 and 1. ELF = 1 corresponds to the perfect localization characteristic of covalent bonds or lone pairs (filled core levels), while ELF = 0 is typical for a vacuum (no electron density) or areas between atomic orbitals. ELF = 0.5 for a homogeneous electron gas, with values of this order indicating regions with bonding of a metallic character. It should be noted that the ELF is useful in distinguishing metallic, covalent, and ionic bonding. The contours of ELF domains for B-C system are shown in Figs. S1– Fig. S3, respectively.39 As seen in them, the high electron localization can be seen in the region between B-B and B-C atoms in B5 C, B4 C, B3 C, and B2 C, indicating strong

covalent B–B and B–C bonding. The ELF of the B–C bonding exceeds 0.9, which is larger than that of B–B bonding. Hence, in these carbides, the B–C bonds have stronger covalency than the B–B bonds. Moreover, the unusual three-center bonding appear in the boron icosahedral of R3c-B4 C, which may lead to some exceptional properties. No B–B bonding were found in the predicted BC, BC2 , BC3 , BC4 , and BC5 . The ELF of the B–C bonding in BC, BC2 , BC3 , BC4 , and BC5 is also very large, indicating a strong B–C bonding in them. The strong bonding may be helpful to their hardness.

IV. CONCLUSION

In conclusion, using the developed PSO technique on crystal structure prediction, we report the possible structures of B-C system, which are mechanically and dynamically stable. Then, the structures, stability, elastic properties, and electronic properties of boron carbide with different boron concentration have been explored by first-principles calculations based on density functional theory. The results show that these boron carbides have high bulk moduli, which can be grouped into low compressible materials. Boron carbides with rich carbon content become more stable than those with high boron content under high pressure. Except B4 C and BC4 , other boron carbides have high shear moduli (above 240 GPa), indicating their strong resistance to shape change at a constant volume. The calculated theoretical Vickers hardnesses of the predicted boron carbides exceed 40 GPa except B4 C and BC4 , indicating that they are potential superhard materials. The DOSs of B5 C and B4 C show that they are p-type semiconductor. The B–B bonding exists only in B5 C, B4 C, B3 C, and B2 C, and their electronic localization function show that the strength of the B–C bonding is stronger than that of the B–B bonding. The phonon calculations have confirmed that all boron carbides are dynamically stable. The superior mechanical properties of B-C system highlight their potential

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Sat, 23 Aug 2014 04:50:45

224704-8

Wang et al.

applications as intrinsically superhard materials. We hope that these calculations will stimulate extensive experimental work on these technologically important boron carbide. ACKNOWLEDGMENTS

This research was sponsored by the National Natural Science Foundation of China (Nos. 21071045 and 51371076), the Program for New Century Excellent Talents in University (No. NCET-10-0132), and the Program for Innovative Research Team (in Science and Technology) in University of Henan Province (No. 13IRTSTHN017). 1 S.

Veprek, J. Vac. Sci. Technol. A 17, 2401 (1999). V. Novikov, S. N. Dub, and V. I. Mal’nev, Sov. J. Superhard Mater. (USRR) 5, 16 (1983). 3 V. L. Solozhenko, O. O. Kurahevych, D. Andrault, Y. L. Godec, and M. Mezouar, Phys. Rev. Lett. 102, 015506 (2009). 4 P. V. Zinin, L. C. Ming, H. A. Ishii, R. Jia, T. Acosta, and E. Hellebrand, J. Appl. Phys. 111, 114905 (2012). 5 M. Calandra and F. Mauri, Phys. Rev. Lett. 101, 016401 (2008). 6 Y. S. Yao, J. S. Tse, and D. D. Klug, Phys. Rev. B 80, 094106 (2009). 7 V. Domnich, Y. Gogotsi, M. Trenary, and T. Tanaka, Appl. Phys. Lett. 81, 3783 (2002). 8 F. L. Pasquale and J. A. Kelber, Appl. Surf. Sci. 258, 2639 (2012). 9 H. Werheit, J. Phys.: Condens. Matter 18, 10655 (2006) 10 V. I. Ivashchenko and V. I. Shevchenko, Phys. Rev. B 80, 235208 (2009). 11 M. M. Li, X. Fan, and W. T. Zheng, J. Phys.: Condens. Matter 25, 425502 (2013). 12 F. Mauri, N. Vast, and C. J. Pickard, Phys. Rev. Lett. 87, 085506 (2001). 13 L. F. Xu, Z. S. Zhao, L. M. Wang, B. Xu, J. L. He, Z. Y. Liu, and Y. J. Tian, J. Phys. Chem. C 114, 22688–22690 (2010). 14 N. Nakae, J. Ishisada, H. Dekura, and K. Shirai, J. Phys.: Conf. Ser. 215, 012116 (2010). 15 Q. Li, H. Wang, Y. J. Tian, Y. Xia, T. Cui, J. L. He, Y. M. Ma, and G. T. Zou, J. Appl. Phys. 108, 023507 (2010). 2 N.

J. Chem. Phys. 140, 224704 (2014) 16 H. Y. Liu, Q. Li, L. Zhu, and Y. M. Ma, Phys. Lett. A 375, 771–774 (2011). 17 A.

S. Mikhaylushkin, X. Zhang, and A. Zunger, Phys. Rev. B 87, 094103 (2013). 18 Y. C. Wang, J. Lv, and Y. M. Ma, Phys. Rev. B 82, 094116 (2010). 19 J. Lv, Y. C. Wang, L. Zhu, and Y. M. Ma, J. Chem. Phys. 137, 084104 (2012). 20 Y. C. Wang, J. Lv, L. Zhu, and Y. M. Ma, Comput. Phys. Commun. 183, 2063–2070 (2012). 21 M. Hu, F. Tian, Z. S. Zhao, Q. Huang, B. Xu, L. M. Wang, H. T. Wang, Y. J. Tian, and J. L. He, J. Phys. Chem. C 116, 24233–24238 (2012). 22 Z. S. Zhao, F. Tian, X. Dong, Q. Li, Q. Q. Wang, H. Wang, X. Zhong, B. Xu, D. L. Yu, J. L. He, H. T. Wang, Y. M. Ma, and Y. J. Tian, J. Am. Chem. Soc. 134, 12362–12365 (2012). 23 G. Kresse and D. Joubert, Phys. Rev. B 59, 1758 (1999). 24 G. Kresse and J. Hafner, Phys. Rev. B 47, 558 (1993). 25 J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). 26 R. Hill, Proc. Phys. Soc. A. 65, 349 (1952). 27 D. M. Ceperley and B. J. Alder, Phys. Rev. Lett. 45, 566 (1980). 28 Z. J. Wu, E. J. Zhao, H. P. Xiang, X. F. Hao, X. J. Liu, and J. Meng, Phys. Rev. B 76, 054115 (2007). 29 T. Rouxel, H. Ji, J. P. Guin, F. Augereau, and B. Rufflé, J. Appl. Phys. 107, 094903 (2010). 30 P. Ravindran, L. P. Fast, B. Johansson, J. Wills, and O. Eriksson, J. Appl. Phys. 84, 4891 (1998). 31 D. B. Sirdeshmukh, L. Sirdeshmukh, and K. G. Sirdeshmukh, Micro- and Macro-Properties of Solides: Thermal, Mechanical and Dielectric Properties, Springer Series in Material Science Vol. 80 (Springer, 2006), p. 315. 32 S. C. Abrahams and F. S. L. Hsu, J. Chem. Phys. 63, 1162 (1975). 33 X. Q. Chen, H. Y. Niu, D. Z. Li, and Y. Y. Li, Intermetallics 19, 1275 (2011). 34 X. Q. Chen, H. Y. Niu, C. Franchini, D. Z. Li, and Y. Y. Li, Phy. Rev. B 84, 121405 (2011). 35 S. F. Pugh, Philos. Mag. 45, 823 (1954). 36 V. Domnich, Y. Gogotsi, and M. Trenary, Mater. Res. Soc. Symp. Proc. 649, Q8.9.1–Q8.9.6 (2001). 37 N. V. Novikov and S. N. Dub, J. Hard Mater. 2, 3–11 (1991). 38 A. D. Becke and K. E. Edgecombe, J. Chem. Phys. 92, 5397–5403 (1990). 39 See supplementary material at http://dx.doi.org/10.1063/1.4882071 for the calculated electronic localization function (ELF) for the studied boron carbides.

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Sat, 23 Aug 2014 04:50:45