- Email: [email protected]
Structural evolution in the crystallization of rapid cooling silver melt
Annals of Physics 354 (2015) 499–510
Contents lists available at ScienceDirect
Annals of Physics journal homepage: www.elsevier.com/locate/aop
Structural evolution in the crystallization of rapid cooling silver melt Z.A. Tian a,b,∗ , K.J. Dong b , A.B. Yu b a
School of Physics and Electronics, Hunan University, Changsha 410082, China
b
Laboratory for Simulation and Modelling of Particulate Systems School of Materials Science and Engineering, University of New South Wales, Sydney, NSW 2052, Australia
graphical
abstract
a
b
c
d
e
f
g
h
i
highlights • • • •
A comprehensive structural analysis is conducted focusing on crystallization. The involved atoms in our analysis are more than 90% for all samples concerned. A series of distinct intermediate states are found in crystallization of silver melt. A novelty icosahedron-saturated state breeds the metastable bcc state.
∗ Corresponding author at: Laboratory for Simulation and Modelling of Particulate Systems School of Materials Science and Engineering, University of New South Wales, Sydney, NSW 2052, Australia. E-mail address: [email protected] (Z.A. Tian). http://dx.doi.org/10.1016/j.aop.2014.12.021 0003-4916/© 2015 Elsevier Inc. All rights reserved.
500
article
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
info
Article history: Received 1 April 2014 Accepted 26 December 2014 Available online 9 January 2015 Keywords: Rapid cooling Molecular dynamics simulation Crystallization Structure evolution Ostwald’s rule of stages ISRO saturation stage
abstract The structural evolution in a rapid cooling process of silver melt has been investigated at different scales by adopting several analysis methods. The results testify Ostwald’s rule of stages and Frank conjecture upon icosahedron with many specific details. In particular, the cluster-scale analysis by a recent developed method called LSCA (the Largest Standard Cluster Analysis) clarified the complex structural evolution occurred in crystallization: different kinds of local clusters (such as ico-like (ico is the abbreviation of icosahedron), ico-bcc like (bcc, body-centred cubic), bcc, bcc-like structures) in turn have their maximal numbers as temperature decreases. And in a rather wide temperature range the icosahedral short-range order (ISRO) demonstrates a saturated stage (where the amount of ico-like structures keeps stable) that breeds metastable bcc clusters. As the precursor of crystallization, after reaching the maximal number bcc clusters finally decrease, resulting in the final solid being a mixture mainly composed of fcc/hcp (face-centred cubic and hexagonal-closed packed) clusters and to a less degree, bcc clusters. This detailed geometric picture for crystallization of liquid metal is believed to be useful to improve the fundamental understanding of liquid–solid phase transition. © 2015 Elsevier Inc. All rights reserved.
1. Introduction The properties of metal materials are significantly influenced by the kinetic pathway of the transformation during crystallization. In 1897 Ostwald [1] proposed the rule of stages that an unstable system could transform into another transient state by the smallest loss of free energy, before finally reaching a stable state. Although theoretical basis for this rule is yet to be established [2], there are increasing evidences supporting it in crystallization with slow kinetics, such as colloids [3–5], proteins [6,7], compositions [8], and some simulations based on ideal models, such as the Lennard-Jones fluids [9–13], or 2D square lattice [14]. It is very interesting to find out if this rule holds for crystallization of realistic metals, which have quick kinetics. An important related question is what intermediate states can be identified for a metal to transform from liquid to crystal? Stranski and Totomanow [15] argued that the phase first formed is the one that has the lowest free-energy barrier of formation rather than the phase that is globally stable under the conditions prevailing. Recent molecular dynamics (MD) studies for Fe, Mo, and V reveal that the solid–liquid interfacial free energy for body-centred cubic (bcc) crystal is lower than that for face-centred cubic (fcc) crystal by 30%–35% [16,17]. Thus bcc structures should probably be formed as an intermediate state from metal liquid to fcc crystal. Some studies have identified the metastable bcc phase at the earlier stage of solidification, but its role is not clearly identified. For example, a simulation for a 106 -particle L-J system reported that the bcc structures form prior to fcc/hcp ones in a cooling process, with its amount being so small that it is believed to play a very limited role in crystallization [18]. On the other hand, the recent MD simulations for silver melt demonstrate that an intermediate bcc phase is indeed the precursor of the final fcc crystal [19], but for such a small system (only containing 500 atoms) the scale effect is inevitable. Therefore, it is necessary to further investigate the role of the metastable bcc structures during crystallization of metals. In addition, to explain the phenomena that metal liquids can be supercooled to 0.2Tm (the melting temperature) without crystallization [20,21], Frank [22] argued that it is icosahedron (ico) that remarkably stabilizes the supercooled state. Since then, many experiments [23–25] and simulations [26–29] verified that the number of icosahedron is indeed increasing with the decrease of tempera-
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
501
ture in cooling processes. Then a question arises: what relationship exists between the icosahedron and the metastable bcc structures (if any)? A key issue responsible for the limited knowledge about the structural evolution at the earlier stage of solidification is due to the lack of effective analysis methods for different structures. For example, in characterizing the local structures of clusters that comprise a central atom and its neighbours, only a few (e.g. fcc, hcp, and bcc) kinds of clusters can be quantified, with their total amount usually less than 10% in most of the previous reports [18,30]. To date, it is not clear about other transitional structures or states formed by the rest (90%) atoms during such phase transitions. In this paper, an MD simulation of a rapid cooling process for silver liquid is conducted, with a large number (50,000) of atoms to eliminate the finite size effect. The so-called largest standard cluster analysis (LSCA) [31] is adopted to quantify the structural evolution in the liquid–solid transition. The method can cover all kinds of local clusters, and over 95% of atoms are shown to be involved in the 26 major clusters during the cooling process concerned. The results convincingly reveal that the cooling system experiences a series of distinct intermediate stages, marked by the (local) maxima in percentage at different temperatures. Particularly, the icosahedron-saturated stage breeds the metastable bcc clusters. It is concluded that the Frank conjecture is valid and the rule of stages is applied to this system. The rest of this paper is organized as follows. In Section 2, we describe the force field, simulation conduction, and the structural analysis methods. Section 3 is devoted to the identification of a hidden stage existing in crystalline phase transition. Section 4 shows the number evolution of major clusters. And the conclusions are presented in Section 5. 2. Computational methods 2.1. Force field The Quantum Sutton–Chen (QSC) many-body potential [32,33] is adopted here to simulate the solidification process of liquid silver. The total energy of a QSC system is thus calculated by Utot =
Ui =
1 i
i
= D·
1 i
2 j= ̸ i
2 i= ̸ j
1/2 i
D · V (rij ) − c · D · ρ
V (rij ) − c ρ
1/2 i
(1)
where rij is the distance between atoms i and j, V (rij ) is the pairwise repulsive potential, and ρi is the local energy density associated with the atom i. V (rij ) and ρi are defined as V (rij ) =
ρi =
n α
(2)
rij
ϕ(rij ) =
j̸=i
α m j̸=i
rij
.
(3)
The parameters of this potential are obtained by fitting experimental properties such as the density, cohesive energy, elastic constant, phonon-dispersion curves, together with quantum corrections through zero-point energy and pressure. In this way, the QSC potential can describe quantitatively the temperature-dependent properties in the melting and phase transformation [33]. Also, it produces accurate values of surface energies, vacancy formation energies, and stacking fault energies for crystal phases [34] and metallic glasses [35]. The values of the parameters for silver are c = 94.948, D = 4.0072 meV, m = 6, n = 11, and α = 4.0691 Å, the same as those used in the previous study [33]. 2.2. MD simulation The rapid cooling of silver liquid is simulated by using the constant-pressure MD techniques [36,37]. The MD simulations are performed in a cubic box with the three dimensional (3D) periodic
502
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
Fig. 1. Topology of the largest standard cluster (LSC) with a central atom. A tDh (truncated decahedron) LSC [2/555, 10/422] composed of a central atom and 12 neighbours (a); A CNS of S555 composed of a bonded reference pair and 5 CNNs (b), together with the topology between these 5 CNSs (c); Another CNS S422 (d), and the topology of its 4 CNNs (e).
boundary conditions. The potential is cut off at 22 atom unit (a.u.), and the simulated system contains 50,000 silver atoms. The equations of motion are integrated through the velocity leapfrog algorithm with a time step of 2.0 fs. The calculations start at 1973 K (the melting point Tm of silver is about 1233.85 K), and the initial configuration is produced based upon the Boltzmann distribution law. The system first runs enough time steps at 1973 K to obtain an equilibrium liquid, determined by the energy change of system. The damped force method (i.e., the Gaussian thermostat) [38,39] is then adopted to continuously decrease the system temperature to 273 K at a cooling rate of 5.0 × 1011 K/s. Every Kelvin, the velocities, positions and other data of every atom in the system are recorded for further analysis. 2.3. Analysis method In experiments, monitoring the time dependence of temperature (T –t curve) is a simple but effective method to detect phase transition. And crystallization can be identified by a level-off stage on the T –t curve, which results from the latent heat releasing. For a simulated cooling system at a certain cooling rate, an abrupt decrease on the E–T (Energy versus Temperature) curve corresponds to crystallization, while a slow and smooth change to a vitrification process. Four structural analysis methods are adopted in this work to investigate the structure at various length scales. The first method is based on the Pair Distribution Function (PDF or g (r )) which examines one-dimensional (1D) structure (of paired atoms). The validity of this simulation has been established by comparing the simulated and measured g (r ) in the previous study [19]. The second method is based on the Angular Distribution Function (ADF) related to three-body (two-dimensional (2D) structure) correlation [40,41], providing a type of overall information about the phase transition, such as the symmetry of microstructures. The third method is based on the so-called bond orientation order (BOO) Ql , which analyses the spherical harmonics and is sensitive to the orientation deviation of bonds (formed between a reference centre and its neighbours) [27]. The last method used for threedimensional (3D) structural analysis is LSCA [31], which can comprehensively quantify all kinds of local clusters beyond the nearest neighbours. The LSCA method is a more recent one and is briefly described here for clarity. As shown in Fig. 1, a cluster in LSCA is considered to be composed of a central atom and its near neighbours within a spherical region with a radius rc . In such a cluster, any pair of atoms are bonded or paired if they separate less than rc . A part of such cluster is called as a common-neighbour-subcluster (CNS, see Fig. 1(b)
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
503
and (d)), which comprises of a reference pair (the centre and one neighbour) and their common near neighbours (CNNs, Fig. 1(c) and (e)). According to the topology of the CNNs, a CNS can be systematically denoted by a CNS-index Sijk, where i is the number of CNNs, j the number of all bonds between CNNs, and k is the number of bonds in the longest continuous chain formed by part or all of the j bonds. The CNSs shown in Fig. 1(b) and (d) are S555 and S422 respectively. And in a 13-atom truncated decahedron (tDh, shown Fig. 1(a)), 10 neighbours comprise 10 S422 with the centre as well as their CNNs, and other 2 neighbours result in 2 S555, thus it is denoted as [2/555, 10/422]. Similarly, ico, fcc, hcp, and bcc clusters are [12/555], [12/421], [6/421, 6/422], and [6/444, 8/666], respectively. And for convenience, the centre of an ico LSC is hereafter referred to as an ico-atom, and so on. Therefore, all atoms in a system can be classified; and based on this basic classification, very large structures in a disordered system also can be readily identified. For example, the nano cluster shown in sFig. 1 (see Appendix) is a bcc crystalline cluster composed of 303 atoms, formed in the cooling process of liquid silver that will be discussed in detail in this paper. In the past, the cluster results depend on the value of rc [26,42], resulting in an uncertainty in structural quantification. This deficiency can be overcome by finding a natural and unique upper limit of rc , described as follows. For each possible reference pair, a too large rc will invalidate CNS-index. For example, as rc is large enough, j > i will be inevitable, which destroys the uniqueness of Sijk. Decreasing rc can remove the longest bond one by one till the CNS-index validity comes back. Thus for a reference pair, a unique maximal rc can be determined, and then such a maximal rc for all CNSs in a cluster can also be uniquely determined. This cluster is the largest one that all CNSs can be denoted by valid CNS-indexes, so it is called as the largest standard cluster (LSC), and the analysis is LSCA (the largest standard cluster analysis). Searching the unique rc for each cluster can be done using a computer with an appropriate program, such as that developed in our previous study [31]. Note that this parameter-free method is general and can be applied to other particulate systems, e.g., recently we have applied it in investigating the structural evolution of the packing of uniform spheres [43]. The four methods may have different advantages or disadvantages. However, their combined usage can lead to the construction of a comprehensive picture about the structural evolution. In the following sections, it will be shown that the overall feature can first be identified by the statistical property of energy, where the E–T curve shows that crystallization occurs at (630, 570) K. The structure evolution of the system can then be explored by other three methods. At the pairwise scale, the PDF curves confirm the outcome of the energy analysis. At the triple-wise scale, the ADF reveals the existence of a hidden stage during the crystallization, which is verified by the local BOO parameter Q4,local . At the cluster scale, LSCA provides more detailed information about local structures and their evolution. Based on the resulting information, the role of the icosahedral short-range order (ISRO) and metastable bcc clusters, as well as the relationship between energy change and structure evolution, are discussed. 3. A hidden stage in crystallization The much higher decreasing rate of E–T curve at T ∈ (625, 570) K (see Fig. 2(a)) reveals that this cooling process results in crystallization. The crystalline characteristic on the g (r ) curves can be recognized at 600 K and become more typical at 570 K (Fig. 2(b)) with the increasing of a sub peak between the first two major peaks. As T < 570 K, g (r ) curves have three extra distinct lower sub peaks within the first three major peaks. Briefly it can be determined that the crystallization commences and ends at about 625 K and 570 K respectively. However, neither specific structures (which should be indicated by discrete spectral PDF) in a sample nor structural evolution involved in successive samples can be identified by such continuous g (r ) curves. As a 2D analysis approach, ADF contour can produce more details about structural evolution during phase transition. Usually ADF is normalized as X (cos θ )dθ = 1, where X (cos θ ) is the probability at cos θ , and θ is the angle formed by two lines from the centre atom to any two neighbours in a LSC. However, typical structural characteristic may be smeared if many ADF curves are plotted in one picture, because the position and height of ADF peaks vary with temperature, caused by drastic structural changes during a phase transition. So the relative normalization is conducted before depicting ADF contours, where the maximum values of all ADF curves are aligned to an identical value of 1. That is, the symbol A(cos θ ) = X (cos θ )/ max{X (cos θ )} is employed to plot ADF contours.
504
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
Fig. 2. Temperature dependence of the potential (E–T curve) in the vicinity of solidification stage (a), and g (r ) curves at several selected temperatures (b).
Fig. 3. ADF contour for all LSCs in the system. For the lower peak in the range of 625–570 K, the distinct drop of its central position indicates that the system has transferred from liquid to a metastable stage; and its width also decreases with 590 K as a turning point. The double peaks centred at ±0.50 are the characteristic of fcc and/or hcp structures.
ADF contours in the temperature range of (770, 445) are depicted in Fig. 3. At a high temperature, the typical two-peak characteristic for liquid is very clear: the higher one is centred at 0.58 covering (0.51, 0.62), and the lower one covering (−0.51, −0.15) is smooth and wider with −0.33 as the centre (see the two dark strips near the right side of the contour box). With the decrease of temperature, the
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
505
Table 1 ADF data for three regular crystalline LSCs. cos θ
−1
fcc hcp bcc
0.09 0.05 0.08
√ − 3/3
−5/6
−1/2 0.36 0.27
0.09 0.26
a
−1/3
0
0.05 0.13
0.18 0.18 0.13
1/3
1/2
√
3/3
0.36 0.36 0.13
0.26
b
Fig. 4. The total ADF X (cos θ) for all LSCs at 590 K (a), and local Q4,local as a function of temperature during phase transition (b).
ADF patterns change evidently: the central position moves from −0.33 to −0.50 for the lower one; and from 0.58 to 0.50 for the higher one. As shown in the middle of Fig. 3, the nucleation can reasonably be regarded to be completed when the temperature decreases to 570 K, because after this temperature the ADF contour no longer distinctly changes in pattern. Comparing the typical characteristic of ADF contours as T < 570 K shown in Fig. 3 (with the distribution of cos θ for perfect crystals listed in Table 1, where vacations mean the ordinary zero probability, due to their discrete spectral ADFs), one can obtain that the final state of this cooling process should be mainly composed of fcc/hcp crystalline structures, because the major peaks at ±0.5 in the ADF contour are corresponding to those for fcc and √ hcp clusters rather than for bcc ones which hold two major speaks at ± 3/3 ≈ ±0.577. However, the ADF contour cannot determine whether the final solid is fcc or hcp crystal, or their mixture which may involve other kinds of clusters. A specific stage at around 590 K can be found in the ADF contour (see the middle of Fig. 3), where the width of the lower peak (centred at around −0.5) does not keep narrowing but broadening temporarily. Is there a hidden stage formed around 590 K? For further information, the ADF curve at 590 K is depicted in Fig. 4(a) that has three peaks centred at −0.5, 0.0, and 0.55 respectively, besides the one at −1. The centre of the highest peak of the total ADF is at about 0.55, which is larger than 0.5 corresponding to the √ major peak for fcc/hcp crystals but less than 0.577 for bcc. With the highest major peak-pair at ±( 3)/3, bcc clusters can move the centre of the highest peak of the total ADF from 0.5 to 0.55; and with the second peak-pair at ±1/3, bcc clusters can also effectively broaden the peak of the total ADF centred at −0.5. Note that other clusters in the supercooled liquid may also produce such √ effects, e.g., ico clusters with the highest major peak-pair at ±( 5)/5. In other words, ADF cannot clarify what is the hidden stage itself. To verify this hidden stage, we calculate the local BOO parameter [44] Q4,local , which is defined
as Q4,local =
N
j =1
4π 9
2 1/2 1 njb j , where nb is the number of nearest neighi=1 Y4m (θi , ϕi ) m=−4 nj
4
b
bours of sphere j, Y4m (θ , ϕ) is a spherical harmonic, and θi and ϕi are the polar and azimuthal angles of bond i. As shown in Fig. 4(b), with the decrease of temperature, Q4,local decreases to a minimum (about 0.1011) at around 590 K, and then increases. The existence of this minimum indicates a hidden stage, but without enough structural details. For example, it is 0.76, 0.191, 0.097, 0.036, and 0 for perfect regular local structures of sc (simple cubic), fcc, hcp, bcc, and ico LSCs, respectively. It is very
506
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
a
b
Fig. 5. The CNS total number (a) and the LSC type number (b) as a function of temperature during phase transition.
difficult to determine the distribution of components from such an overall statistic parameter that involves many kinds of possible local structures. Therefore a further analysis is necessary to reveal the structure of this hidden stage. 4. Structural evolution during phase transition LSCA can provide more detailed structural information. As the building blocks of LSC, totally 23 kinds of CNSs have been identified in 325 samples (one per degree from 773 K to 449 K), which is much less than that in sphere packing [43]. The total number of CNSs is probably associated with many factors (which needs further investigation in the future). However, the distinct maximum at around 590 K in Fig. 5(a) must correspond to a critical structural change at the CNS length scale. Generally, around each atom, a unique LSC can be identified by LSCA, thus there are always 50,000 LSCs in each sample of the simulated system. The type number varies remarkably with temperature. As shown in Fig. 5(b), the LSC type number decreases with the decrease of temperature. This suggests that the structures at the LSC scale are generally more ordered at lower temperatures. To reveal the predominant structural characteristic, one choice is to investigate the first N types (Top-N) of LSCs in number percentage for each sample. With N = 10 totally 26 (out of 131, 972) kinds of LSCs are picked out from the 325 samples considered, and for each sample more than 95% atoms are involved in these top-10 LSCs. Fig. 6 shows their %–T curves (the number percentage as a function of temperature) during the phase transition, among which Fig. 6(a) has once been reported in our previous paper as a preliminary result [31]. Although different from each other in structure, the 26 dominant LSCs consist of only 7 CNSs: S421, S422, S433, S444, S544, S555, and S666. Except for S433 and S544 that are believed to be the deformed S444 and S555 respectively, other 5 CNSs can comprise LSCs of regular structures, such as fcc [12/421], hcp [6/421,6/422], bcc [6/444, 8/666] crystals, as well as non-crystal icosahedron [12/555] and truncated decahedron [10/422, 2/555]. In terms of the number of CNSs in a LSC, all the 26 LSCs are here classified into 5 groups as follows, where for conciseness, the sum of S544 and S555 is represented by S5, and that of S433, S444 and S666 is S46. (1) (2) (3) (4) (5)
stable LSCs, including fcc and hcp ones that composed of S421 and S422 (Fig. 6(a)). bcc LSC, composed of 6 S444 and 8 S666 (Fig. 6(a)). ico-like LSCs, including ico and other 9 LSCs (Fig. 6(b) and (c)), in which S5 > S46. ico-bcc LSCs, including 5 LSCs shown in Fig. 6(d)–(f), where S5 ≈ S46. bcc-like LSCs, including 9 LSCs shown in Fig. 6(g)–(i), where S5 < S46.
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
a
b
c
d
e
f
g
h
i
507
Fig. 6. The number percentage of 26 LSCs as a function of temperature (%–T curve) in phase transition: (a) for 3 crystalline (fcc, hcp and bcc) LSCs; (b) and (c) for 9 ico-like LSCs (including icosahedron, S5 > S46); (d)–(f) for 5 ico-bcc LSCs (S5 ≈ S46); and (g)–(i) for 8 bcc-like LSCs (S5 < S46). S5 is the sum of S544 and S555, and S46 is the sum of S433, S444 and S666.
Interestingly, different LSCs in a same group are not only similar in structure, but also have a similar evolution during phase transition, as shown in Fig. 6. Generally, there are two evolutionary patterns: monotonic increase and of a single-peak. Only the stable (fcc and hcp) LSCs are in monotonic increase, indicating that the final solid is a mixture mainly composed of fcc and hcp structures and hence the phase transition is crystallization. This is consistent with the result obtained from energy, PDF and ADF analyses. Furthermore, the growth rate of fcc being always higher than hcp is also in line with that fcc is more stable than hcp by subtle energy difference [45,46]. The evolution of all other 24 LSCs is of a single peak; however, there are differences in details. First, the peaks are at different temperatures. For example, 4 ico-bcc (Fig. 6(d)–(f)), bcc (Fig. 6(a)), and one bcc-like ([4/444, 4/555, 6/666], Fig. 6(g)) LSCs have their maxima at about 610, 590, and 570 K, respectively. Second, some %–T curves present dogleg patterns around certain temperatures in the decrease stage, e.g., those for the 3 ico-bcc LSCs at 590 K (Fig. 6(d) and (f)), indicating the subtle changes in their decrease rates. Third, some peaks are very sharp, such as that for bcc LSCs (Fig. 6(a)) or 4 ico-bcc ones (Fig. 6(d)–(f)), while others are not so distinct, such as those for ico-like (Fig. 6(b) and (c)) and bcc-like (Fig. 6(h) and (i)) ones. To clarify the major structural characteristic, we depict the evolution for the three groups (ico-like, bcc-like and ico-bcc) that include multi LSCs, together with that for all the 26 LSCs. As for groups, the peaks on the plots in Fig. 7 are much sharper than their corresponding individuals shown in Fig. 6. Notably, for total ico-like LSCs the curve demonstrates a distinct level-off (saturation) stage covering (660, 630) K (see Fig. 7(a)). Fig. 7(b) and (c) show three critical temperatures of 610, 590, and 570 K with the abrupt changes on the plots. Fig. 7(d) not only highlights the critical temperature of 590 K, but also reveals that the present quantification is quite comprehensive, as the total amount of the 26 LSCs is always more than 25% (at 775 K), and close to 50% at the end of the phase transition (570 K). Correspondingly, atoms involved in the 26 LSCs are over 90% and 99.7%, respectively (see the inset in Fig. 7(d)). Figs. 6 and 7 give a clear picture about the phase transition: it is indeed achieved step by step, with ico-like, ico-bcc, bcc, bcc-like LSCs obtain their maxima in order. When temperature is high (i.e. T > 660 K), fcc and hcp clusters are very few (Fig. 6(a)), whereas other 24 (including bcc) LSCs exist in super-cooled liquid with rather large amounts and their amount increases slowly with the decrease of temperature. First, the ico-like LSCs have its maxima around 17% at about 660 K, followed by a level-off (saturation) stage ended at 630 K (see Fig. 7(a)), during which the increase of bcc and bcc-like LSCs accelerates abruptly at about 645 K (Figs. 6(a) and 7(c)). Second, the ico-bcc LSCs have their maximum about 13% at 610 K (Fig. 7(b)). The next critical temperature is 590 K where bcc LSCs have their maximum about 4.5%. Finally at 570 K, except fcc and hcp LSCs (see Fig. 6(a)), all others
508
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
a
b
c
d
Fig. 7. %–T curves for three groups of LSCs: ico-like (a), ico-bcc (b), bcc-like (c), and all 26 LSCs (d). The inset in (d) is for the atoms involved in all the 26 LSCs.
are decreasing marked with the beginning of the decrease in number of the bcc-like LSCs (Fig. 7(c)), and the ico LSCs decrease to very few and almost no longer decrease (Fig. 7(a)). Therefore the phase transition ends at 570 K. Only fcc/hcp LSCs keep increasing in the whole phase transition (Fig. 6(a)). The above structural evolution clarifies the specific stage around 590 K, which is aforementioned by the ADF contours (see Fig. 3), the BOO parameter Q4 (Fig. 4), and the CNS number (Fig. 5). First, it is indeed a special stage, evidenced by the facts: three groups (ico-bcc, bcc, and bcc-like) and all the 26 LSCs have critical changes around 590 K. Specifically, the bcc LSCs have their maxima (Fig. 6(a)); the decreasing rate of ico-bcc LSCs distinctly decreases (Fig. 7(b)); and the increasing rate of bcc-like LSCs and of all the 26 major LSCs reduces obviously (see Fig. 7(c) and (d)). Second, the major structures essentially have typical bcc characteristic. Besides the bcc LSCs, the structures of ico-bcc and bcc-like LSCs are also very close to bcc ones with a sizable fraction of S444 and S666. Thus this stage can be called the metastable bcc stage, of which the bcc-related LSCs increase in the first half while decrease in the second half. So far, we confirm that the phase transition proceeds via several steps with 660, 645, 630, 610, 590, and 570 K as the critical temperatures in terms of the evolution of the number percentage of LSCs. At this point, the relationship between energy change and structural evolution can be clearly phrased. The strong stabilization of ISRO for supercooled liquid decreases the energy difference between the LSCs that are inherent in liquids and their parent phase of supercooled liquid. Thus, as shown in Fig. 2(a), as bcc LSCs begin to remarkably increase at about 645 K, the slope of energy curve does not have obvious changes. So do other translational LSCs at 610 and 590 K. On the contrary, for the crystalline LSCs that are incompatible with disordered liquids, things are different. At 625 and 570 K, where the nucleation of fcc/hcp clusters begins and ends, the slope of energy curve changes remarkably, indicating that the free energy of fcc/hcp is much lower than bcc ones. This is consistent with the fact that the most stable phase of silver is not bcc but fcc crystal. In addition, the fcc and hcp almost simultaneously nucleate, indicating that the free energy difference between them is very small. It is also in good agreement with the recent simulation and experiment [45,46]. 5. Conclusions A rapid cooling process of a melt containing 50,000 silver atoms has been simulated by employing the quantum Sutton–Chen many-body potential. The energy analysis and g(r) plots reveal that the
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
509
liquid–solid transition is crystallization occurred between 625 and 570 K. The ADF contour and local BOO parameter Q4,local further demonstrate a transitional state during the crystallization, but the specific structures of the transitional stage cannot be clarified by the energy, g(r) or ADF analysis. A recently developed structural analysis method, i.e. LSCA, has been adopted to quantify the detailed structural evolution during this crystallization. Without any pre-set parameters, LSCA determines the LSC around each atom in a system by a simple topological criterion, giving a comprehensive and deterministic quantification of clusters including 10–20 atoms. Totally there are 131,972 clusters identified for the samples examined. The evolution of the major 26 local clusters involving over 95% atoms reveals that Ostwald’s rule of stages is valid in the solidification of silver liquid under rapid cooling. The phase transition proceeds step by step with ico-like, ico-bcc, bcc, bcc-like clusters obtaining their maximum numbers at T = (660, 630), 610, 590, and 570 K sequentially, and finally growing into a fcc/hcp mixed phase with less bcc structures. The most striking result is that the ISRO saturation stage remarkably declines the energy difference between bcc structures and its parent phase. Consistent with the Frank conjecture, the inherent ISRO in liquid metals becomes stronger with the increase of supercooling, and finally results in a saturation stage that covers (660, 630) K, where all the major 26 clusters are increasing and competing. And before it comes to an end the bcc clusters increase in a very high speed that is much bigger than that for all other clusters. These findings provide a more detailed and clearer geometric picture for crystallization of liquid metal. They can help improve the fundamental understanding of the related systems too. Acknowledgements The authors are grateful to the Australian Research Council (grant No: DE120100960) for the partial financial support of this work. Z. A. Tian is also grateful to the University of New South Wales for a Vice Chancellor’s postdoctoral fellowship (SIR50/MATSCI/PS22949). Appendix A. Supplementary data Supplementary material related to this article can be found online at http://dx.doi.org/10.1016/j. aop.2014.12.021. 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]
W. Ostwald, Z. Phys. Chem. 22 (1897) 289. P.R. ten Wolde, D. Frenkel, Science 277 (1997) 1975. D.F. Stefan Auer, Nature 409 (2001) 1020. V.J. Anderson, H.N.W. Lekkerkerker, Nature 416 (2002) 811. E. Sanz, C. Valeriani, D. Frenkel, M. Dijkstra, Phys. Rev. Lett. 99 (2007) 055501. P.G. Vekilov, Cryst. Growth Des. 4 (2004) 671. S.H. Shin, S. Chung, B. Sanii, L.R. Comollie, C.R. Bertozzi, J.J.D. Yoreoa, Proc. Natl. Acad. Sci. USA 109 (2012) 12968. S.Y. Chung, Y.M. Kim, J. Kim, Y.J. Kim, Nat. Phys. 5 (2009) 68. P.R. ten Wolde, M.J. Ruiz-Montero, D. Frenkel, Phys. Rev. Lett. 75 (1995) 2714. P.R. ten Wolde, M.J. Ruiz-Montero, D. Frenkel, J. Chem. Phys. 104 (1996) 9932. D. Moroni, P.R. ten Wolde, P.G. Bolhuis, Phys. Rev. Lett. 94 (2005) 235703. J.F. Lutsko, G. Nicolis, Phys. Rev. Lett. 96 (2006) 046102. C. Desgranges, J. Delhommelle, J. Am. Chem. Soc. 128 (2006) 10368. L.O. Hedges, S. Whitelam, J. Chem. Phys. 135 (2011) 164902. I.N. Stranski, Z. Totomanow, Z. Phys. Chem. 163 (1933) 339. D.Y. Sun, M. Asta, J.J. Hoyt, M.I. Mendelev, D.J. Srolovitz, Phys. Rev. B 69 (2004) 020102. J.J. Hoyt, M. Asta, D.Y. Sun, Phil. Mag. 86 (2006) 3651. W.C. Swope, H.C. Andersen, Phys. Rev. B 41 (1990) 7042. Z.A. Tian, R.S. Liu, C.X. Zheng, H.R. Liu, Z.Y. Hou, P. Peng, J. Phys. Chem. A 112 (2008) 12326. D. Turnbull, J. Chem. Phys. 20 (1952) 411. D. Turnbull, J. Appl. Phys. 21 (1950) 804. F.C. Frank, Proc. R. Soc. Lond. Ser. A 215 (1952) 43. M. Celino, V. Rosato, A. Di Cicco, A. Trapananti, C. Massobrio, Phys. Rev. B 75 (2007) 174210. T. Schenk, D. Holland-Moritz, V. Simonet, R. Bellissent, D.M. Herlach, Phys. Rev. Lett. 89 (2002) 075507.
510
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
[25] G.W. Lee, A.K. Gangopadhyay, K.F. Kelton, R.W. Hyers, T.J. Rathz, J.R. Rogers, D.S. Robinson, Phys. Rev. Lett. 93 (2004) 037802. [26] P. Ganesh, M. Widom, Phys. Rev. B 74 (2006) 134205. [27] P.J. Steinhardt, D.R. Nelson, M. Ronchetti, Phys. Rev. B 28 (1983) 784. [28] S. Nos, F. Yonezawa, J. Chem. Phys. 84 (1985) 1803. [29] H. Jónsson, H.C. Andersen, Phys. Rev. Lett. 60 (1988) 2295. [30] R.S. Liu, K.J. Dong, Z.A. Tian, H.R. Liu, P. Peng, A.B. Yu, J. Phys.: Condens. Matter 19 (2007) 196103. [31] Z.A. Tian, R.S. Liu, K.J. Dong, A.B. Yu, Europhys. Lett. 96 (2011) 36001. [32] A. Sutton, J. Chen, Phil. Mag. Lett. 61 (1990) 139. [33] Y. Qi, T. Cagin, Y. Kimura, W.A. Goddard, Phys. Rev. B 59 (1999) 3527. [34] P. Xu, T. Cagin, W.A. Goddard, J. Chem. Phys. 123 (2005) 104506. [35] H.J. Lee, T. Cagin, W.L. Johnson, W.A. Goddard, J. Chem. Phys. 119 (2003) 9858. [36] H.C. Andersen, J. Chem. Phys. 72 (1980) 2384. [37] D. Brown, J.H.R. Clarke, Mol. Phys. 51 (1984) 1243. [38] W.G. Hoover, A.J.C. Ladd, B. Moran, Phys. Rev. Lett. 48 (1982) 1818. [39] D.J. Evans, J. Chem. Phys. 78 (1983) 1818. [40] R.D. Mountain, P.K. Basu, J. Chem. Phys. 78 (1983) 7318. [41] D.H. Li, R.A. Moore, S. Wang, J. Chem. Phys. 89 (1988) 4309. [42] D.Q. Yu, M. Chen, X.J. Han, Phys. Rev. E 72 (2005) 051202. [43] Z.A. Tian, K.J. Dong, A.B. Yu, Phys. Rev. E 89 (2014) 032202. [44] A.R. Kansal, S. Torquato, F.H. Stillinger, Phys. Rev. E 66 (2002) 041109. [45] A.D. Bruce, A.N. Jackson, G.J. Ackland, N.B. Wilding, Phys. Rev. E 61 (2000) 906. [46] J.X. Zhu, M. Li, R. Rogers, W. Meyer, R.H. Ottewill, W.B. Russell, P.M. Chaikin, Nature 387 (1997) 883.
Contents lists available at ScienceDirect
Annals of Physics journal homepage: www.elsevier.com/locate/aop
Structural evolution in the crystallization of rapid cooling silver melt Z.A. Tian a,b,∗ , K.J. Dong b , A.B. Yu b a
School of Physics and Electronics, Hunan University, Changsha 410082, China
b
Laboratory for Simulation and Modelling of Particulate Systems School of Materials Science and Engineering, University of New South Wales, Sydney, NSW 2052, Australia
graphical
abstract
a
b
c
d
e
f
g
h
i
highlights • • • •
A comprehensive structural analysis is conducted focusing on crystallization. The involved atoms in our analysis are more than 90% for all samples concerned. A series of distinct intermediate states are found in crystallization of silver melt. A novelty icosahedron-saturated state breeds the metastable bcc state.
∗ Corresponding author at: Laboratory for Simulation and Modelling of Particulate Systems School of Materials Science and Engineering, University of New South Wales, Sydney, NSW 2052, Australia. E-mail address: [email protected] (Z.A. Tian). http://dx.doi.org/10.1016/j.aop.2014.12.021 0003-4916/© 2015 Elsevier Inc. All rights reserved.
500
article
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
info
Article history: Received 1 April 2014 Accepted 26 December 2014 Available online 9 January 2015 Keywords: Rapid cooling Molecular dynamics simulation Crystallization Structure evolution Ostwald’s rule of stages ISRO saturation stage
abstract The structural evolution in a rapid cooling process of silver melt has been investigated at different scales by adopting several analysis methods. The results testify Ostwald’s rule of stages and Frank conjecture upon icosahedron with many specific details. In particular, the cluster-scale analysis by a recent developed method called LSCA (the Largest Standard Cluster Analysis) clarified the complex structural evolution occurred in crystallization: different kinds of local clusters (such as ico-like (ico is the abbreviation of icosahedron), ico-bcc like (bcc, body-centred cubic), bcc, bcc-like structures) in turn have their maximal numbers as temperature decreases. And in a rather wide temperature range the icosahedral short-range order (ISRO) demonstrates a saturated stage (where the amount of ico-like structures keeps stable) that breeds metastable bcc clusters. As the precursor of crystallization, after reaching the maximal number bcc clusters finally decrease, resulting in the final solid being a mixture mainly composed of fcc/hcp (face-centred cubic and hexagonal-closed packed) clusters and to a less degree, bcc clusters. This detailed geometric picture for crystallization of liquid metal is believed to be useful to improve the fundamental understanding of liquid–solid phase transition. © 2015 Elsevier Inc. All rights reserved.
1. Introduction The properties of metal materials are significantly influenced by the kinetic pathway of the transformation during crystallization. In 1897 Ostwald [1] proposed the rule of stages that an unstable system could transform into another transient state by the smallest loss of free energy, before finally reaching a stable state. Although theoretical basis for this rule is yet to be established [2], there are increasing evidences supporting it in crystallization with slow kinetics, such as colloids [3–5], proteins [6,7], compositions [8], and some simulations based on ideal models, such as the Lennard-Jones fluids [9–13], or 2D square lattice [14]. It is very interesting to find out if this rule holds for crystallization of realistic metals, which have quick kinetics. An important related question is what intermediate states can be identified for a metal to transform from liquid to crystal? Stranski and Totomanow [15] argued that the phase first formed is the one that has the lowest free-energy barrier of formation rather than the phase that is globally stable under the conditions prevailing. Recent molecular dynamics (MD) studies for Fe, Mo, and V reveal that the solid–liquid interfacial free energy for body-centred cubic (bcc) crystal is lower than that for face-centred cubic (fcc) crystal by 30%–35% [16,17]. Thus bcc structures should probably be formed as an intermediate state from metal liquid to fcc crystal. Some studies have identified the metastable bcc phase at the earlier stage of solidification, but its role is not clearly identified. For example, a simulation for a 106 -particle L-J system reported that the bcc structures form prior to fcc/hcp ones in a cooling process, with its amount being so small that it is believed to play a very limited role in crystallization [18]. On the other hand, the recent MD simulations for silver melt demonstrate that an intermediate bcc phase is indeed the precursor of the final fcc crystal [19], but for such a small system (only containing 500 atoms) the scale effect is inevitable. Therefore, it is necessary to further investigate the role of the metastable bcc structures during crystallization of metals. In addition, to explain the phenomena that metal liquids can be supercooled to 0.2Tm (the melting temperature) without crystallization [20,21], Frank [22] argued that it is icosahedron (ico) that remarkably stabilizes the supercooled state. Since then, many experiments [23–25] and simulations [26–29] verified that the number of icosahedron is indeed increasing with the decrease of tempera-
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
501
ture in cooling processes. Then a question arises: what relationship exists between the icosahedron and the metastable bcc structures (if any)? A key issue responsible for the limited knowledge about the structural evolution at the earlier stage of solidification is due to the lack of effective analysis methods for different structures. For example, in characterizing the local structures of clusters that comprise a central atom and its neighbours, only a few (e.g. fcc, hcp, and bcc) kinds of clusters can be quantified, with their total amount usually less than 10% in most of the previous reports [18,30]. To date, it is not clear about other transitional structures or states formed by the rest (90%) atoms during such phase transitions. In this paper, an MD simulation of a rapid cooling process for silver liquid is conducted, with a large number (50,000) of atoms to eliminate the finite size effect. The so-called largest standard cluster analysis (LSCA) [31] is adopted to quantify the structural evolution in the liquid–solid transition. The method can cover all kinds of local clusters, and over 95% of atoms are shown to be involved in the 26 major clusters during the cooling process concerned. The results convincingly reveal that the cooling system experiences a series of distinct intermediate stages, marked by the (local) maxima in percentage at different temperatures. Particularly, the icosahedron-saturated stage breeds the metastable bcc clusters. It is concluded that the Frank conjecture is valid and the rule of stages is applied to this system. The rest of this paper is organized as follows. In Section 2, we describe the force field, simulation conduction, and the structural analysis methods. Section 3 is devoted to the identification of a hidden stage existing in crystalline phase transition. Section 4 shows the number evolution of major clusters. And the conclusions are presented in Section 5. 2. Computational methods 2.1. Force field The Quantum Sutton–Chen (QSC) many-body potential [32,33] is adopted here to simulate the solidification process of liquid silver. The total energy of a QSC system is thus calculated by Utot =
Ui =
1 i
i
= D·
1 i
2 j= ̸ i
2 i= ̸ j
1/2 i
D · V (rij ) − c · D · ρ
V (rij ) − c ρ
1/2 i
(1)
where rij is the distance between atoms i and j, V (rij ) is the pairwise repulsive potential, and ρi is the local energy density associated with the atom i. V (rij ) and ρi are defined as V (rij ) =
ρi =
n α
(2)
rij
ϕ(rij ) =
j̸=i
α m j̸=i
rij
.
(3)
The parameters of this potential are obtained by fitting experimental properties such as the density, cohesive energy, elastic constant, phonon-dispersion curves, together with quantum corrections through zero-point energy and pressure. In this way, the QSC potential can describe quantitatively the temperature-dependent properties in the melting and phase transformation [33]. Also, it produces accurate values of surface energies, vacancy formation energies, and stacking fault energies for crystal phases [34] and metallic glasses [35]. The values of the parameters for silver are c = 94.948, D = 4.0072 meV, m = 6, n = 11, and α = 4.0691 Å, the same as those used in the previous study [33]. 2.2. MD simulation The rapid cooling of silver liquid is simulated by using the constant-pressure MD techniques [36,37]. The MD simulations are performed in a cubic box with the three dimensional (3D) periodic
502
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
Fig. 1. Topology of the largest standard cluster (LSC) with a central atom. A tDh (truncated decahedron) LSC [2/555, 10/422] composed of a central atom and 12 neighbours (a); A CNS of S555 composed of a bonded reference pair and 5 CNNs (b), together with the topology between these 5 CNSs (c); Another CNS S422 (d), and the topology of its 4 CNNs (e).
boundary conditions. The potential is cut off at 22 atom unit (a.u.), and the simulated system contains 50,000 silver atoms. The equations of motion are integrated through the velocity leapfrog algorithm with a time step of 2.0 fs. The calculations start at 1973 K (the melting point Tm of silver is about 1233.85 K), and the initial configuration is produced based upon the Boltzmann distribution law. The system first runs enough time steps at 1973 K to obtain an equilibrium liquid, determined by the energy change of system. The damped force method (i.e., the Gaussian thermostat) [38,39] is then adopted to continuously decrease the system temperature to 273 K at a cooling rate of 5.0 × 1011 K/s. Every Kelvin, the velocities, positions and other data of every atom in the system are recorded for further analysis. 2.3. Analysis method In experiments, monitoring the time dependence of temperature (T –t curve) is a simple but effective method to detect phase transition. And crystallization can be identified by a level-off stage on the T –t curve, which results from the latent heat releasing. For a simulated cooling system at a certain cooling rate, an abrupt decrease on the E–T (Energy versus Temperature) curve corresponds to crystallization, while a slow and smooth change to a vitrification process. Four structural analysis methods are adopted in this work to investigate the structure at various length scales. The first method is based on the Pair Distribution Function (PDF or g (r )) which examines one-dimensional (1D) structure (of paired atoms). The validity of this simulation has been established by comparing the simulated and measured g (r ) in the previous study [19]. The second method is based on the Angular Distribution Function (ADF) related to three-body (two-dimensional (2D) structure) correlation [40,41], providing a type of overall information about the phase transition, such as the symmetry of microstructures. The third method is based on the so-called bond orientation order (BOO) Ql , which analyses the spherical harmonics and is sensitive to the orientation deviation of bonds (formed between a reference centre and its neighbours) [27]. The last method used for threedimensional (3D) structural analysis is LSCA [31], which can comprehensively quantify all kinds of local clusters beyond the nearest neighbours. The LSCA method is a more recent one and is briefly described here for clarity. As shown in Fig. 1, a cluster in LSCA is considered to be composed of a central atom and its near neighbours within a spherical region with a radius rc . In such a cluster, any pair of atoms are bonded or paired if they separate less than rc . A part of such cluster is called as a common-neighbour-subcluster (CNS, see Fig. 1(b)
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
503
and (d)), which comprises of a reference pair (the centre and one neighbour) and their common near neighbours (CNNs, Fig. 1(c) and (e)). According to the topology of the CNNs, a CNS can be systematically denoted by a CNS-index Sijk, where i is the number of CNNs, j the number of all bonds between CNNs, and k is the number of bonds in the longest continuous chain formed by part or all of the j bonds. The CNSs shown in Fig. 1(b) and (d) are S555 and S422 respectively. And in a 13-atom truncated decahedron (tDh, shown Fig. 1(a)), 10 neighbours comprise 10 S422 with the centre as well as their CNNs, and other 2 neighbours result in 2 S555, thus it is denoted as [2/555, 10/422]. Similarly, ico, fcc, hcp, and bcc clusters are [12/555], [12/421], [6/421, 6/422], and [6/444, 8/666], respectively. And for convenience, the centre of an ico LSC is hereafter referred to as an ico-atom, and so on. Therefore, all atoms in a system can be classified; and based on this basic classification, very large structures in a disordered system also can be readily identified. For example, the nano cluster shown in sFig. 1 (see Appendix) is a bcc crystalline cluster composed of 303 atoms, formed in the cooling process of liquid silver that will be discussed in detail in this paper. In the past, the cluster results depend on the value of rc [26,42], resulting in an uncertainty in structural quantification. This deficiency can be overcome by finding a natural and unique upper limit of rc , described as follows. For each possible reference pair, a too large rc will invalidate CNS-index. For example, as rc is large enough, j > i will be inevitable, which destroys the uniqueness of Sijk. Decreasing rc can remove the longest bond one by one till the CNS-index validity comes back. Thus for a reference pair, a unique maximal rc can be determined, and then such a maximal rc for all CNSs in a cluster can also be uniquely determined. This cluster is the largest one that all CNSs can be denoted by valid CNS-indexes, so it is called as the largest standard cluster (LSC), and the analysis is LSCA (the largest standard cluster analysis). Searching the unique rc for each cluster can be done using a computer with an appropriate program, such as that developed in our previous study [31]. Note that this parameter-free method is general and can be applied to other particulate systems, e.g., recently we have applied it in investigating the structural evolution of the packing of uniform spheres [43]. The four methods may have different advantages or disadvantages. However, their combined usage can lead to the construction of a comprehensive picture about the structural evolution. In the following sections, it will be shown that the overall feature can first be identified by the statistical property of energy, where the E–T curve shows that crystallization occurs at (630, 570) K. The structure evolution of the system can then be explored by other three methods. At the pairwise scale, the PDF curves confirm the outcome of the energy analysis. At the triple-wise scale, the ADF reveals the existence of a hidden stage during the crystallization, which is verified by the local BOO parameter Q4,local . At the cluster scale, LSCA provides more detailed information about local structures and their evolution. Based on the resulting information, the role of the icosahedral short-range order (ISRO) and metastable bcc clusters, as well as the relationship between energy change and structure evolution, are discussed. 3. A hidden stage in crystallization The much higher decreasing rate of E–T curve at T ∈ (625, 570) K (see Fig. 2(a)) reveals that this cooling process results in crystallization. The crystalline characteristic on the g (r ) curves can be recognized at 600 K and become more typical at 570 K (Fig. 2(b)) with the increasing of a sub peak between the first two major peaks. As T < 570 K, g (r ) curves have three extra distinct lower sub peaks within the first three major peaks. Briefly it can be determined that the crystallization commences and ends at about 625 K and 570 K respectively. However, neither specific structures (which should be indicated by discrete spectral PDF) in a sample nor structural evolution involved in successive samples can be identified by such continuous g (r ) curves. As a 2D analysis approach, ADF contour can produce more details about structural evolution during phase transition. Usually ADF is normalized as X (cos θ )dθ = 1, where X (cos θ ) is the probability at cos θ , and θ is the angle formed by two lines from the centre atom to any two neighbours in a LSC. However, typical structural characteristic may be smeared if many ADF curves are plotted in one picture, because the position and height of ADF peaks vary with temperature, caused by drastic structural changes during a phase transition. So the relative normalization is conducted before depicting ADF contours, where the maximum values of all ADF curves are aligned to an identical value of 1. That is, the symbol A(cos θ ) = X (cos θ )/ max{X (cos θ )} is employed to plot ADF contours.
504
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
Fig. 2. Temperature dependence of the potential (E–T curve) in the vicinity of solidification stage (a), and g (r ) curves at several selected temperatures (b).
Fig. 3. ADF contour for all LSCs in the system. For the lower peak in the range of 625–570 K, the distinct drop of its central position indicates that the system has transferred from liquid to a metastable stage; and its width also decreases with 590 K as a turning point. The double peaks centred at ±0.50 are the characteristic of fcc and/or hcp structures.
ADF contours in the temperature range of (770, 445) are depicted in Fig. 3. At a high temperature, the typical two-peak characteristic for liquid is very clear: the higher one is centred at 0.58 covering (0.51, 0.62), and the lower one covering (−0.51, −0.15) is smooth and wider with −0.33 as the centre (see the two dark strips near the right side of the contour box). With the decrease of temperature, the
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
505
Table 1 ADF data for three regular crystalline LSCs. cos θ
−1
fcc hcp bcc
0.09 0.05 0.08
√ − 3/3
−5/6
−1/2 0.36 0.27
0.09 0.26
a
−1/3
0
0.05 0.13
0.18 0.18 0.13
1/3
1/2
√
3/3
0.36 0.36 0.13
0.26
b
Fig. 4. The total ADF X (cos θ) for all LSCs at 590 K (a), and local Q4,local as a function of temperature during phase transition (b).
ADF patterns change evidently: the central position moves from −0.33 to −0.50 for the lower one; and from 0.58 to 0.50 for the higher one. As shown in the middle of Fig. 3, the nucleation can reasonably be regarded to be completed when the temperature decreases to 570 K, because after this temperature the ADF contour no longer distinctly changes in pattern. Comparing the typical characteristic of ADF contours as T < 570 K shown in Fig. 3 (with the distribution of cos θ for perfect crystals listed in Table 1, where vacations mean the ordinary zero probability, due to their discrete spectral ADFs), one can obtain that the final state of this cooling process should be mainly composed of fcc/hcp crystalline structures, because the major peaks at ±0.5 in the ADF contour are corresponding to those for fcc and √ hcp clusters rather than for bcc ones which hold two major speaks at ± 3/3 ≈ ±0.577. However, the ADF contour cannot determine whether the final solid is fcc or hcp crystal, or their mixture which may involve other kinds of clusters. A specific stage at around 590 K can be found in the ADF contour (see the middle of Fig. 3), where the width of the lower peak (centred at around −0.5) does not keep narrowing but broadening temporarily. Is there a hidden stage formed around 590 K? For further information, the ADF curve at 590 K is depicted in Fig. 4(a) that has three peaks centred at −0.5, 0.0, and 0.55 respectively, besides the one at −1. The centre of the highest peak of the total ADF is at about 0.55, which is larger than 0.5 corresponding to the √ major peak for fcc/hcp crystals but less than 0.577 for bcc. With the highest major peak-pair at ±( 3)/3, bcc clusters can move the centre of the highest peak of the total ADF from 0.5 to 0.55; and with the second peak-pair at ±1/3, bcc clusters can also effectively broaden the peak of the total ADF centred at −0.5. Note that other clusters in the supercooled liquid may also produce such √ effects, e.g., ico clusters with the highest major peak-pair at ±( 5)/5. In other words, ADF cannot clarify what is the hidden stage itself. To verify this hidden stage, we calculate the local BOO parameter [44] Q4,local , which is defined
as Q4,local =
N
j =1
4π 9
2 1/2 1 njb j , where nb is the number of nearest neighi=1 Y4m (θi , ϕi ) m=−4 nj
4
b
bours of sphere j, Y4m (θ , ϕ) is a spherical harmonic, and θi and ϕi are the polar and azimuthal angles of bond i. As shown in Fig. 4(b), with the decrease of temperature, Q4,local decreases to a minimum (about 0.1011) at around 590 K, and then increases. The existence of this minimum indicates a hidden stage, but without enough structural details. For example, it is 0.76, 0.191, 0.097, 0.036, and 0 for perfect regular local structures of sc (simple cubic), fcc, hcp, bcc, and ico LSCs, respectively. It is very
506
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
a
b
Fig. 5. The CNS total number (a) and the LSC type number (b) as a function of temperature during phase transition.
difficult to determine the distribution of components from such an overall statistic parameter that involves many kinds of possible local structures. Therefore a further analysis is necessary to reveal the structure of this hidden stage. 4. Structural evolution during phase transition LSCA can provide more detailed structural information. As the building blocks of LSC, totally 23 kinds of CNSs have been identified in 325 samples (one per degree from 773 K to 449 K), which is much less than that in sphere packing [43]. The total number of CNSs is probably associated with many factors (which needs further investigation in the future). However, the distinct maximum at around 590 K in Fig. 5(a) must correspond to a critical structural change at the CNS length scale. Generally, around each atom, a unique LSC can be identified by LSCA, thus there are always 50,000 LSCs in each sample of the simulated system. The type number varies remarkably with temperature. As shown in Fig. 5(b), the LSC type number decreases with the decrease of temperature. This suggests that the structures at the LSC scale are generally more ordered at lower temperatures. To reveal the predominant structural characteristic, one choice is to investigate the first N types (Top-N) of LSCs in number percentage for each sample. With N = 10 totally 26 (out of 131, 972) kinds of LSCs are picked out from the 325 samples considered, and for each sample more than 95% atoms are involved in these top-10 LSCs. Fig. 6 shows their %–T curves (the number percentage as a function of temperature) during the phase transition, among which Fig. 6(a) has once been reported in our previous paper as a preliminary result [31]. Although different from each other in structure, the 26 dominant LSCs consist of only 7 CNSs: S421, S422, S433, S444, S544, S555, and S666. Except for S433 and S544 that are believed to be the deformed S444 and S555 respectively, other 5 CNSs can comprise LSCs of regular structures, such as fcc [12/421], hcp [6/421,6/422], bcc [6/444, 8/666] crystals, as well as non-crystal icosahedron [12/555] and truncated decahedron [10/422, 2/555]. In terms of the number of CNSs in a LSC, all the 26 LSCs are here classified into 5 groups as follows, where for conciseness, the sum of S544 and S555 is represented by S5, and that of S433, S444 and S666 is S46. (1) (2) (3) (4) (5)
stable LSCs, including fcc and hcp ones that composed of S421 and S422 (Fig. 6(a)). bcc LSC, composed of 6 S444 and 8 S666 (Fig. 6(a)). ico-like LSCs, including ico and other 9 LSCs (Fig. 6(b) and (c)), in which S5 > S46. ico-bcc LSCs, including 5 LSCs shown in Fig. 6(d)–(f), where S5 ≈ S46. bcc-like LSCs, including 9 LSCs shown in Fig. 6(g)–(i), where S5 < S46.
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
a
b
c
d
e
f
g
h
i
507
Fig. 6. The number percentage of 26 LSCs as a function of temperature (%–T curve) in phase transition: (a) for 3 crystalline (fcc, hcp and bcc) LSCs; (b) and (c) for 9 ico-like LSCs (including icosahedron, S5 > S46); (d)–(f) for 5 ico-bcc LSCs (S5 ≈ S46); and (g)–(i) for 8 bcc-like LSCs (S5 < S46). S5 is the sum of S544 and S555, and S46 is the sum of S433, S444 and S666.
Interestingly, different LSCs in a same group are not only similar in structure, but also have a similar evolution during phase transition, as shown in Fig. 6. Generally, there are two evolutionary patterns: monotonic increase and of a single-peak. Only the stable (fcc and hcp) LSCs are in monotonic increase, indicating that the final solid is a mixture mainly composed of fcc and hcp structures and hence the phase transition is crystallization. This is consistent with the result obtained from energy, PDF and ADF analyses. Furthermore, the growth rate of fcc being always higher than hcp is also in line with that fcc is more stable than hcp by subtle energy difference [45,46]. The evolution of all other 24 LSCs is of a single peak; however, there are differences in details. First, the peaks are at different temperatures. For example, 4 ico-bcc (Fig. 6(d)–(f)), bcc (Fig. 6(a)), and one bcc-like ([4/444, 4/555, 6/666], Fig. 6(g)) LSCs have their maxima at about 610, 590, and 570 K, respectively. Second, some %–T curves present dogleg patterns around certain temperatures in the decrease stage, e.g., those for the 3 ico-bcc LSCs at 590 K (Fig. 6(d) and (f)), indicating the subtle changes in their decrease rates. Third, some peaks are very sharp, such as that for bcc LSCs (Fig. 6(a)) or 4 ico-bcc ones (Fig. 6(d)–(f)), while others are not so distinct, such as those for ico-like (Fig. 6(b) and (c)) and bcc-like (Fig. 6(h) and (i)) ones. To clarify the major structural characteristic, we depict the evolution for the three groups (ico-like, bcc-like and ico-bcc) that include multi LSCs, together with that for all the 26 LSCs. As for groups, the peaks on the plots in Fig. 7 are much sharper than their corresponding individuals shown in Fig. 6. Notably, for total ico-like LSCs the curve demonstrates a distinct level-off (saturation) stage covering (660, 630) K (see Fig. 7(a)). Fig. 7(b) and (c) show three critical temperatures of 610, 590, and 570 K with the abrupt changes on the plots. Fig. 7(d) not only highlights the critical temperature of 590 K, but also reveals that the present quantification is quite comprehensive, as the total amount of the 26 LSCs is always more than 25% (at 775 K), and close to 50% at the end of the phase transition (570 K). Correspondingly, atoms involved in the 26 LSCs are over 90% and 99.7%, respectively (see the inset in Fig. 7(d)). Figs. 6 and 7 give a clear picture about the phase transition: it is indeed achieved step by step, with ico-like, ico-bcc, bcc, bcc-like LSCs obtain their maxima in order. When temperature is high (i.e. T > 660 K), fcc and hcp clusters are very few (Fig. 6(a)), whereas other 24 (including bcc) LSCs exist in super-cooled liquid with rather large amounts and their amount increases slowly with the decrease of temperature. First, the ico-like LSCs have its maxima around 17% at about 660 K, followed by a level-off (saturation) stage ended at 630 K (see Fig. 7(a)), during which the increase of bcc and bcc-like LSCs accelerates abruptly at about 645 K (Figs. 6(a) and 7(c)). Second, the ico-bcc LSCs have their maximum about 13% at 610 K (Fig. 7(b)). The next critical temperature is 590 K where bcc LSCs have their maximum about 4.5%. Finally at 570 K, except fcc and hcp LSCs (see Fig. 6(a)), all others
508
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
a
b
c
d
Fig. 7. %–T curves for three groups of LSCs: ico-like (a), ico-bcc (b), bcc-like (c), and all 26 LSCs (d). The inset in (d) is for the atoms involved in all the 26 LSCs.
are decreasing marked with the beginning of the decrease in number of the bcc-like LSCs (Fig. 7(c)), and the ico LSCs decrease to very few and almost no longer decrease (Fig. 7(a)). Therefore the phase transition ends at 570 K. Only fcc/hcp LSCs keep increasing in the whole phase transition (Fig. 6(a)). The above structural evolution clarifies the specific stage around 590 K, which is aforementioned by the ADF contours (see Fig. 3), the BOO parameter Q4 (Fig. 4), and the CNS number (Fig. 5). First, it is indeed a special stage, evidenced by the facts: three groups (ico-bcc, bcc, and bcc-like) and all the 26 LSCs have critical changes around 590 K. Specifically, the bcc LSCs have their maxima (Fig. 6(a)); the decreasing rate of ico-bcc LSCs distinctly decreases (Fig. 7(b)); and the increasing rate of bcc-like LSCs and of all the 26 major LSCs reduces obviously (see Fig. 7(c) and (d)). Second, the major structures essentially have typical bcc characteristic. Besides the bcc LSCs, the structures of ico-bcc and bcc-like LSCs are also very close to bcc ones with a sizable fraction of S444 and S666. Thus this stage can be called the metastable bcc stage, of which the bcc-related LSCs increase in the first half while decrease in the second half. So far, we confirm that the phase transition proceeds via several steps with 660, 645, 630, 610, 590, and 570 K as the critical temperatures in terms of the evolution of the number percentage of LSCs. At this point, the relationship between energy change and structural evolution can be clearly phrased. The strong stabilization of ISRO for supercooled liquid decreases the energy difference between the LSCs that are inherent in liquids and their parent phase of supercooled liquid. Thus, as shown in Fig. 2(a), as bcc LSCs begin to remarkably increase at about 645 K, the slope of energy curve does not have obvious changes. So do other translational LSCs at 610 and 590 K. On the contrary, for the crystalline LSCs that are incompatible with disordered liquids, things are different. At 625 and 570 K, where the nucleation of fcc/hcp clusters begins and ends, the slope of energy curve changes remarkably, indicating that the free energy of fcc/hcp is much lower than bcc ones. This is consistent with the fact that the most stable phase of silver is not bcc but fcc crystal. In addition, the fcc and hcp almost simultaneously nucleate, indicating that the free energy difference between them is very small. It is also in good agreement with the recent simulation and experiment [45,46]. 5. Conclusions A rapid cooling process of a melt containing 50,000 silver atoms has been simulated by employing the quantum Sutton–Chen many-body potential. The energy analysis and g(r) plots reveal that the
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
509
liquid–solid transition is crystallization occurred between 625 and 570 K. The ADF contour and local BOO parameter Q4,local further demonstrate a transitional state during the crystallization, but the specific structures of the transitional stage cannot be clarified by the energy, g(r) or ADF analysis. A recently developed structural analysis method, i.e. LSCA, has been adopted to quantify the detailed structural evolution during this crystallization. Without any pre-set parameters, LSCA determines the LSC around each atom in a system by a simple topological criterion, giving a comprehensive and deterministic quantification of clusters including 10–20 atoms. Totally there are 131,972 clusters identified for the samples examined. The evolution of the major 26 local clusters involving over 95% atoms reveals that Ostwald’s rule of stages is valid in the solidification of silver liquid under rapid cooling. The phase transition proceeds step by step with ico-like, ico-bcc, bcc, bcc-like clusters obtaining their maximum numbers at T = (660, 630), 610, 590, and 570 K sequentially, and finally growing into a fcc/hcp mixed phase with less bcc structures. The most striking result is that the ISRO saturation stage remarkably declines the energy difference between bcc structures and its parent phase. Consistent with the Frank conjecture, the inherent ISRO in liquid metals becomes stronger with the increase of supercooling, and finally results in a saturation stage that covers (660, 630) K, where all the major 26 clusters are increasing and competing. And before it comes to an end the bcc clusters increase in a very high speed that is much bigger than that for all other clusters. These findings provide a more detailed and clearer geometric picture for crystallization of liquid metal. They can help improve the fundamental understanding of the related systems too. Acknowledgements The authors are grateful to the Australian Research Council (grant No: DE120100960) for the partial financial support of this work. Z. A. Tian is also grateful to the University of New South Wales for a Vice Chancellor’s postdoctoral fellowship (SIR50/MATSCI/PS22949). Appendix A. Supplementary data Supplementary material related to this article can be found online at http://dx.doi.org/10.1016/j. aop.2014.12.021. 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]
W. Ostwald, Z. Phys. Chem. 22 (1897) 289. P.R. ten Wolde, D. Frenkel, Science 277 (1997) 1975. D.F. Stefan Auer, Nature 409 (2001) 1020. V.J. Anderson, H.N.W. Lekkerkerker, Nature 416 (2002) 811. E. Sanz, C. Valeriani, D. Frenkel, M. Dijkstra, Phys. Rev. Lett. 99 (2007) 055501. P.G. Vekilov, Cryst. Growth Des. 4 (2004) 671. S.H. Shin, S. Chung, B. Sanii, L.R. Comollie, C.R. Bertozzi, J.J.D. Yoreoa, Proc. Natl. Acad. Sci. USA 109 (2012) 12968. S.Y. Chung, Y.M. Kim, J. Kim, Y.J. Kim, Nat. Phys. 5 (2009) 68. P.R. ten Wolde, M.J. Ruiz-Montero, D. Frenkel, Phys. Rev. Lett. 75 (1995) 2714. P.R. ten Wolde, M.J. Ruiz-Montero, D. Frenkel, J. Chem. Phys. 104 (1996) 9932. D. Moroni, P.R. ten Wolde, P.G. Bolhuis, Phys. Rev. Lett. 94 (2005) 235703. J.F. Lutsko, G. Nicolis, Phys. Rev. Lett. 96 (2006) 046102. C. Desgranges, J. Delhommelle, J. Am. Chem. Soc. 128 (2006) 10368. L.O. Hedges, S. Whitelam, J. Chem. Phys. 135 (2011) 164902. I.N. Stranski, Z. Totomanow, Z. Phys. Chem. 163 (1933) 339. D.Y. Sun, M. Asta, J.J. Hoyt, M.I. Mendelev, D.J. Srolovitz, Phys. Rev. B 69 (2004) 020102. J.J. Hoyt, M. Asta, D.Y. Sun, Phil. Mag. 86 (2006) 3651. W.C. Swope, H.C. Andersen, Phys. Rev. B 41 (1990) 7042. Z.A. Tian, R.S. Liu, C.X. Zheng, H.R. Liu, Z.Y. Hou, P. Peng, J. Phys. Chem. A 112 (2008) 12326. D. Turnbull, J. Chem. Phys. 20 (1952) 411. D. Turnbull, J. Appl. Phys. 21 (1950) 804. F.C. Frank, Proc. R. Soc. Lond. Ser. A 215 (1952) 43. M. Celino, V. Rosato, A. Di Cicco, A. Trapananti, C. Massobrio, Phys. Rev. B 75 (2007) 174210. T. Schenk, D. Holland-Moritz, V. Simonet, R. Bellissent, D.M. Herlach, Phys. Rev. Lett. 89 (2002) 075507.
510
Z.A. Tian et al. / Annals of Physics 354 (2015) 499–510
[25] G.W. Lee, A.K. Gangopadhyay, K.F. Kelton, R.W. Hyers, T.J. Rathz, J.R. Rogers, D.S. Robinson, Phys. Rev. Lett. 93 (2004) 037802. [26] P. Ganesh, M. Widom, Phys. Rev. B 74 (2006) 134205. [27] P.J. Steinhardt, D.R. Nelson, M. Ronchetti, Phys. Rev. B 28 (1983) 784. [28] S. Nos, F. Yonezawa, J. Chem. Phys. 84 (1985) 1803. [29] H. Jónsson, H.C. Andersen, Phys. Rev. Lett. 60 (1988) 2295. [30] R.S. Liu, K.J. Dong, Z.A. Tian, H.R. Liu, P. Peng, A.B. Yu, J. Phys.: Condens. Matter 19 (2007) 196103. [31] Z.A. Tian, R.S. Liu, K.J. Dong, A.B. Yu, Europhys. Lett. 96 (2011) 36001. [32] A. Sutton, J. Chen, Phil. Mag. Lett. 61 (1990) 139. [33] Y. Qi, T. Cagin, Y. Kimura, W.A. Goddard, Phys. Rev. B 59 (1999) 3527. [34] P. Xu, T. Cagin, W.A. Goddard, J. Chem. Phys. 123 (2005) 104506. [35] H.J. Lee, T. Cagin, W.L. Johnson, W.A. Goddard, J. Chem. Phys. 119 (2003) 9858. [36] H.C. Andersen, J. Chem. Phys. 72 (1980) 2384. [37] D. Brown, J.H.R. Clarke, Mol. Phys. 51 (1984) 1243. [38] W.G. Hoover, A.J.C. Ladd, B. Moran, Phys. Rev. Lett. 48 (1982) 1818. [39] D.J. Evans, J. Chem. Phys. 78 (1983) 1818. [40] R.D. Mountain, P.K. Basu, J. Chem. Phys. 78 (1983) 7318. [41] D.H. Li, R.A. Moore, S. Wang, J. Chem. Phys. 89 (1988) 4309. [42] D.Q. Yu, M. Chen, X.J. Han, Phys. Rev. E 72 (2005) 051202. [43] Z.A. Tian, K.J. Dong, A.B. Yu, Phys. Rev. E 89 (2014) 032202. [44] A.R. Kansal, S. Torquato, F.H. Stillinger, Phys. Rev. E 66 (2002) 041109. [45] A.D. Bruce, A.N. Jackson, G.J. Ackland, N.B. Wilding, Phys. Rev. E 61 (2000) 906. [46] J.X. Zhu, M. Li, R. Rogers, W. Meyer, R.H. Ottewill, W.B. Russell, P.M. Chaikin, Nature 387 (1997) 883.