Defects in boron carbide: First-principles calculations and CALPHAD modeling

Available online at www.sciencedirect.com

Acta Materialia 60 (2012) 7207–7215 www.elsevier.com/locate/actamat

Defects in boron carbide: First-principles calculations and CALPHAD modeling Arkapol Saengdeejing ⇑, James E. Saal, Venkateswara Rao Manga, Zi-Kui Liu The Pennsylvania State University, University Park, PA 16802-5007, USA Received 7 August 2012; received in revised form 4 September 2012; accepted 7 September 2012 Available online 13 October 2012

Abstract The energetics of defects in B4+xC boron carbide and b-boron are studied through first-principles calculations, the supercell phonon approach and the Debye–Gru¨neisen model. It is found that suitable sublattice models for b-boron and B4+xC are B101(B,C)4 and B11(B,C) (B,C,Va) (B,Va) (B,C,Va), respectively. The thermodynamic properties of B4+xC, b-boron, liquid and graphite are modeled using the CALPHAD approach based on the thermochemical data from first-principles calculations and experimental phase equilibrium data in the literature. The concentrations of various defects are then predicted as a function of carbon composition and temperature. Ó 2012 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. Keywords: CALPHAD; First-principles calculations; Debye–Gru¨neisen model; Boron carbide; b-Boron

1. Introduction Boron carbide, usually referred to as B4+xC (x = 0, . . . , 7), exists as a single phase between 10 and 20 at.% C [1] and melts at temperatures higher than 2700 K. It is the third hardest naturally occurring material, just below diamond and cubic boron nitride. The combination of high melting temperature, outstanding hardness, high elastic modulus, low density and very high resistance to chemical agents [2] makes B4+xC a prime material for lightweight ballistic armor, high-temperature thermoelectric conversion and refractory applications. The nuclear industry has also used B4+xC as shielding and in control rods because of its ability to absorb neutrons without forming any long-life radionuclides. The large solubility of carbon in B4+xC is due to its intriguing crystal structure and the chemical similarity between carbon and boron. Large vacancy concentrations have also been detected [3]. The interactions between boron, carbon and vacancy are essential in understanding the behavior of the phase. ⇑ Corresponding author. Tel.: +81 22 795 4830.

E-mail address: [email protected] (A. Saengdeejing).

B4+xC has the R3m structure type [4–6], closely related to the a-boron structure, as shown in Fig. 1. It consists of a 12-atom icosahedron and a 3-atom chain distributed over four Wyckoff positions: equatorial (6h) sites along the center of the icosahedron, polar (6h) sites that form the two triangular caps of the icosahedron, the two ends of the 3atom chain (2c), and the center of the chain (1b) [3]. The atom at each chain end is bonded covalently to an atom along the equatorial girth of three different icosahedra. Each icosahedron is connected to six other icosahedra through one of the three polar sites on the top and bottom of the icosahedron. Kaufman et al. [7] modeled B4+xC with a simple solution model using the pair potential method to estimate the enthalpy of formation of solution phases. Kasper et al. [8] used the associate model, (B12, B11) (CBC,CBB,BVaB), to model its defects and thermodynamic properties based on phase-equilibrium information and very limited thermochemical data, which is clearly an oversimplification for B4+xC. Recently, Saal et al. [9] proposed that mixing occurs in the icosahedron for C-rich compositions and in the chain for B-rich compositions based on the crystal structures, enthalpy of formation and infrared phonon modes predicted from first-principles

1359-6454/$36.00 Ó 2012 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.actamat.2012.09.029

7208

A. Saengdeejing et al. / Acta Materialia 60 (2012) 7207–7215

Fig. 1. Rhombohedral primitive cell of B4+xC, drawing by VESTA 3 [43], showing four different Wyckoff positions with bonds to the adjacent atoms in the cell: equatorial (1), polar (2), chain end (3) and chain center (4).

calculations and available experimental data. A more comprehensive model is thus needed to represent the complex defect structures in B4+xC and the corresponding thermodynamic properties. The present work aims to develop a more complete understanding of defect mechanisms through first-principles calculations and improve the modeling of thermodynamic properties of the compound and the B–C system [10].

composition of the eutectic reaction of graphite, B4+xC, and liquid was measured by Kieffer et al. [17] as 26 at.% C, about 3% lower than the data by Elliott [1], and the liquidus temperatures reported by Kieffer et al. [17] from 0 to 16 at.% C are higher than those noted by Elliott [1]. The solubility of boron in graphite was reported by Lowell [18] for a wide temperature range with the highest value of 2587 K at around 2.35 at.% C. Kuhlmann et al. [19] found that carbon is soluble at 1 at.% C in boron and proposed that the solubility limit of carbon in b-boron could be as high as 3.75 at.% C. Using density functional theory, Shang et al. [20] showed that a-boron is more stable than bboron at low temperatures in the defect-free structure and predicted the a to b transition temperature at 1388 K. From experiment, a-boron is not stable at low temperatures. Widom and Mihalkovic [21] calculated the energies of defect structures in b-boron and reported that it is possible to make b-boron more stable than a-boron by introducing partially occupied sites into the b-boron structure. 3. Sublattice models One sublattice model has been adopted for the liquid and graphite phases. Their Gibbs energy functions are represented by the following equation: Gm ¼ xB o GB þ xC o GC þ RT ðxB ln xB þ xC ln xC Þ þ ex DGmix ð1Þ where xB and xC are the mole fractions of B and C, and GB and oGC are the Gibbs energies of pure B and C. The third term represents the ideal mixing between species in the sublattice. exDGmix is the Gibbs energy due to non-ideal atomic interactions and is defined using the Redlich–Kister polynomial [22]: X ex k DGmix ¼ xB xC LB;C ðxB  xC Þk ð2Þ o

