Molecular dynamics simulation of homogeneous nucleation of KBr cluster

Journal of Physics and Chemistry of Solids 66 (2005) 1732–1738 www.elsevier.com/locate/jpcs

Molecular dynamics simulation of homogeneous nucleation of KBr cluster Xiaolei Zhua,b,*, Kai Chenc b

a Department of Chemistry, Nanjing University of Technology, Nanjing 210009, People’s Republic of China Coordination Chemistry Institute, State Key Laboratory of Coordination Chemistry, Nanjing University, Nanjing 210093, People’s Republic of China c Department of Chemistry, Nanjing Normal of University, Nanjing 210097, People’s Republic of China

Received 8 February 2005; revised 24 May 2005; accepted 27 June 2005

Abstract We perform molecular dynamics simulations to study the homogeneous nucleation in the freezing of molten potassium bromide clusters. The nucleation rates tend to decrease with increasing cluster size and temperature. The solid–liquid interfacial free energy ssl of 42.4– 52.3 mJ/m2 is close to the values predicted by Turnbull’s relation and comparable to the experimental observation by Buckle and Ubbelohde. It is interesting to find that there is no cluster size effect on the critical nucleus size. Critical nucleus sizes inferred from classical nucleation theory are of 6.5–20.7 KCBrK ionic pairs in the temperature range of 400–600 K. The critical nucleus size at bulk MD freezing temperature obtained by extrapolation is about 45 KCBrK ionic pairs, which is comparable to the experimental value of NaCl. q 2005 Elsevier Ltd. All rights reserved. Keywords: A. Inorganic compound; B. Crystal growth; D. Phase transitions

1. Introduction Since the phase transition and nucleation are microscopic processes, the experimental study of these processes is difficult. Despite the obvious importance of phase transition, nucleation, and nonequilibrium growth of the new phases in science and technology, there is not much knowledge about the details of molecular processes in such transformations. On the other hand, computer simulation has become an important and respected research tool in recent years, complementing the transitional approaches of theory and experiment in the areas of physics, chemistry, biology, astronomy, and geology, etc. Computer simulation provides a direct route from the microscopic details of a system to macroscopic properties of experimental interest. A great many more investigations of the phase transitions and the nucleations in clusters, mostly concerning melting and freezing, have been carried out by computer simulation [1–11] than by experiment [12–15]. Although Bartell and * Corresponding author. Tel.: C86 258 358 7507. E-mail address: [email protected] (X. Zhu).

0022-3697/$ - see front matter q 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.jpcs.2005.06.009

Dibble [12–15] have contributed significantly to the experimental study of phase transitions and nucleations in molecular clusters, almost all the information we have about the phase transitions and the nucleations in ionic clusters is from computer simulation and theoretical analysis. The molecular dynamics (MD) simulations of melting and freezing processes on (KBr)32 have been first performed by Rose and Berry [2,3]. Then, the MD studies for NaCl [4], KI [6,7,9], RbCl [8], and NaBr [10] clusters have been reported by Bartell, Huang, Zhu, and Deng. The nucleation rates and interfacial free energies are derived from the MD results and classic nucleation theory. The orders of interfacial free energies for KI and NaBr clusters are similar to Buckle and Ubbelohde’s results [16]. Recently, a largescale MD simulation in supercooled NaCl has been carried out by Koishi, Yasuoka, and Ebisuzaki [11]. In order to estimate the critical nucleus size, the special purpose computer MDGRAPE-2 is used for the MD simulations of large (NaCl)6912 and (NaCl)62,500.Critical nucleus size, nucleation time lag, and nucleation rate are directly calculated from the results of quenching MD simulations without using classic nucleation theory. In the present work, the critical nucleus size of KBr at bulk MD freezing temperature is first estimated based on the MD simulations of small (KBr)108, KBr)256, (KBr)500, and (KBr)864 clusters

X. Zhu, K. Chen / Journal of Physics and Chemistry of Solids 66 (2005) 1732–1738

and classic nucleation relations. It is interesting to find in the MD simulations of KBr clusters that there is no cluster size effect on the critical nucleus size and the critical nucleus size linearly depends on the freezing temperature.

