J. Phys. Chem. B 2005, 109, 14627-14631
14627
Thermal Vacancy Effects on the Thermodynamic Properties and Stability of Fullerites Viatcheslav I. Zubov* and Igor V. Zubov Department of Theoretical Physics, Peoples’ Friendship UniVersity, 117198 Moscow, Russia ReceiVed: March 4, 2005; In Final Form: May 3, 2005
We study the influence of thermal vacancies on equilibrium thermodynamic properties of the high-temperature phase of fullerites taking into account the strong anharmonicity of the lattice vibrations. Treating a crystal with point defects as a quasi-multicomponent system and using the correlative method of the unsymmetrized self-consistent field and the Girifalco interaction potential for the molecular subsystem, we have obtained the vacancy formation parameters for the C60 fullerite. We also take into account the divacancies. The influence of the lattice defects on the specific heats of fullerites is negligible, since a dominant contribution to them is given by intramolecular degrees of freedom. We have calculated the contributions of vacancies to various thermodynamic properties, the volume thermal expansion coefficient, the isothermal bulk modulus, and the components of the isothermal elastic tensor, depending on the temperature and pressure. Near the estimated triple point, these contributions run to more than 10% and increase still further at a metastable region. We also discuss the influence of defects on the thermodynamic stability of fullerites.
1. Introduction Fullerenes as the third molecular form of carbon were discovered in 1985.1 Since the elaboration of effective methods for fabrication, separation, and deep purification of these molecules C20+2n,2 various physical and chemical properties of fullerenes and fullerites, i.e., crystals formed by such molecules, have been investigated (structural, electrical, magnetic, and optical properties, phase transitions, polymerizarion, etc.).3-12 The most widespread fullerite, C60, and the next to it, C70, are the most completely studied. In particular, at low temperatures the molecules in their lattices are orientationally ordered, whereas at high temperatures they rotate rather freely in a facecentered cubic (fcc) lattice (with a little of the hexagonal close packed (hcp) phase in C70). Similar behavior is to be expected also from other fullerites. The C60 fullerite doped with some alkali atoms (fulleride) at temperatures in the range 20-35 K exhibits superconductivity. On the basis of the correlative method of the unsymmetrized self-consistent field (CUSF)13-17 and the Girifalco intermolecular potential,3,18-20 we study the thermodynamic properties of the high-temperature modification of fullerites, considering them as perfect crystals. In our previous paper,21 using this approach we have investigated the temperature dependence of the saturated vapor pressures of the family of fullerites from C36 to C96, a complete set of their equilibruim elastic and thermal properties along the sublimation curves, their thermodynamic stability, and its loss. However, it is known22-28 that at high temperatures thermodynamic properties of crystals are susceptible to the influence of point defects, particularly to that of thermal vacancies. Although at low temperatures, the quantity of such defects is vanishingly small, their concentration and especially their derivatives increase rapidly with the temperature, and near the triple point they give considerable contributions to the properties expressed in terms of the second derivatives of the free energy. * Author to whom correspondence should be addressed. E-mail:
[email protected].
In the present paper, we calculate the parameters of vacancy formation, their concentration, and their effects on equilibrium thermodynamic properties of the C60 fullerite. 2. Basic Formulas: Parameters of the Vacancy Formation We consider a pure crystal with short-range interactions containing point defects at the temperature T ) Θ/k and pressure P. Its equilibrium state is governed by the Gibbs free energy
G(Θ,P,N;{Ni}) ) G0 +
∑i Nigfi - Θ ln W
(1)
Here, N is the number of atoms or molecules, G0(Θ,P,N) is the free energy of a hypothetical perfect crystal, Ni (,N) is the number of defects of the ith kind (vacancies and their complexes), W(N,{Ni}) is the statistical weight, i.e., the number of ways the {Ni} defects can be placed in the lattice of
L)N+
∑ νiNi
(2)
sites where νi is the number of lattice points entered into the ith defect, and gfi is the Gibbs energy change on forming a single ith defect (without the entropy of mixing) or the free energy of a defect formation. It defines the entropy of the formation
sfi )
() ∂gfi ∂T
P
(3)
its enthalpy
hfi ) gfi + Tsfi
(4)
volume
Vfi
()
∂gfi ) ∂P
10.1021/jp0511303 CCC: $30.25 © 2005 American Chemical Society Published on Web 07/06/2005
T
(5)
14628 J. Phys. Chem. B, Vol. 109, No. 30, 2005
Zubov and Zubov
and its other parameters. Minimizing G with respect to Ni at constant N, Θ, and P give the equilibrium values of Ni. Considering only monovacancies in crystals with short-range interactions and utilizing the concept of the binary (quasi-binary) mixtures29,30 enables one to express the free energy of vacancy formation in terms of that of the hypothetical perfect crystal.31 Here, we reproduce the main formulas for the case of the pairwise central forces. The Gibbs free energy and the enthalpy of vacancy formation are of the form
gf1
) -f ′ + Θ +
P(2Vf1
- V)
(7)
where V is the volume of the unit cell and f′(Θ,V) and uj are the nonkinetic parts of the specific (per atom or molecule) latticevibration Helmholtz free energy and internal energy, respectively. In the CUSF for strongly anharmonic crystals,13-17 f′ is written as
f ′ ) f ′0 + f2 + fH
(8)
Here, f′0 is the zeroth approximation that takes into consideration the strong anharmonicity whereas f 2 and fH are corrections of the perturbation theory that refine the influence of the intermolecular correlations and contain the anharmonicity of higher orders. In its turn
uj ) uj0 - Θ2
(
∂ f +f ∂Θ Θ
)
H
(9)
Including into the zeroth approximation terms up to the fourth order, we have
f′0 )
() (
[( )
)
K0 5Θ β 2 Θ 5β 2 1 12π2Θ - X+ - Θ ln 2 24 X 4 6X V K4
uj0 )
3/4
D1.5
5β (X + 6X )] (10)
K0 (3 + β)Θ + 2 4
(11)
where
K2l )
1
∑ Zk∇2lΦ(Rk)
2l + 1 kg1
l ) 0, 1, 2
(12)
where Φ(r) is the intermolecular potential, Zk are the coordinational numbers, Rk are the radii of the coordinational spheres, β(K2x3/ΘK4) is the solution of the transcendental equation
D-2.5(X + 5β/6X) β ) 3X D-1.5(X + 5β/6X)
(13)
in which Dν are the parabolic cylinder functions. The perturbation theory takes into account the anharmonic terms up to the sixth order. In eqs 8 and 9, we do not take into consideration the quantum corrections because they are negligible for fullerites. Substituting eq 6 into eq 5, we have the differential equation for the vacancy formation volume
( )
∂Vf1 2P ∂P
Θ
BT ) - V
(∂P∂V)
Θ
(15)
is the isothermal bulk modulus of the hypothetical perfect crystal, BpT. Due to the small compressibility of crystals, the approximate solution of eq 12 is
(
Vf1 ≈ V 1 +
(6)
hfV ) - u + P(VfV - V)
2
where
)
2P Θ - p 3BpT BT
(16)
The entropy of the vacancy formation can be calculated from eq 3 or using the thermodynamic relation
sf1 ) k
hf1 - gf1 Θ
(17)
The influence of the monovacancies on the thermodynamic properties of crystals in terms of hfV (or of gfV) and its derivatives. In particular, its effects on the volume thermal expansion coefficient and the isothermal compressibility κT ) 1/BT are given by the formulas24,27
c1 Vf1 hf1 R≈R + V kT2
(18)
c1(Vf1)2 V kT
(19)
p
κT ≈ κpT +
where c1 is the vacancy concentration that in this approximation is
c1 )
N1 N + N1
) e-g1f/Θ
(20)
Here, we do not write the expression for the specific heats because a dominant contribution to them is given by the intramolecular degrees of freedom. As for their lattice parts,31 the vacancies have a pronounced effect on CP but not on CV, and hence we can say that isochoric processes are scarcely affected by vacancies. The calculation of vacancy contributions to the components of elastic tensors CFσ is somewhat more complicated. Here, we take into account that the shear moduli C44 and C11-C12 correspond to isochoric processes. Therefore, contributions of vacancies to them are expected, as well as that to the isochoric specific heat, to be negligible. As this takes place, the contribution to C12 is equal to that to C11. Next, it is evident from the known relation BT ) (CT11 + 2CT12)/3 that it is equal to the contribution to BT T T,p - C1β ) BT - BpT β ) 1, 2 C1β
(21)
At high temperatures, some vacancies unite into bonded complexes. This causes the concentration of vacant lattice points and their contributions to thermodynamic properties to increase in the vicinity of the triple point of a crystal and especially at its metastable superheated state (between the triple and spinodal points). 3. Monovacancies and Divacancies
+
Vf1
Θ )VBT
(14)
Let us return to eq 1. In calculating W, one needs to take into account that two defects cannot be nearest neighbors; this
Thermal Vacancy Effects on Fullerites
J. Phys. Chem. B, Vol. 109, No. 30, 2005 14629
is a new defect. Such a calculation can be performed using the modified Lidiard-Howard method.26 We shall take into consideration only vacancies (ν1 ) 1) and divacancies (ν2 ) 2), then N1-1
[L - (Z + 1)l] ∏ (Z/2)N l)0
2
W(N,N1,N2) )
N 1!
N2!
N2-1
∏ (L - 2m)
(22)
m)0
where Z is the first coordination number. Minimizing eq 1 with respect to N1 and N2 and taking into account eq 2 gives the equilibrium concentration of vacancies and divacancies
N1
c1 )
≈ e-g1 /Θ(1 - z e-g1 /Θ) f
N + N1 + 2N2
f
(
(23) Figure 1. Gibbs free energy (1, 2, 3) and enthalpy (4, 5, 6) of the vacancy formation in C60 fullerite in equilibrium with the vapor (1, 4), at 1 kbar (2, 5), and at 2 kbar (3, 6).
)
2gf2 gb2 f z z + c2 ≈ e-g2 /Θ ) exp 2 2 Θ Θ
(24)
where ) is the divacancy binding Gibbs free energy. The total concentration of the vacant lattice points is gb2
2gf1
gf2
ct ) c1 + 2c2 ≈ e-g1 /Θ[1 + Z e-g1 /Θ(e-g2 /Θ - 1)] (25) f
f
b
) 0 this formula gives the It is easily seen that when vacancy concentration in the absence of the interactions between them, eq 20. One can evaluate the divacancy binding free energy taking into consideration that each monovacancy is surrounded by Z nearest neighbors with one broken bond whereas a divacancy has 2Z - 1 such neighbors. Therefore, for crystals with shortrange forces, it is reasonable to suppose that gb2
gb2 ≈ gf1/Z
(26)
and hence
ct ≈ e-g1 /Θ[1 + Z e-g1 /Θ(e-g1 /ZΘ - 1)] f
f
f
(27)
The thermal expansion coefficient and isothermal bulk modulus of a crystal with vacancies and divacancies can be obtained from eqs 18 and 19 by replacing c1 with ct. 4. Outline of Calculations and Results At the low-temperature, orientationally ordered phase of the fullerites, the vacancy concentration and its influence on the thermodynamic properties are negligible. We investigate their high-temperature phases in which the molecules rotate rather freely. Intermolecular forces at orientationally disordered phases of fullerenes are described by the potential
ΦG(r) ) -R
(
)
1 2 1 + + s(s - 1)3 s(s + 1)3 s4 1 1 1 + (28) β s(s - 1)9 s(s + 1)9 s10
(
vacancies in the most abundant fullerite, C60. In this case, a ) 3.55 × 10-8 cm, R ) 7.494 × 10-14 erg, and β ) 1.3595 × 10-16 erg.3 At first, we have solved the equation of state
)
where s ) r/2a, r is the distance between the centers of the molecules, and R is the radius of their hard core. It has been proposed by Girifalco for C603 and then generalized for the family of fullerenes from C28 to C96.18-20 Here, we study the
P)-
[
]
a 1 dK0 βΘ dK2 (3 - β)Θ dK4 + P2 + PH + + 3V 2 da 2K2 da 4K4 da (29)
(in which a is the lattice parameter) within the perfect crystal approximation, eqs 12 and 13, at various temperatures, from equilibrium with the orientationally ordered phase up to the spinodal. In eq 29, P is the fixed external pressure or the saturated vapor pressure; in the latter case, this equation adds to the condition of phase equilibrium with the gas.17 Then, we have calculated thermodynamic properties, including parameters of vacancies, along the lower branches that are thermodynamically stable. We have performed our calculations in equilibrium with the saturated vapor and at 1 and 2 kbar. The volume of the vacancy formation, eq 16, is very nearly equal to that of the unit cell. The entropy of the vacancy formation, eq 17, depends almost not at all on the pressure and decreases with increasing temperature from 8.8 to 9.0 at the equilibrium curve with lowtemperature modification to 5.7 at the estimated32 melting curve and to 4.17-4.50 at the spinodal curve (in units of the Boltzmann constant). Figure 1 demonstrates the Gibbs free energy Gf1 ) Ngf1 and the entalpy Hf1 ) Nhf1 of the vacancy formation (N is the Avogadro number). For equilibrium with the saturated vapor, we show them up to the spinodal point, and in what follows we shall discuss the influence of vacancies on the thermodynamic stability. In other cases, we restrict our consideration to the estimated32 melting points. One can see that both quantities increase with the pressure and decrease with increasing temperature. The dependence of the logarithm of the vacancy concentration on the inverse of temperature is demonstrated in Figure 2. It is close to linear, the slope increasing with the pressure, which reflects the pressure dependence of Gf1. Note that a fraction of the divacancies in the concentration of the vacant lattice points becomes perceptible in proximity to the estimated triple point, reaches the greatest value at the spinodal point in equilibrium with the vapor where it is about 10%, and
14630 J. Phys. Chem. B, Vol. 109, No. 30, 2005
Zubov and Zubov
Figure 4. Thermal expansion coefficient R, in 10-5 K-1, in equilibrium with the vapor (1) and at 1 kbar (2) and the Gru¨neisen parameter γ (3), (a) without and (b) with vacancies. Figure 2. Logarithm of the vacancy concentration, eq 20, versus the inverse of temperature in equilibrium with the vapor (1), at 1 kbar (2), and at 2 kbar (3). Ttr is the estimated triple point; other vertical lines show the spinodal points.
Figure 5. Solution of eq 29 at high temperatures taking into account the vacancies: (1) the stable (and metastable) branch, BT > 0, (2) the absolutely unstable branch, BT < 0, and the vacancy concentration, (3) eq 20 and (4) eq 27.
fullerites, it reverses sign in the metastable region), and therefore, it is meaningless to discuss the relative contribution from it. In Figure 4, we show the thermal expansion coefficient R and the Gru¨neisen parameter γ ) VBTR/ClV, where ClV is the lattice part of the isochoric specific heat. Near the triple point, the contributions of the vacancies to both properties are about 20%, a decrease of the Gru¨neisen parameter changing to an increase.
Figure 3. Elastic properties of the C60 fullerite (a) without and (b) with vacancies, in equilibrium with the vapor (odd numbers) and at 1 kbar (even numbers): CT11 (1, 2), CT12 (3, 4), and BT (5, 6).
decreases rapidly with the increase in pressure; for more details, see below. Figure 3 displays elastic properties of the C60 fullerite, from the equilibrium curve with the low-temperature phase to the estimated melting curve, (a) in the perfect crystal approximation and (b) with consideration for vacancies. Their concentration and effects on thermodynamic properties decrease very steeply with increasing pressure. (The pressure “forces out” the vacancies of a crystal.) Therefore, we do not show here the results for 2 kbar. Near the estimated triple point, the contributions of the vacancies comprise 10% in the modulus of the uniaxial compression CT11 and 25% in the isothermal bulk modulus BT. The component CT12 is not positively definite (in the case of
In the vicinity of the spinodal points TS(P), the vacancy concentration c1(T) as well as ct(T) are free of a singularity, as contrasted with the isothermal compressibility κT and the thermal expansion coefficient R, which tend to infinity. Consequently, the loss of stability of the solid phase is caused by the lattice vibrations rather than by the vacancies (the instability of the subsystem of the lattice vibration but not of that of the vacancies). Nevertheless, the vacancies have effects on the effective force coefficients, eq 12, and hence reduce the spinodal temperature TS. To study this matter, we have solved also eq 29, changing K2l f (1 - ct)K2l, in equilibrium with the vapor. In such a case, we obtained TS = 1844 K, instead of 1917 K for the perfect C60 crystal;17 i.e., the vacancies reduce the spinodal temperature (of the system crystal-vapor) by about 4%. In Figure 5, we show the high-temperature part of the temperature dependence of the lattice parameter and the vacancy concentration. The latter increases from 2 × 10-3 at the estimated triple point to 2.15 × 10-2, where about 10% of vacancies unite into divacancies.
Thermal Vacancy Effects on Fullerites 5. Discussion We have finished the investigation of the influence of thermal vacancies on the thermodynamic properties and stability of the C60 fullerite depending on the temperature and the pressure. The impact of the vacancies on the thermodynamic stability is studied as well, especially the effects of the vacancies on the force constants. One can compare the present results with those for simple van der Waals crystals (solidified noble gases)31 and ensure that in C60 the effect is of great intensity. The influence of vacancies on the properties of the pressure decreases very rapidly with increasing pressure. Through the use of the Girifalco potential with corresponding parameters,20 it is easy to calculate the properties for other fullerites. The results can be used for the calculation of other properties of fullerites in which the vacancy contribution is proportional of the vacancy concentration Ct as well. References and Notes (1) Kroto, H. W.; Heath, J. R.; O’Brien, S. C.; Curl, R. F.; Smalley, R. E. Nature 1985, 318, 162. (2) Kra¨tschmer, W.; Lamb, L. D. Fostinoupolos, K.; Huffman, R. D. Nature 1990, 347, 354. (3) Girifalco, L. A. J. Phys. Chem. 1992, 96, 858. (4) Tanigaki, K.; Hirosawa, I.; Ebbesen, T. W.; Mizukishimakawa, J.; Kubo, Y.; Tsai, J. S.; Kuroshima, S. Nature 1992, 356, 419. (5) The Fullerenes; Kroto, H. W., Fischer, J. E., Cox D. E., Eds.; Pergamon Press: Oxford, U.K., New York, 1993. (6) Eletskii, A. V.; Smirnov, B. M. Usp. Fiz. Nauk 1995, 165, 977 (in Russian). (7) Dresselhaus, M. S.; Dresselhaus, G.; Eklund, P. C. Science of Fullerenes and Carbon Nanotubes; Academic: San Diego, CA, 1995. (8) Buntar, V.; Sauerzopf, F. M.; Weber, H. W. Aust. J. Phys. 1997, 50, 329.
J. Phys. Chem. B, Vol. 109, No. 30, 2005 14631 (9) Fullerenes and Related Structures; Hirsh, A., Ed.; Springer-Verlag, Berlin, Germany, 1999. (10) Fartaria R. P. S.; Fernandes, F. M. S. S.; Freitas, F. M. M. J. Phys. Chem. B 2002, 106, 10227. (11) Fernandes, F. M. S. S.; Freitas, F. M. M.; Fartaria R. P. S. J. Phys. Chem. B 2004, 108, 9251. (12) Licenkov, S. V., Chernozatonskii, L. A.; Stankevich, I. V. Fiz. TVerd. Tela 2004, 46, 2238 (in Russian). (13) Zubov, V. I. Phys. Status Solidi B 1978, 87, 385, 88, 43. (14) Yukalov, V. I.; Zubov, V. I. Fortschr. Phys. 1983, 31, 627. (15) Zubov, V. I.; Sanchez, J. F.; Tretiakov; N. P.; Yusef, A. E. Int. J. Mod. Phys. B 1995, 9, 803. (16) Zubov, V. I.; Tretiakov; N. P.; Sanchez Ortiz, J. F.; Caparica, A. A. Phys. ReV. B 1996, 53, 12080. (17) Zubov, V. I.; Sanchez Ortiz, J. F.; Teixeira Rabelo, J. N.; Zubov, I. V. Phys. ReV. B 1997, 55, 6747. (18) Kniaz, K.; Girifalco, L. A.; Fischer, J. E. J. Phys. Chem. 1995, 99, 16804. (19) Abramo, M. C.; Caccamo, C. J. Phys. Chem. Solids 1996, 57, 1751. (20) Zubov, V. I. Mol. Mater. 2000, 13, 385. (21) Zubov, V. I.; Zubov, I. V.; Teixeira Rabelo, J. N. J. Phys. Chem. B 2003, 107, 10458. (22) Hill, T. L. Statistical Mechanics; McGraw-Hill, New York, 1960. (23) Korpiun, R.; Coufal, H. J. Phys. Status Solidi A 1971, 6, 187. (24) Shocknecht, W, E,; Simmons, R. O. AIP Conf. Proc. 1971, 3, 169. (25) Aksienov, V. L. Fiz. TVerd. Tela 1972, 14, 1986. (26) Girifalco, L. A. Statistical Physics of Materials; Wiley: New York, 1973. (27) Chadwick, A. V.; Glyde, H. R. Rare Gas Solids; Academic Press: London, New York, 1977; Vol. II, p 1152. (28) Card, D. N.; Jacobs, P. W. M. Mol. Phys. 1977, 34, 1. (29) Magalinskii, V. B. IzV. VUZOV, Fiz. 1976, No. 1, 85. (30) Magalinskii, V. B. IzV. VUZOV, Fiz. 1977, No. 9, 17. (31) Zubov, V. I. Phys. Status Solidi B 1980, 101, 95. (32) Zubov, V. I.; Rodrigues, C. G.; Zubov, I. V. Phys. Status Solidi B 2002, 233, 151.