2. Experimental data in the literature There are five phases in the B–C system: liquid, a-boron, b-boron, B4+xC and graphite. Limited experimental information is reported in the literature. Froment et al. [11] measured boron activity across the B4+xC phase at 2300 K using dense graphite Knudsen-cell mass spectrometry. They also measured the boron activity at 673 K using a solid-state potentiometric cell [12]. The heat capacity was measured by Matsui et al. [13] at BxC (x = 3, 4, 5) in both single- and two-phase regions. Elliott [1] investigated the phase equilibria using X-ray, metallographic and thermal analysis methods of sintered and arc-melted specimens. A high-purity (>99.6 at.%) of boron was used to determine the carbon solubility limits in B4+xC phase as 9.1 at.% C at the B-rich side and 19.8 at.% at the C-rich side at 2300 K. The B4+xC solubility limit has been reported by others as well [14–16]. Bouchacourt and Thevenot [14] studied hot-pressed boron carbide from 0 to 60 at.% C and determined the carbon solubility limit in B4+xC to be 8.8 and 20 at.% C at the B- and C-rich ends, respectively. Beauvy [15] reported 21.6 ± 0.8 at.% C solubility at the C-rich end and the congruent melting point of B4+xC to be about 18.2 at.% C at 2723 K. The most recent determination of solubility limit in B4+xC was reported by Schwetz and Karduck [16] to be about 18.8 at.% C from electron probe microanalysis. The liquid

k¼0

where kLB,C is the kth-order interaction parameter for mixing. The interaction parameters are expressed as a linear function of temperature by kL = a + bT. It was observed experimentally in B4+xC that vacancies can replace the chain center atoms [3], and carbon can only replace a single boron at a polar site in each B12 icosahedron [23]. When all three atoms in the chain are replaced by vacancies, the structure is equivalent to a-boron. Consequently, Saal et al. [9] proposed the (B)11(B,C) (B,C,Va) (B,Va) (B,C,Va) sublattice model to represent the B4+xC solution phase. The first sublattice represents the B in the icosahedron, the second sublattice is the carbon and boron mixing in the polar site of icosahedrons, the third and fifth sublattices are the two ends of the chain with B, C, and vacancy (Va), and the fourth sublattice is the center of the chain with B and Va. With this sublattice model, the number of end-members, i.e. structures that have only one species in each sublattice, is 36. From those 36 endmembers, 24 of them are distinguishable based on the

A. Saengdeejing et al. / Acta Materialia 60 (2012) 7207–7215

symmetry of the chain end site. One of the end-members, B11CCBC or CCBC for short, represents the B4+xC structure at 20 at.% C. Based on the above model, the range of carbon content in B4+xC is from 0 to 21.43 at.% C. The Gibbs energy function for the B4+xC solution phase per mole of formula, i.e. 5 + x moles of atoms, is thus written as:   þ xC Ggraphite GmB4þx C ¼ ð5 þ xÞ xB GbB B C XXXX IV V abcd þ y IIa y III  T Df S abcd Þ b y c y d ðDf H a

þ RT

b

c

i

j

XX

d

y ij ln y ij þ ex DGmix

ð3Þ

where a, b, c and d are the species in sublattices II, III, IV and V, respectively, y ij is the site fraction of species j in sublattice i, and DfHabcd and DfSabcd are the enthalpy and entropy of formation of the end-member abcd with respect to b-boron and graphite. Again, exDGmix is the excess Gibbs energy due to non-ideal interactions similar to Eq. (2) and can be written as follows: X XXX  k ex IV V DGmix ¼ y III y IIB y IIC k LB:B;C:b:c:d y IIB  y IIC b yc yd c

b

d

XXX V þ y IIa y IV c yd a

c

(

k

d

X  III  III k III k  y III B y C LB:a:B;C:c:d y B  y C þ

X

k

 III  III k III k y III B y Va LB:a:B;Va:c:d y B  y Va

k

þ

X

III k y III C y Va LB:a:C;Va:c:d

 III k y C  y III Va

)

k

X XXX  IV  V IV k IV k þ y IIa y III y IV b yd B y Va LB:a:b:B;Va:d y B  y Va a

b

d

a

b

c

k

( X XXX  k IV þ y IIa y III y VB y VC k LB:a:b:c:B;C y VB  y VC b yc þ

X

k

 k y VB y VVa k LB:a:b:c:B;Va y VB  y VVa

k

þ

X

y VC y VVa k LB:a:b:c:C;Va



y VC

k  y VVa

) þ ...

ð4Þ

k

Eq. (4) only shows the excess Gibbs energy for single sublattice mixing. Correlated simultaneous mixing in more than one sublattice is not shown here due to the complexity of the formalism and space limitations but it is included in the current model. The b-boron crystal structure is very complex and has the same space group as B4+xC, R 3m, and the unit cell can be considered as a face-centered cubic-like structure consisting of a B84 quasisphere with B10  B  B10 chains located in the octahedral interstices [24]. According to Werheit [25], carbon can replace one of the B atoms in each of the four B12 icosahedra of the unit cell, so that a maximum of four C atoms can be placed in a single unit