1733

heating is continued at 1060 K to generate 150 melted clusters each with 3000 or more time steps than the previous one. Then, we examine nucleation and crystallization of (KBr)256, (KBr)500, (KBr)864 clusters by quenching a thoroughly melted cluster to a certain temperature (400, 450, 500, 550, or 600 K) for a time period of 160–230 ps.

2. Method

Table 1 Potential parameters of KBr in molecular dynamics simulation

Aij (J/mol) ˚) sij (A

KC–KC

KC–BrK

BrK–BrK

25.628 2.926

20.358 3.179

15.268 3.432

3.1. Thermodynamic properties The various diagnostic criteria are often used to identify phase transitions in literatures [4,7,21,22] including total energies, Lindemann indices, pair correlation functions, diffusion coefficients, Voronoi polyhedra, and bond-order parameter. All these criteria give actually the same melting temperatures. Thus, the melting points of clusters are estimated using total energy curves during slow heating runs (Fig. 1(a)). The melting temperatures of (a) –540 –560 Etot(kJ/mol)

which consists of Coulomb and repulsive interactions. Here, rij is the interionic distance between ions i and j, ˚ [19]. The value of KBr sij Z si C sj , and rZ0.335 A parameters Aij and sij listed in Table 1 are taken from Tosi and Fumi [20]. Clusters with free boundaries are chosen in order to avoid the interference with nucleation caused by periodic boundary conditions. Starting from the initial configuration at the bath temperature of 298 K, the temperature is maintained by rescaling each of the first 10,000 time steps, step size always being 8 fs. After this period of equilibrium, the temperature rescaling is switched off and constant energy MD simulations are performed for 10,000 time steps. Then, the process is repeated. Slow heating stages begin at 400 K, with temperature step of 20 K. Every stage is first simulated at constant temperature for 4000 time steps and then for 6000 time steps at constant energy. Runs are continued until the bath temperatures reach 100 K above the melting temperature of the cluster to ensure complete melting. In order to acquire more accurate estimate of melting temperatures of clusters, additional runs are carried out in the vicinity of the melting transition. The temperature step is reduced to 5 K. In slow cooling processes, each configuration is cooled to 410 K in decrements of 20 K, spend 4000 time steps at constant temperature and 6000 time steps at constant energy. To carry out nucleation analysis, the solid clusters are quickly heated from 400 to 1060 K using 20,000 time steps for (KBr)256, (KBr)500, (KBr)864. For each cluster size,

3. Results and discussion

–580 –600 –620 –640 200

400

600

800

1000

1200

T(K) (b) –540

–560 Etot(kJ/mol)

Molecular dynamics simulations are carried out on 108-, 256-, 500-, and 864- ion pair clusters of potassium bromide using the modified version of the program MDIONS [17]. The initial configuration of a cluster is constructed based on a face˚ centered-cubic (fcc) cell with unit cell parameter of 6.5966 A [18]. The external shape of a cluster is taken as a cube because the cubic cluster of alkali halide has lowest energy [7]. In MD simulation, we use the Born–Mayer–Huggins’ potential functions [19,20] to describe the interionic interactions, i.e.
 zi zj e 2 sij Krij Fij ðrij Þ Z C Aij exp (1) rij r

–580

–600

–620

400

500

600

700

800

900

1000 1100

T (K) Fig. 1. Temperature dependence of the total energy during heating (a) and cooling (b). Open squares for (KBr)108; open triangles for (KBr)256; filled triangles for (KBr)500; open circles for (KBr)864.

X. Zhu, K. Chen / Journal of Physics and Chemistry of Solids 66 (2005) 1732–1738

(a)

-608

Etot (kJ/mol)

-610 -612 -614 -616 -618 -620 -622 -624

0

50

100

150

200

150

200

150

200

time(ps) (b)

0.5 0.4

Q6

(KBr)108, (KBr)256, (KBr)500, and (KBr)864 derived from Fig. 1(a) are 923, 946, 957, 968 K, respectively. It appears that there is an approximately linear relation between the melting points of (KBr)n clusters and the reciprocals of cluster radii. The bulk melting point of 1012(3) K inferred from intercept is in good agreement with the experimental value of 1007 K [23]. As shown in Fig. 1(a), the total energy curves for heating process are obtained from MD simulations. The heats of fusion for (KBr)n clusters can be estimated from the differences between the total energy curves of solid clusters and liquid clusters. The heats of fusion for (KBr)108, (KBr)256, (KBr)500, (KBr)864 are 18.3, 19.4, 21.0, and 21.7 kJ molK1, respectively. When the heats of fusion of (KBr) n are plotted against NK1/3, the intercept (25.0(7) kJ molK1) corresponds to the bulk heat of fusion, which is close to the experimental heat of fusion of 25.5 kJ/ mol [24]. There might exist differences in local density between the nuclei and the remained parts. However, the difference in local density is expected to be small, contrary to a large difference in density between gas and liquid. On the other hand, in MD simulation, the structure of a cluster is distorted by thermal vibration, especial for the surface of a cluster. Therefore, the computed density of a cluster is noisy and rough. Fig. 2 shows the average number density of (KBr)864 as functions of axial distance (r) (the density is averaged over three axial directions) and time for the quenching run of (KBr)864 from 1060 to 550 K. It is not difficult to see that the densities of (KBr)864 clusters in core are between those of bulk liquid and solid KBr as expected. In order to examine the time evolution of the mean density of (KBr)864 cluster during quenching run, the densities in Fig. 2 are further averaged over the axial ˚ . The derived densities are plotted distance of 6.6–16.5 A against time in Fig 3(c). It demonstrates that the average density of cluster increases with time evolution of nucleation and crystal growth, but it is not sensitive enough for identifying nucleation time.

0.3 0.2 0.1 0.0 0

50

100 time(ps)

(c) 0.0130

D(A-3)

1734

0.0125

0.0120

0.0115 0

50

100 time(ps)

Fig. 3. Time dependence of the total energy (a), global order Parameter Q6 (b), and average number density (c) during the quenching run of (KBr)864 from 1060 to 550 K.

3.2. Nucleation time

D(A–3)

0.015

bulk solid(298K)

0.010 bulk liquid(945K) 0.005 0.000 5

10

15

20 r(A)

25

30

35

Fig. 2. The number density as functions of the dxial distance (r) and time (21.0–207.9 ps, with time interval of 21.0 ps) during the quenching run of (KBr)864 from 1060 to 550 K.

Clearly, crystallization information must exhibit the changes of local atomic arrangement and energy. Therefore, the nucleation time during quenching run can be identified by total energy (Fig. 3(a)) and bond order parameter (Fig. 3(b)). In bond order parameter method [22], the vectors rij connecting contiguous atoms i and j are defined as ‘bonds’ to distinguish the liquid atom and solid atom within a given radius rcut of i in cluster. The orientation of a bond rij for a given reference system is constrained by the spherical harmonics functions Ylm ðqij ; fij Þ h Ylm ðrij Þ, where qij and Fij are the polar and azimuthal angles of vectors rij in the reference system. The local structures around ion i can be

X. Zhu, K. Chen / Journal of Physics and Chemistry of Solids 66 (2005) 1732–1738

represented by qlm ðiÞ Z

Nb ðiÞ 1 X Y ðr Þ Nb ðiÞ jZ1 lm ij

and solid phases is ignored, the nucleation energy barrier DG* for a spherical nucleus is represented in terms of (2)

The second-order invariant is expressed in terms of " #1=2 l 4p X 2 ql ðiÞ h jq ðiÞ j : (3) 2l C 1 mZKl lm The global order parameter Ql can be derived from !1=2 l 4p X 2  Ql Z jQ j (4) 2l C 1 mZKl lm

DG * Z 16ps3sl =3DG2v

0.8

Q lm Z iZ1 N P

(5)

450K 500K 550K 600K

0 In(Nn/N0

Nb ðiÞqlm ðiÞ

(9)

where ssl is the interfacial free energy between solid and liquid, DGv is the free energy change of freezing per unit volume of solid. For slow cooling runs (Fig. 1(b)), the self-diffusion coefficients can be estimated from the slope of the mean-

where N P

1735

–0.8 –1.6 –2.4 –3.2

Nb ðiÞ

iZ1

–4

In MD simulation, the structure of cluster will be modified by thermal vibration. Fortunately, we can examine the distribution of the global order parameter Ql The time of nucleation in any quenching run can be estimated as time at which the global order parameter Q6 exhibit steadily increasing with time evolution as shown in Fig. 3(b).

(a) (KBr)

–4.8 0

20

40

60

80

1.6

3.3. Nucleation rate

400K 450K 500K 550K 600K

(6)

Where N0 represents the number of clusters, J is nucleation rate, Vc is approximately taken as the volume of the cluster, t0 is the time lag, and tn is the nucleation time (a time when the nth nucleation event has taken place). Nn Z N0 KnC 1. If the ln(Nn/N0) is plotted against the time tn, as exhibited in Fig. 4, the nucleation rates (J) can be inferred from the slopes of fitting lines. According to the theory of homogeneous nucleation, the nucleation rate can be expressed as J Z A expðKDG =kB TÞ

(7)

where A is a prefactor and DG* is the nucleation energy barrier. In nucleation analysis, classical nucleation theory (CNT) developed by Turnbull and co-worker is taken into account. Based on CNT, the prefactor is formulated by A Z 16ð3=4pÞ1=3 ðssl =kB TÞ1=2 D=Vm2=3 Dr 2

(8)

where D is the coefficient of diffusion in the liquid, Vm is the volume of a molecule in solid, and Dr is the molecular jump distance from liquid to the solid, which can be taken as Vm1=3 . On the other hand, if the strain energy between the liquid

In(Nn/N0

0 –0.8 –1.6 –2.4 –3.2

(b) (KBr) –4 0

50

100

150

500

200

250

time(ps) 2 1

400K 450K 500K 550K

0 In(Nn/N0

Nn Z N0 eKJVc ðtnKt0 Þ

100

time(ps)

0.8

It is well known that the number of unfrozen clusters, Nn, follows first-order rate law at a given temperature, i.e.

256

–1 –2 –3 –4

(c) (KBr)

864

–5 0

20

40

60

80

100

time(ps) Fig. 4. Plots of ln(Nn/N0) vs time of nucleation.

120

1736

X. Zhu, K. Chen / Journal of Physics and Chemistry of Solids 66 (2005) 1732–1738

Table 2 Nucleation rates (!10-36, mK3 s-1), growth rates (m sK1), interfacial free energies (mJ m-2), and sizes of critical nuclei T (K)

J

G

s

Nc

cal1a

cal2b

0.7

(KBr)256 450 500 550 600

7.41(40) 6.30(22) 4.21(22) 2.79(11)

29.5 30.8 30.3 28.5

42.4 43.5 44.4 44.4

41.1 42.8 44.4 46.1

9.5 11.7 15.0 19.0

400 450 500 550 600

2.30(17) 1.47(6) 1.39(6) 1.04(9) 0.37(3)

37.2 39.7 40.8 39.9 36.9

45.0 48.2 49.3 50.0 51.2

46.8 48.3 49.7 51.2 52.6

6.5 9.1 11.5 14.8 20.7

400 450 500 550

1.00(3) 0.82(4) 0.63(3) 0.56(5)

40.4 42.9 43.3 42.2

47.6 50.3 52.0 52.3

49.9 51.3 52.7 54.0

6.5 8.7 11.6 14.7

(KBr)500

(KBr)864

a b

From this work. From Turnbull relation.

square displacement by DZ

1 dhr 2 ðrÞi 6 dt

0.8

(10)

It is not surprised to note that the coefficients of BrK ion are slightly larger than those of KK ion for (KBr)256, (KBr)500, and (KBr)864. Arrheniius parameters for the averaged coefficients of KK and BrK diffusion in supercooled liquid potassium bromide cluster (for DZ D0 expðKEa =RTÞ) are given in Table 3, which are useful for nucleation analysis because the coefficient D is included in the prefactor of nucleation equation. Since 30 independent quenching runs are performed on KBr liquid clusters for each freezing temperature and each cluster size, the nucleation rates derived from the nucleation MD runs are relatively reliable. It is not difficult to see from Table 2 that the nucleation rates increase with the decreasing of temperature and cluster size. Note that for each cluster size, the free energy barriers to nucleation DG* linearly depend on the freezing temperatures as shown in Fig. 5. In this case, the factor exp(KDG*/kBT) in Eq. (7) becomes a constant in the temperature range studied. Thus, the nucleation rates will depend on prefactors (A). Since increasing temperatures prevail over the increased cost in an interfacial free energies and diffusion coefficients in the temperature ranges are studied, the nucleation rates slowly decrease with increasing temperature as exhibited in Table 2, which is consistent with the nucleation theory. On the other hand, the nucleation rates increase with decreasing cluster size as shown in Table 2. Several factors can account for such tendency. First, as shown in Table 3, the diffusion coefficients are greater in smaller clusters, which are included in the nucleation prefactor (A). Second, the nucleation starts near the surface of a cluster

∆G*

Cluster

0.9

0.6 0.5 0.4 0.3 0.2 350

400

450

500

550

600

650

T(K) Fig. 5. Temperature dependence of nucleation energy barrier. Open squares for (KBr)256; open circles for (KBr)500; open triangles for (KBr)864.

for alkali halide [4,6,7,9,10]. Smaller cluster has larger surface-to-volume ratio, which will result in more nucleation sites. Third, the size effect of nucleation rate is also related to the size effects of heat capacity and interfacial free energy. In addition, the Laplace pressure [9] can have effect on nucleation rate. 3.4. Growth rate Kelton and Greer [25] proposed the average sizedependent growth rate G using reaction rate theory, i.e. 16D GðRÞ Z 2 Dr







3Vm Vm 2ssl DGv K sinh 4p 2kB T R

(11)

where D is the diffusion coefficient in the supercooled liquid, Dr can be considered as Vm1=3 , Vm is the volume per molecule in solid phase, DGv is the free energy decrease per unit volume, and ssl is the interfacial free energy per unit area. Based on Eq. (11), the growth rates G are estimated for different temperature and cluster size as shown in Table 2. Table 3 Thermodynamic and physical properties used in the nucleation analysis Properties

Value

Reference

DHfus (kJ/mol), (KBr)256 DHfus (kJ/mol), (KBr)500 DHfus (kJ/mol), (KBr)864 Cp(l)–Cp(s) (J mol-1 KK1), (KBr)256 Cp(l)–Cp(s) (J mol-1KK1), (KBr)500 Cp(l)-Cp(s) (Jmol-1 KK1), (KBr)864 D (m2/s, !108), (KBr)256

19.4 21.0 21.7 10.75

This work This work This work This work

9.46

This work

8.95

This work

D(m2/s, !108), (KBr)500 D(m2/s, !108), (KBr)864 V (m3 mol-1)

2.18 exp(K9700/RT) 7.54 exp(K18, 270/RT) 11.48 exp(K21, 680/RT) 43.2!10K6

This work This work This work [23]

X. Zhu, K. Chen / Journal of Physics and Chemistry of Solids 66 (2005) 1732–1738

3.5. Mononuclear and polynuclear growth

50

3.6. Interfacial free energy It is well known that the value of solid-liquid interfacial free energy ssl is very difficult to measure directly. Fortunately, the interfacial free energy can be inferred from the MD nucleation rate in terms of classical relations. The results in Table 2 demonstrate that the interfacial free energies tend to slowly increase with increasing temperature and cluster size. An empirical estimate of ssl, suggested by Turnbull [27], can be predicted by kT DHfus ðNA V 2 Þ1=3

(12)

where NA is Avogadro’s number and V is the molar volume. kT is taken as 0.32. As shown in Table 2, our magnitude of the interfacial free energy, 42.4–52.3 mJ/m2, is reasonable when compared to Turnbull’s estimation of 41.1–54.0 mJ/m2 in the temperature range studied. The ssl of bulk KBr derived from experiment by Buckle and Ubbelohde [16] is 57.1 mJ/m2 at 800 K. This value is comparable to our result with consideration of the small temperature and cluster size dependences of ssl.

40

Nc

Kashchiev [26] has suggested a criterion for whether the crystallization in a given volume is mononuclear growth, or polynuclear growth, that is, if the ratio G=JVc4=3 is obviously greater than one, the freezing is mononuclear. It can be comprehended that if the growth of the first critical nucleus is so fast that there is not enough time for the second nucleus to grow, the crystallization is mononuclear. Based on the results in Table 2, we can estimate the ratio G=JVc4=3 for different temperature and cluster size. It is not difficult to see that all of the ratios are less than one. Thus, it implies that polynuclear growth may occur in the freezing of (KBr)256, (KBr)500, and (KBr)864. Such behavior is observed in the simulations. However, most quenching runs for smaller clusters lead to single crystals. In fact, during crystal growth process, sometimes, it is found that several independent nuclei grow up, match each other, and finally yield a single crystal.

ssl Z

1737

R=0.9829

30

20

10 400

500

600

700

800

900

1000 1100

Tf(K) Fig. 6. Correlation of critical nucleus size with freezing temperature.

to determine critical nucleus size. It is interesting to note that there is no cluster size effect on the sizes of critical nuclei as shown in Table 2. If sizes of critical nuclei for all cluster sizes are plotted against freezing temperatures, there is a good linear relation between them, as presented in Fig. 6. When all critical nucleus size for three cluster sizes are taken into account and extrapolated to the MD freezing point of bulk KBr (1012 K), we obtain the result 45 KCBrK ionic pairs (90 ions). Buckle and Ubbelohde [16,28,29] obtained liquid–solid interfacial free energy and critical sizes of some alkali halides from the time lag and the free energy of a bulk. The estimated critical nucleus size of NaCl is 99 ions at the temperature of 912 K, which is consistent with our result. 3.8. Crystallite structure Voronoi polyhedra method [21] is often used to analysis local arrangement of atoms around an atom in condensed system. Frenkel and co-workers [30] have studied the crystal nucleation in a Lennard–Jones system. This system has a stable face-centered-cubic (fcc) phase before melting. They found that the precritical nuclei are body-centered-cubic (bcc) ordered by using Voronoi polyhedra analysis. As they grow to their critical size, they become more fcc ordered in the core. In our Voronoi 450 400

3.7. Critical nucleus size

total

350

R*

ZK2ssl =DGv

(13)

300 NVP

In nucleation process, there is decay and growth balance for embryos of critical nuclei. Although, we can estimate the sizes of critical nuclei by identifying the number of ions satisfying fcc environments using bond order parameter method and Voronoi polyhedra technique, it is relatively rough because the structure of a cluster is distorted by thermal oscillation in MD simulation. Therefore, we prefer to use the classical expression

250 200

(0363)

150

(0444)

100 50 0

(0282) 0

50

100

150

200

250

time(ps) Fig. 7. Time evolution of the number of Voronoi polyhedra.

1738

X. Zhu, K. Chen / Journal of Physics and Chemistry of Solids 66 (2005) 1732–1738

polyhedra analysis, this behavior is not found. Fig. 7 illustrates that the precritical, critical, and postcritical nuclei have fcc structure. On the other hand, (0363)- and (0444)-polyhedra play an important role for the nucleation and crystal growth of (KBr)864.

4. Conclusions We carry out the MD simulations in (KBr)108, (KBr)256, (KBr)500, and (KBr)864 systems in order to estimate the critical nucleus size at the MD freezing temperature of bulk KBr. The thermodynamic properties including freezing temperature, heat capacity, heat of fusion, diffusion coefficient, and averaged number density of KBr cluster are inferred from the MD results. It is observed that the nucleation rates slowly decrease with increasing temperature, which is consistent with nucleation theory. On the other hand, the nucleation rates increase with decreasing cluster size. It is mainly due to the fact that the diffusion coefficients are greater in smaller clusters and alkali halides tend to nucleate near the surface of a cluster. The solid– liquid interfacial free energy of KBr clusters have the same magnitude as these of predicted by Turnbull’s relation and by Buckle and Ubbelohde. Critical nucleus size of potassium bromide derived from MD simulation together with classical nucleation theory is in good agreement with the experimental value of NaCl. The nucleation rates and critical nucleus sizes of KBr estimated from present MD simulations will be helpful for experimental nucleation studies of KBr clusters.

References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21]

[22]

[23] [24] [25] [26]

Acknowledgements This work is supported by a grants from Personnel Bureau of Nanjing of China (Project No. 2004103TSNB443) and Jiangsu Science & Technology Department of China (Project No. BK2005118).

[27] [28] [29] [30]

W.C. Swope, H.C. Andersen, Phys. Rev. B41 (1990) 7042. J.P. Rose, R.S. Berry, J. Chem. Phys. 98 (1993) 2346. J.P. Rose, R.S. Berry, J. Chem. Phys. 98 (1993) 3262. J. Huang, X. Zhu, L.S. Bartell, J. Phys. Chem. A102 (1998) 2708. L.S. Bartell, J. Huang, J. Phys. Chem. A102 (1998) 8722. X.L. Zhu, X.Z. You, et al., Chem. Phys. 250 (1999) 303. X.L. Zhu, X.Z. You, Chem. Phys. 269 (2001) 243. H. Deng, J. Huang, J. Solid State Chem. 159 (2001) 10. J. Huang, L.S. Bartell, J. Phys. Chem. A106 (2002) 2404. X.L. Zhu, J. Mol. Struct. (Theochem) 680 (2004) 137. K. Takahiro,Y, K. Takahiro, Y. Kengji, E. Tyoshikazu, J. Chem. Phys. 119 (2003) 11298. L.S. Bartell, T.S. Dibble, J. Am. Chem. Soc. 112 (1990) 890. L.S. Bartell, T.S. Dibble, J. Phys.Chem. 95 (1991) 1159. L.S. Bartell, T.S. Dibble, J. Phys.Chem. 96 (1992) 8203. L.S. Bartell, T.S. Dibble, J. Phys.Chem. 96 (1992) 2317. E.R. Buckle, A.R. Ubbelohde, Proc. R. Soc. London, Ser. A 261 (1961) 197. N. Anastasiou, D. Fincham, MDIONS, CCP5 Program Library, SERC DARESbury Laboratory, Daresbury, UK. Swanson, Tatge, JC Fel, Reports, NBS 1950. M.P. Tosi, F. Fumi, J. Phys. Chem. Solid 25 (1964) 45. M.L. Huggins, E. Mayer,J, J. Chem. Phys. 1 (1933) 3233. M. Tanemura, Y. Hiwatari, H. Matsuda, T. Ogawa, N. Ogita, A. Heda, Prog., Theor. Phys. 58 (1977) 1079; J.N. Cape, J.L. Finney, L.V. Woodcook, J. Chem Phys. 1 (1981) 2366. P.J. Steinhard, D.R. Nelson, M. Ronchett, Phys. Rev. B28 (1983) 784; Y. Chushak, L.S. Bartell, J. Phys Chem. A104 (2000) 9328. R.C. Weast, CRC Handbook of Chemistry and Physics, 63rd ed., CRC Press, Boca Raton, FL, 1982. 9 ed.. I. Barin, O. Knacke, Thermochemical Properties of Inorganic Substances, Springer-Verlag, Berlin, 1973. K.F. Kelton, A.L. Greer, J. Non-Cryst. Solids 79 (1986) 295. D. Kashchiev, D. Verdoes, G.M. van Rosmalen, J. Cryst. Growth 110 (1991) 373. D. Turnbull,, J. Appl. Phys. 21 (1950) 1022. E.R. Buckle, A.R. Ubbelohde, Proc. R. Soc. London, Ser. A 529 (1960) 325. E.R. Buckle, Proc. R. Soc. London, Ser. A 261 (1961) 189. Pieter Rein ten Wolde, Maria J. Ruiz-Montero, Daan Frenkel, J. Chem. Phys. 104 (1996) 9932.