7209

cell, corresponding to 3.75 at.% C solubility limit. Therefore, the (B)101(B,C)4 sublattice model is selected to describe the solubility of the carbon in b-boron. Its Gibbs energy function per mole of formula can be represented by the following equation:    GmB4þx C ¼ y IIB 105GbB þ 4Ggraphite þ y IIC 101GbB B B C  þ Df H B101 C4  Df S B101 C4  T   ð5Þ þ 4RT y IIB ln y IIB þ y IIC ln y IIC þ ex DGmix where y IIB and y IIC are the site fractions of B and C in the second sublattice, and exDGmix is the excess Gibbs energy of mixing similar to Eq. (2). 4. First-principles calculations and the Debye–Gru¨neisen model For given atomic structures and atomic numbers such as the end-members, first-principles calculations based on density functional theory are performed to obtain the total energy of the structure, which is used to calculate its enthalpy of formation. In the present work, the projector augmented-wave (PAW) method [26] within the generalized gradient approximation (GGA) [27] as implemented in VASP [26,28,29] is used. The Brillouin zone integration is carried out using the Methfessel–Paxton order 1 smearing method [30]. An energy cutoff of 520 eV is employed. The k-mesh is set be 8  8  8, i.e. about 5000 k-points per reciprocal atom, and tested to ensure convergence of the total energy to a precision of better than 0.1 meV atom1 for all the calculations. A final self-consistent static calculation with the tetrahedron method with Blo¨chl corrections [31] is performed to obtain accurate total energy and electronic density of states (DOS). The free energy of all end-members of B4+xC and bboron as shown in Eqs. (3) and (4) can be predicted as a function of temperature and volume as follows [10]: F ðV ; T Þ ¼ E0 ðV Þ þ F vib ðV ; T Þ ð6Þ where E0 is the total energy at 0 K, and Fvib is the free energy contribution from the vibrations of atoms. Thermal electronic excitations can also be considered but are not necessary since B4+xC is an insulator. The equilibrium volume at a given temperature is obtained by minimizing the free energy with respect to volume. E0(V) is fitted to the four-parameter Birch–Murnaghan equation of states as suggested by Shang et al. [32]: E ¼ a þ bV 2=3 þ cV 4=3 þ dV 2 ð7Þ where a, b, c and d are parameters evaluated from the total energy calculated at a series of fixed volumes from which the bulk modulus and equilibrium volume at 0 K can be calculated. The vibrational contribution, Fvib, can be obtained from either the phonon properties [33] or the Debye–Gru¨neisen model [34] with the former being more accurate but computationally more expensive. Shang et al. [32] systematically compared the results from these two approaches

7210

A. Saengdeejing et al. / Acta Materialia 60 (2012) 7207–7215

and found out that the Debye model gives reasonably good results. Furthermore, some structures in the present work are mechanically unstable at zero Kelvin, and the Debye model is currently the only practical way to estimate their thermodynamic properties at finite temperatures. In the Debye model, Fvib is defined as: F vib ¼ Evib  TS vib

ð8Þ

where Evib is the vibrational energy, and Svib the vibrational entropy. They are given as follows [34]:   HD Evib ¼ E0K þ 3k B TD ð9Þ T      HD 4 HD ð10Þ  ln 1  e T S vib ¼ 3k B D 3 T where D(x) is the Debye function: Z x t3 DðxÞ ¼ 3x3 dt t 0 e 1

ð11Þ

E0K is the 0 K zero point energy, generated from fluctuations at the quantum level due to the Heisenberg uncertainty principle and defined by: 9 E0K ¼ k B HD 8

ð12Þ

HD is the characteristic Debye temperature calculated from the following equation [34]:  12  c rB V0 ð13Þ HD ¼ sA M V where s is a scaling parameter, A a constant equal to 231.104 with the bulk modulus, B, in GPa and the inter˚ ,M the average atomic mass, V the atomic radius, r, in A volume, V0 the equilibrium volume at 0 K, and c the Gru¨neisen parameter. 5. Results and discussion 5.1. First-principles The enthalpy and entropy of formation of the 24 distinguishable end-members in the sublattice model of B4+xC shown in Eq. (3) are obtained from first-principles calculations and the Debye–Gru¨neisen model. They can be found in Table 1 as the parameters in the Gibbs energies of the end-members. Saal et al. [9] performed the first-principles calculations for all the end-members of B4+xC at 0 K using equilibrium volumes and PAW-GGA pseudopotentials. It was predicted that four end-members have negative enthalpies of formation with the rest having positive values. One of those is the stable structure, CCBC end-member, of B4+xC at 20 at.% C with x = 0. The other three are BCBC, BCVC and CCVC. The current enthalpies of formation agree with the previously predicted values but have slightly lower values due to different reference state for enthalpy of formation calculations. In Saal et al., the reference state is

a-boron, while the current work uses b-boron as the reference state for boron. Table 2 lists the equilibrium volume (V0), bulk modulus (B)  0 and first derivative of the bulk modulus to the pressure BP for all 24 end-members.The equilibrium volumes for ˚3 all of the end-members range from 6.627 to 8.277 A 1 atom and the bulk modulus is about 193 GPa for most end-members, except BVBV and CVBV at 130 and 149 GPa, respectively. There are two parameters in the Debye–Gru¨neisen model in Eq. (13): the Gru¨neisen parameter, c, and the scaling parameter, s. The Gru¨neisen parameter describes the anharmonic effects in the vibrating lattice and is calculated from the following equation: 2 V @ 2 P =@V 2 c¼  3 2 @P =@V

ð14Þ

Moruzzi et al. [34] suggested that the scaling parameter is closely related to the anisotropy of bulk modulus and evaluated this parameter to be 0.617 from Eq. (13) with known Debye temperatures for body-centered cubic nonmagnetic elements. In the present work, the Debye temperatures of the four mechanically stable end-members are obtained from phonon calculations with the harmonic approximation, and their scaling parameters calculated from Eq. (13) are listed in Table 3. It can be seen that the temperature difference between the Debye temperature from phonon calculations and the Debye–Gru¨neisen model is less than 50 K for the four stable end-members of B4+xC if the scaling parameter of a-boron is used. Therefore, the scaling parameter of a-boron is adopted for all endmembers of B4+xC in the present work. Fig. 2 shows the calculated heat capacity of the CCBC end-member from the Debye–Gru¨neisen model, phonon calculations, the present CALPHAD (CALculation of PHAse Diagram) modeling and experimental data [13] at 20 at.% C. It can be seen that the heat capacity with the Neumann–Kopp approximation, i.e. with the heat capacity of formation being zero, in the present CALPHAD approach is very close to the prediction from the phonon calculations. The heat capacity from the Debye–Gru¨neisen model is lower than those of both phonon and present CALPHAD data at high temperatures, but agrees very well at low temperatures. Since the only data from the Debye– Gru¨neisen model to be used in the CALPHAD model is the entropy of formation at 300 K, the Debye–Gru¨neisen model predictions should be suitable for the CALPHAD model. The enthalpy and entropy of formation of b-boron endmembers, 1289 J mol1 and 0.95 J mol1 K1, respectively, with the B101(B,C)4 sublattice model are obtained in the same way. Furthermore, the enthalpy of mixing of the carbon in the structure is obtained through first-principles calculations by replacing one or three of the B atoms in the B12 icosahedra with C, resulting in 432 and 960 J mol1, respectively.

A. Saengdeejing et al. / Acta Materialia 60 (2012) 7207–7215 Table 1 Gibbs energy of end-members of the phases in the B–C system. Phases

Endmembers

Functions

Liquid [35] 0

GLiquid B

GLiquid B

0

GLiquid C

GLiquid C

Graphite [35] 0

GGraphite B

GGraphite B

0

GGraphite C

GGraphite C

b-boron [35] 0

Gbboron B101 B4

105Gbboron B

0

Gbboron B101 C 4

101Gbboron þ 4GGraphite  135340  11:84T B C

0

C GBB4þx 11 BBBB

15Gbboron þ 265910  15:119T B

B4+xC

C GBB4þx 11 BBBC 0 B4þx C GB11 BBBVa 0 B4þx C GB11 BBVaB 0 B4þx C GB11 BBVaC 0 B4þx C GB11 BBVaVa 0 B4þx C GB11 BCBB 0 B4þx C GB11 BCBC 0 B4þx C GB11 BCBVa 0 B4þx C GB11 BCVaB 0 B4þx C GB11 BCVaC 0 B4þx C GB11 BCVaVa 0 B4þx C GB11 BVaBB 0 B4þx C GB11 BVaBC 0 B4þx C GB11 BVaBVa 0 B4þx C GB11 BVaVaB 0 B4þx C GB11 BVaVaC 0 B4þx C GB11 BVaVaVa 0 B4þx C GB11 CBBB 0 B4þx C GB11 CBBC 0 B4þx C GB11 CBBVa 0 B4þx C GB11 CBVaB 0 B4þx C GB11 CBVaC 0 B4þx C GB11 CBVaVa 0 B4þx C GB11 CCBB 0 B4þx C GB11 CCBC 0 B4þx C GB11 CCBVa 0 B4þx C GB11 CCVaB 0 B4þx C GB11 CCVaC 0 B4þx C GB11 CCVaVa 0 B4þx C GB11 CVaBB 0 B4þx C GB11 CVaBC 0 B4þx C GB11 CVaBVa 0 B4þx C GB11 CVaVaB 0 B4þx C GB11 CVaVaC 0 B4þx C GB11 CVaVaVa 0

14Gbboron þ GGraphite  12237  6:419T B C 14Gbboron þ 257344  14:111T B 14Gbboron þ 140063  6:969T B 13Gbboron þ GGraphite þ 60140  2:728T B C 13Gbboron þ 244469 þ 7:345T B 14Gbboron þ GGraphite  12237  6:419T B C 13Gbboron þ 2GGraphite  147927 þ 2:371T B C 13Gbboron þ GGraphite þ 58121 þ 6:167T B C 13Gbboron þ GGraphite þ 60140  2:728T B C

7211

Table 2 ˚ 3 atom1), bulk modulus, B (GPa), first Equilibrium volume, V0 (A derivative of bulk modulus to the pressure, B0P , calculated from first principles by fitting to the four-parameter Birch–Murnaghan equation of state. End-members

V0

B

B0P

BBBB BBBC BBBV BBVB BCVB BBVV BCBC BCBV BCVC BCVV BVBV BVVV CBBB CBBC CBBV CBVB CCVB CBVV CCBC CCBV CCVC CCVV CVBV CVVV

7.7508 7.3112 7.7508 8.2770 7.2737 7.4694 7.4209 8.0140 7.7872 7.9076 7.1150 7.2576 7.3001 7.3153 8.0132 8.0133 8.0818 7.2737 7.2582 7.7534 7.4852 7.5960 8.0814 6.6273

174 197 174 186 234 216 198 196 206 172 130 201 202 210 192 192 149 234 203 207 221 178 149 202

3.4107 3.6043 3.4107 3.5980 3.4326 3.6436 3.8186 3.4365 3.3045 3.7948 6.1291 4.0648 3.9324 3.4055 3.3896 3.3845 3.4511 3.4326 4.0898 3.5264 3.5986 3.7316 3.4594 3.5214

12Gbboron þ 2GGraphite  74481 þ 0:888T B C 12Gbboron þ GGraphite þ 69297  14:854T B C

5.2. CALPHAD

14Gbboron þ 257344  14:111T B 13Gbboron þ GGraphite þ 58121 þ 6:167T B C 13Gbboron þ 527223  69:534T B 13Gbboron þ 244469 þ 7:345T B 12Gbboron þ GGraphite þ 69297  14:854T B C 12Gbboron B 14Gbboron þ GGraphite þ 233739  5:932T B C 13Gbboron þ 2GGraphite þ 50656 þ 0:041T B C 13Gbboron þ GGraphite þ 131253  4:182T B C 13Gbboron þ GGraphite þ 131247  4:197T B C 12Gbboron þ 2GGraphite þ 56858 þ 0:097T B C 12Gbboron þ GGraphite þ 237615  24:937T B C 13Gbboron þ 2GGraphite þ 50656 þ 0:041T B C 12Gbboron þ 3GGraphite  187507 þ 8:824T B C 12Gbboron þ GGraphite þ 43261  5:808T B C 12Gbboron þ 2GGraphite þ 56858 þ 0:097T B C 11Gbboron þ 3GGraphite  39106  4:027T B C 11Gbboron þ 2GGraphite þ 155584  13:199T B C 13Gbboron þ GGraphite þ 131253  4:182T B C 12Gbboron þ GGraphite þ 43261  5:808T B C 12Gbboron þ GGraphite þ 237587  25:015T B C 12Gbboron þ GGraphite þ 237615  24:937T B C 11Gbboron þ GGraphite þ 155584  13:199T B C 11Gbboron þ GGraphite þ 122155  5:778T B C

The thermodynamic descriptions for elemental B and C are taken from the SGTE Pure4 database [35]. The B4+xC sublattice model has already been discussed in the previous section along with the enthalpy and entropy of formation for all 36 end-members. Table 1 lists the molar Gibbs energy for all the end-members, with their enthalpies and entropies of formation taken directly from the first-principles calculations. The excess Gibbs energy in Eq. (3) remains to be evaluated with a total of 630 possible regular interaction parameters. Taking into account the symmetry of the three-atom chain and limiting the choice to single sublattice mixing, the number of such parameters decreases to 108. The interaction parameters are chosen from experimental observations on the site occupancy of B4+xC [3,19,23,36– 40]. Accordingly, seven interaction parameters, listed in Table 4, are selected for optimization of the B4+xC Gibbs energy description. B C The first two interaction parameters, 0 LB4þx and 11 BðB;CÞBC 0 B4þx C LB11 BCBðB;CÞ , represent the mixing between boron and carbon in the chain end positions and have the same value. This mixing affects the transition from the BBBC/BCBB structure at the B-rich end of the B4+xC to the higher carbon concentrations by adding more carbon to the chainB C B C end positions. 0 LB4þx and 1 LB4þx are the regular 11 ðB;CÞCBC 11 ðB; CÞCBC and subregular interaction parameters, respectively, that

7212

A. Saengdeejing et al. / Acta Materialia 60 (2012) 7207–7215

Table 3 Scaling parameter (s) for the Debye–Gru¨neisen model that reproduced the Debye temperature (Hph) from phonon calculations for various structures. The last two columns are the calculated Debye temperature from the Debye–Gru¨neisen model using a-boron scaling parameter (HD) and the difference (DH) between HD and Hph, Gru¨neisen parameter at equilibrium volume (c0). s

Hph

HD

jDHj

c0

a-Boron CCBC BCBC BBBC BCVC

0.969 0.989 0.930 1.000 0.981

1372 1461 1349 1369 1390

1372 1435 1390 1326 1366

0 26 41 43 24

2.2809 1.4508 1.6551 1.6355 1.4856

Heat capacity, J/mol-atom-K

Structure

30 25 20

and vacancy in the chain center, reflecting a possible short-range ordering tendency. For the solubility of carbon in b-boron, the first-principles data from the previous section are used for the evaluation of the B101(B,C)4 sublattice model parameters. For boron solubility in graphite, the SGTE Pure4 database [35] contains the Gibbs energy description of boron in the graphite structure, and the interaction parameters are evaluated using the experimental data of the solubility limit of boron in graphite from Lowell [18]. The interaction parameters in liquid are evaluated using the experimental liquidus, peritectic, and eutectic reactions from Elliot [1]. The model parameters of the B–C system are shown in Tables 1 and 4.

Phonon Debye CALPHAD Exp.

15 10 5 0

0

200

400

600 800 1000 1200 1400 Temperature, K

Fig. 2. Heat capacity of CCBC end-member calculated from quasiharmonic approximation, Debye–Gru¨neisen model and experimental data [13].

Table 4 List of interaction parameters of the B–C system. Phase

Parameters

Functions

0 Liquid LðB;CÞ 1 Liquid LðB;CÞ

57102 + 1.181T

0 Graphite LðB;CÞ

12344 + 18.185T

0 bboron LB101 ðB;CÞ 4

48353  406.690T

C 0 B4þx C LB11 BðB;CÞBC and 0 LBB4þx 11 BCBðB;CÞ 0 B4þx C LB11 ðB;CÞCBC 1 B4þx C LB11 ðB;CÞCBC 0 B4þx 4C LB11 CCðB;VaÞC 1 B4þx C LB11 CCðB;VaÞC 0 B4þx C LB11 ðB;CÞCðB;VaÞC

517218 + 18.533T

Liquid 21188  0.718T

Graphite

b-Boron

B4+xC

413746 + 11.887T 114798 510458  0.530T 84694 + 21.902T 331138

denote the boron and carbon mixing in the icosahedra B C B C position. 0 LB4þx and 1 LB4þx are the interaction 11 CCðB;VaÞC 11 CCðB;VaÞC parameters for the boron and vacancy interaction in the chain center. The last interaction parameter, 0 B4þx C LB11 ðB;CÞCðB;VaÞC , depicts the simultaneous mixing between boron and carbon in the icosahedra position and boron

Fig. 3. Boron activity at 2300 K (a) and 637 K (b) in the B4+xC region calculated from this work compared with the Kasper et al. [8] model (dotted lines) and experimental data by Froment and co-workers [11,12].

A. Saengdeejing et al. / Acta Materialia 60 (2012) 7207–7215

Fig. 4. Calculated B–C phase diagram with experimental data [1,16,18] and from the model parameters reported by Kasper et al. [8] (dotted lines).

Fig. 3 compares the boron activity calculated at 2300 and 673 K from this work, the CALPHAD model from Kasper et al. [8] and experiment [11], respectively. It can be seen that the calculated boron activities from the present work agree with the experimental data much better than those by Kasper et al. [8]. The calculated phase diagram agrees with most of the experimental data [1,18] except the maximum solubility limit of the carbon in B4+xC as shown in Fig. 4. Table 5 compares the three invariant reactions calculated from the model with experiments. The carbon solubility limit in the C-rich B4+xC has been reported by four authors [1,14–16], but the results show considerable scatter. Elliott [1] reported values of 19.8 at.% C using X-ray diffraction. Bouchacourt and Thevenot [14] gave a value of 20 at.% C using the same technique. Beauvy’s value [15] is the highest at 21.6 ± 0.8 at.% C. The most recent value, determined using electron probe microanalysis by Schwetz and Karduck [16], is 18.8 at.% C. The maximum carbon solubility predicted from the present model is 18.7 at.% C at 2300 K and agrees with the experimental determination from Schwetz and Karduck [16]. The predicted minimum carbon solubility in B4+xC agrees very well with the experimental data. The model accurately reproduces the solubility of boron in graphite. The solubility limit of the carbon in b-boron is unknown experimentally. Werheit et al. [25] showed that

7213

Fig. 5. Calculated B–C phase diagram on the B-rich side with experimental data [1,16,17] and from the model parameters reported by Kasper et al. [8] (dotted lines).

carbon is still soluble in b-boron at 1 at.% C. According to the data from first-principles calculations, the temperature-independent interaction parameter is too small, resulting in almost no carbon solubility in b-boron phase. The sublattice model for b-boron phase, B93(B,C)12, in Kasper’s model [8] has the value of interaction parameter about 1637 kJ, compared to 48 kJ from this model. The larger negative value of interaction parameter results in the solubility of carbon in the b-boron phase. First-principles calculations on different carbon configurations in b-boron structure have been performed. The results are still not be able to reproduce the carbon solubility in b-boron phase. A fairly large negative temperature-dependent term is used to achieve the carbon solubility in b-boron phase. It is not common to have such a negative value in the temperature-dependent term in the interaction parameter because the temperature-independent term is fixed, but for the data from first-principles the negative temperature-dependent term is required to achieve the carbon solubility in bboron. The B-rich side of the phase diagram is shown in Fig. 5. The peritectic temperature is about 26 K higher than experiment and the composition of the congruent melting of B4+xC is 3 K lower than the value reported by Beauvy [15].

Table 5 Calculated invariant and congruent reactions compared with experiments (in parenthesis). Reactions

Liquid B4+xC Liquid B4+xC + Graphite Liquid + B4+xC b-boron

at.%C Liquid

B4+xC

– 29.0 (29.0) 0.09 (0.10)

– 18.8 (19.9) 8.9 (9.5)

Temperature, K

Ref.

2720 (2723) 2648 (2648) 2408 (2348)

[15] [1] [1]

7214

A. Saengdeejing et al. / Acta Materialia 60 (2012) 7207–7215

Fig. 6. Predicted site occupancy of the B11(B,C) (B,C,Va) (B,Va) (B,C,Va) sublattice at 2300 K. The letter in the legend refers to the elements in the sublattice, while the number refers to the sublattice.

As the carbon concentration increases, carbon begins to enter the icosahedral sites while the carbon occupancy in the chain ends continues to increase. When the carbon concentration reaches around 14 at.%, a new defect mechanism emerges with some of the chain center sites becoming vacant. The carbon occupancy in the chain ends reaches 100% with the carbon concentration at about 18 at.%, while the occupancy of carbon in the icosahedron site and the vacancy in the chain center continue to increase as the carbon concentration increases, reaching 85% and 30%, respectively, at the C-rich solubility limit. The presence of vacancies in the chain center positions was discussed by Kwei et al. [3], where data suggested that a large fraction of boron in the chain center can be replaced by vacancies, up to 17% based on neutron powder diffraction data. Fig. 7 shows the calculated site occupancy of vacancies in the chain center together with experimental data from Kwei et al. The possible explanation of the discrepancy between experimental and predicted value might come from the temperature effect. The experimental value is obtained from the sample at room temperature while our predicted value is calculated at high temperatures. The predicted vacancy concentration at the chain center position also shows the trend of decreasing as the temperature is decreased. At the C-rich end of B4+xC, the structure is dominated by the B11C icosahedral structure while the inter-icosahedral chain comprises a mixture of C–B–C and C–Va–C chains. 6. Summary

Fig. 7. Predicted B4+xC site occupancy of vacancies in the chain center position as a function of temperature compared with experimental data [3].

Fig. 6 presents the site occupancies of various species in B4+xC at 2300 K as a function of carbon concentration. It is observed that at low carbon concentrations, carbon mainly enters the ends of the chain, i.e. the 3rd and 5th sublattices, with the icosahedron and chain center sites fully occupied by boron. This prediction agrees with the experimental observations by Raman spectroscopy [41] and acoustic measurements [42] that the C–B–B or C–B– C inter-icosahedra chain without carbon in the icosahedra is the main structure at this composition.

The energetics of defects in B4+xC boron carbide and bboron are systematically studied with their enthalpy and entropy of formation predicted from first-principles calculations, resulting in a two-sublattice model, B101(B,C)4, for b-boron, and a five-sublattice model, B11(B,C) (B,C,Va) (B,Va) (B,C,Va), for B4+xC. These energetics enable us to accurately model the thermodynamic properties of these two phases along with liquid and graphite in the B–C system using the CALPHAD approach. The calculated phase diagram reproduces experimental data very well. The calculated boron activities in B4+xC at high temperatures show better agreement with experimental data than previous modeling in the literature. The defect structure evolution in B4+xC with carbon concentrations are predicted from the model and show good agreement with the observations reported in the literature. Acknowledgements This work is funded by the National Science Foundation (NSF) through Focused Research Group (FRG) Grant DMR-0514592 led by Prof. David Larbalestier. First-principles calculations were carried out on the CyberSTAR and the LION clusters at the Pennsylvania State University supported in part by the Materials Simulation Center and the Graduate Education and Research Services

A. Saengdeejing et al. / Acta Materialia 60 (2012) 7207–7215

at the Pennsylvania State University. We would also like to thank Dr. Shunli Shang, Dr. Dongwon Shin, Dr. Raymundo Arroyave, and Ms. Sara Prins for their helps in firstprinciples calculations and CALPHAD modeling and stimulating discussions. References [1] Elliott RP. The boron-carbon system. Tech. rep. U.S. at. energy comm. Chicago: Illinois Inst. of Technol; 1961. [2] Thevenot F. Boron carbide – a comprehensive review. J Eur Ceram Soc 1990;6(4):205–25. [3] Kwei GH, Morosin B. Structures of the boron-rich boron carbides from neutron powder diffraction: implications for the nature of the inter-icosahedral chains. J Phys Chem 1996;100(19):8031–9. [4] Zhdanov GS, Sevast’yanov NG. Crystal structure of boron carbide (b4c). Dokl Akad Nauk SSSR, Ser A 1941;32:432–4. [5] Clark HK, Hoard JL. The crystal structure of boron carbide. J Am Chem Soc 1943;65:2115–9. [6] Lundstroem T, Bolmgren H. On refractory boron compounds crystallizing in b4c-type structures. Jpn J Appl Phys Ser 1994;10:1–4. [7] Kaufman L, Uhrenius B, Birnie D, Taylor K. Coupled pair potential, thermochemical and phase diagram data for transition metal binary systems-vii. CALPHAD: Comput Coupl Phase Diagr Thermochem 1984;8(1):25–66. [8] Kasper B, Seifert HJ, Kussmaul A, Lukas HL, Aldinger F. Development of a thermodynamic data set for the system B–C–N– O–Si. In: Materialwissenschaftliche Grundlagen, symposium 7, Werkstoffwoche ’96, Stuttgart 1996, 1997. p. 623–8. [9] Saal JE, Shang S, Liu ZK. The structural evolution of boron carbide via ab-initio calculations. Appl Phys Lett 2007;91(23). [10] Liu ZK. First-principles calculations and CALPHAD modeling of thermodynamics. J Phase Equilib Diffus 2009;30(5):517–34. [11] Froment AK, Chatillon C, Colin M. Determination of boron and carbon activities in the non-stoichiometric b4xc compound by the multiple knudsen-cell mass spectrometric method. Rev Int Hautes Temp 1991;27(4):141–58. [12] Froment AK, Chatillon C, Fouletier J, Fouletier M. Carbon activity measurements in boron carbides using a solid state potentiometric cell. J Nucl Mater 1992;188:280–4. [13] Matsui T, Arita Y, Naito K, Imai H. High temperature heat capacities and electrical conductivities of boron carbides. J Nucl Mater 1991;186(1):7–12. [14] Bouchacourt M, Thevenot F. Analytical investigations in the boron– carbon system. J Less-Common Met 1981;82:219–26. [15] Beauvy M. Stoichiometric limits of carbon-rich boron carbide phases. J Less-Common Met 1983;90(2):169–75. [16] Schwetz KA, Karduck P. Investigations in the boron–carbon system with the aid of electron probe microanalysis. J Less-Common Met 1991;175(1):1–11. [17] Kieffer R, Gugel E, Leimer G, Ettmayer P. Boron–carbon system. Ber Dtsch Keram Ges 1971;48(9):385–9. [18] Lowell CE. Solid solution of boron in graphite. J Am Ceram Soc 1967;50(3):142–4. [19] Kuhlmann U, Werheit H, Schwetz KA. Distribution of carbon atoms on the boron carbide structure elements. J Alloys Compd 1962;189(2):249–58.

7215

[20] Shang SL, Wang Y, Arroyave R, Liu ZK. Phase stability in alphaand beta-rhombohedral boron. Phys Rev B 2007;75(9):092101. [21] Widom M, Mihalkovicˇ M. Symmetry-broken crystal structure of elemental boron at low temperature. Phys Rev B 2008;77:064113. [22] Redlich O, Kister AT. Algebraic representation of thermodynamic properties and the classification of solutions. Ind Eng Chem 1948;40(2):345–8. [23] Emin D. Structure and single-phase regime of boron carbides. Phys Rev B 1988;38(9):6041–55. [24] Naslain R. Crystal chemistry of boron and of some boron-rich phases; preparation of boron modifications. Boron Refract. Borides. New York (NY): Springer-Verlag; 1977. [25] Werheit H, Kuhlmann U, Laux M, Lundstroem T. Structural and electronic properties of carbon-doped b-rhombohedral boron. Phys Status Solidi B 1993;179(2):489–511. [26] Blochl PE. Projector augmented-wave method. Phys Rev B 1994;50(24):17953–79. [27] Perdew JP, Wang Y. Accurate and simple analytic representation of the electron-gas correlation-energy. Phys Rev B 1992;45(23):13244–9. [28] Kresse, G. VASP, Vienna ab initio simulation package, 2004. . [29] Kresse G, Joubert D. From ultrasoft pseudopotentials to the projector augmented-wave method. Phys Rev B 1999;59(3):1758. [30] Methfessel M, Paxton AT. High-precision sampling for brillouin-zone integration in metals. Phys Rev B 1989;40(6):3616–21. [31] Blochl PE, Jepsen O, Andersen OK. Improved tetrahedron method for brillouin-zone integrations. Phys Rev B 1994;49(23):16223–33. [32] Shang SL, Wang Y, Kim D, Liu ZK. First-principles thermodynamics from phonon and Debye model: application to Ni and Ni3Al. Comput Mater Sci 2010;47(4):1040–8. [33] Wei S, Chou MY. Ab-initio calculation of force constants and full phonon dispersions. Phys Rev Lett 1992;69(19):2799–802. [34] Moruzzi VL, Janak JF, Schwarz K. Calculated thermal-properties of metals. Phys Rev B 1988;37(2):790–9. [35] Dinsdale AT. SGTE data for pure elements. CALPHAD: Comput Coupl Phase Diagr Thermochem 1991;15(4):317–425. [36] Bylander DM, Kleinman L. Structure of boron carbide (b13c2). Phys Rev B 1991;43(2):1487–91. [37] Aselage TL, Tissot RG. Lattice constants of boron carbides. J Am Ceram Soc 1992;75(8):2207–12. [38] Lazzari R, Vast N, Besson JM, Baroni S, Dal Corso A. Atomic structure and vibrational properties of icosahedral b4c boron carbide. Phys Rev Lett 1999;83(16):3230–3. [39] Mauri F, Vast N, Pickard CJ. Atomic structure of icosahedral b4c boron carbide from a first principles analysis of NMR spectra. Phys Rev Lett 2001;87(8):085506/1–4. [40] Feng Y, Seidler GT, Cross JO, Macrander AT, Rehr JJ. Role of inversion symmetry and multipole effects in nonresonant X-ray raman scattering from icosahedral b4c. Phys Rev B 2004;69(12):125402/1–8. [41] Tallant DR, Aselage TL, Campbell AN, Emin D. Boron carbide structure by Raman spectroscopy. Phys Rev B 1989;40(8):5649–56. [42] Gieske JH, Aselage TL, Emin D. Elastic properties of boron carbides. AIP Conf Proc 1991;231:376–9 [Boron-Rich Solids]. [43] Momma K, Izumi F. VESTA3 for three-dimensional visualization of crystal, volumetric and morphology data. J Appl Crystallogr 2011;44(6):1272–6.