Molecular dynamics investigation of the local structure in iron melts and its role in crystal nucleation during rapid solidification

Qi Zhang a, Jincheng Wang *a, Sai Tang b, Yujian Wang a, Junjie Li a, Wenquan Zhou a and Zhijun Wang a
aState Key Laboratory of Solidification Processing, Northwestern Polytechnical University, Xi’an 710072, P. R. China. E-mail: jchwang@nwpu.edu.cn
bPowder Metallurgy Research Institute, Central South University, Changsha, 410083, P. R. China

Received 6th September 2018 , Accepted 19th November 2018

First published on 20th November 2018


A comprehensive investigation on local structures in iron melts and their role in nucleation under various cooling rates was performed by means of large-scale molecular dynamics simulations. The embedded atoms method (EAM) was adopted to describe the interactions between iron atoms. Connections between short-range order (SRO), medium-range order (MRO), and crystalline nucleation from iron melts were constructed using several structural analysis techniques, including the radial distribution function, common neighbor analysis method, the Voronoi tessellation, and bond order analysis. The simulation results showed that abundant types of atomic clusters with SRO, mainly including the icosahedral-like (ICO-like) and fcc-like clusters, were predominant in undercooled iron melts. The obtained microstructures were determined by the competition between the ICO-like and crystal-like configurations. There existed a critical cooling rate, below which the fcc-like configurations gain the advantage upon cooling and where crystallization could take place; otherwise, the ICO-like configurations are favored and the glass phases could be obtained. Furthermore, it was proved that the crystal nucleation could be divided into three stages: first, a fluctuation and competition between crystal-like and ICO-like clusters in undercooled melts; second, the formation and growth of MRO clusters via the transformation of atomic configurations from ICO-like to crystal-like; finally, the nucleation of bcc nuclei from the core of steady MRO clusters. This process agrees with the Ostwald's step rule and the findings from other investigations. Based on the analysis of the compositional origin of MRO clusters, we further found that the MRO clusters were mainly composed of fcc-like instead of ICO-like configurations, indicating a negative role of ICO-like configurations in crystal nucleation.


1. Introduction

Metals and alloys are considered to be among the most fundamental and important materials in modern industries and in daily life. Usually, solid metallic materials are obtained from liquid–solid phase transformations originating from crystal nucleation and a subsequent growth process in undercooled melts. Generally, dynamically and locally structural heterogeneity exists in undercooled melts due to rapid cooling conditions, which could directly influence the nucleation from melts or even avoid crystallization to form an amorphous phase (AP). Furthermore, this controls the morphology of phases, which further determines the physical, chemical, and mechanical properties of the fabricated materials. Therefore, in order to control and predict solidification microstructures with a high accuracy, it is necessary to understand the mechanisms of nucleation and the subsequent microstructural evolution on an atomic scale. Despite its importance, the atomic mechanisms of nucleation and subsequent microstructural evolution are still poorly understood. On the one hand, dendritic growth and subsequent grain growth generally happens on a diffusive time scale and meso-spatio scale, which are much larger than those of crystal nucleation. Therefore, this makes the whole process of solidification a very complex physical phenomenon spanning a broad range of spatiotemporal scales. On the other hand, it is hard to in situ observe the process of nucleation and growth in undercooled melts directly at an atomic level in experiment.

Computational simulations can overcome these drawbacks and have contributed to shedding light on the complex phenomena that occur during solidification. In 1957, Alder and Wainwright1 studied the crystallization process of a hard-sphere system using molecular dynamics (MD) simulations. Since then, numerous simulations have been performed to study the solidification of various systems, such as Lennard-Jones particles, pure metals, and binary alloys. The issues investigated mainly concerned classical nucleation theory (CNT), including the solid–liquid interfacial energy and its anisotropy, the nucleation barrier, and quantitative calculations of the nucleation rates.2–5 However, one rigorous challenge of CNT is that nucleation can be a multi-stage process. In recent years, lots of theoretical,6,7 experimental,8–11 and numerical12–14 studies have indicated that the nucleation of crystallites appear to occur through at least two distinct stages in various systems, due to the formation of metastable precursors with energy much closer to undercooled melts and subsequent nucleation of the most stable crystalline phase, which follows Ostwald's step rule, which greatly deviates from the classical picture described by CNT. Ten Wolde et al.15 found that the precursors were bcc-like in the Lennard-Jones (L-J) system, and later transformed into fcc crystals from the inner of the preordered clusters, while they were still surrounded by bcc shells. Similar results were also obtained by Desgranges16 and Wang.17 By contrast, Swope and Andersen18 found that both bcc and fcc structural nuclei formed at the early stage of crystallization in an L-J system. However, only the fcc nuclei had the chance to grow, while the bcc nuclei finally disappear. Meanwhile, some other researchers declared that non-equilibrium liquids often transform into AP at the early stage of nucleation, and homogeneous crystal nucleation subsequently takes place.12,19 By means of phase field crystal (PFC) simulation, Tang et al.19 suggested that the first formed AP were composed of clusters with short-range order (SRO) and medium-range order (MRO). Also, they suggested that the crystal nucleation begins with the formation of MRO clusters structurally similar to the subsequently nucleated crystal, while SRO clusters play little role in nucleation. Benefitting by the recent progress in high-performance computing, Shibuta et al.20 studied the solidification of iron melts in a billion-atoms system via molecular dynamics simulations. They suggested that the local accumulation of ICO-like structures in undercooled melts near the previously formed grains leads to local heterogeneity in the distribution of grains, which could further promote the formation of new nuclei. Also, Kurtuldu et al.21,22 confirmed experimentally that fcc crystals can form in supercooled liquids by heteroepitaxial growth from an ICO quasicrystals (iQC) template. They suggested that the ICO order is responsible for the nucleation of the fcc phase.

However, there are also abundant works suggesting that ICO configurations could suppress the occurrence of nucleation. Actually, over 60 years ago, Charles Frank23 suggested that the formation of energetically favored fivefold symmetric icosahedra in supercooled liquids might suppress crystallization, leading to vitrification. Inaddition, Leocmach and Tanaka24 proposed that ICO short-range order (ISRO) plays an important role in preventing the crystallization from happening, thus leading to supercooling and glass transition. According to Jade Taffs,25 a fivefold symmetry has both kinetic and thermodynamic effects in suppressing crystallization in hard sphere systems, and the increase in the fivefold symmetry could suppress the kinetics of the homogenous nucleation of fcc crystals.

From the above statements, it is clear that the local atomic clusters in undercooled melts not only provide atoms that are indispensable to the nucleation of crystals, but also determine the pathways of nucleation. Among these local structures, the icosahedral order is ubiquitous in the structures of melts and is the archetypal model of amorphous structures.24,26–28 However, the understanding of the local structures of undercooled melts and their effects on nucleation is still far from complete. In particular, the role of fivefold symmetric ISRO clusters in the nucleation mechanism is still confusing. Since the cooling condition is a critical factor greatly influencing the local structures of metallic melts during solidification, it is of significance to investigate the evolution of atomic clusters under different cooling conditions. In order to obtain a comprehensive understanding of the nucleation mechanism of metallic melts, in this work, we investigated the entire process of microstructural evolution during solidification, from SRO clusters to MRO clusters in iron melts and ultimately the crystalline state.

2. Methods

2.1. Molecular dynamics simulation

The interaction potentials play a determinant role on the accuracy of MD computations. In this work, the embedded atoms method (EAM) potential proposed by Mendelev et al.29 was adopted to describe the interactions between iron atoms. The total energy U is given by:
 
image file: c8cp05654d-t1.tif(1)
where the subscripts i and j denote different atoms, N is the number of atoms in the system, Fi is the embedding energy to embed an atom i into a local site with electron density ρi, and φij(rij) is the pair potential interaction between atom i and atom j. In addition, rij is the distance between atoms i and j. The electron density ρi is given by,
 
image file: c8cp05654d-t2.tif(2)
where ψ(rij) is the contribution to the electron density of atom i from its neighboring atom j. This EAM potential is one of the most established potentials for body-centered-cubic (bcc) metals and by using this, the MD simulations can accurately predict the behaviors of a solid–liquid interface as well as the elastic and thermodynamic properties and even the melting properties of iron.

The simulations were performed with the large-scale atomic/molecular massively parallel simulator (LAMMPS) code30 developed by Sandia National Laboratory and distributed as an open source code. The simulation cell was a cubic box (34.44 × 34.44 × 34.44 nm3) with periodic boundary conditions in the x, y, and z directions, containing 34[thin space (1/6-em)]56[thin space (1/6-em)]000 atoms initially arranged as an ideal bcc crystalline structure. The lattice constant was aF = 2.87 Å. All the MD simulations were carried out in an isothermal–isobaric (NPT) ensemble. The Nosé–Hoover method31,32 was applied to control the temperature and pressure of the modeling system. The velocity-Verlet algorithm was used to integrate the Newton's motion equations with a time step of 2 fs in all the simulations. The master list distance cutoff in our simulations was 7.3 Å.

First, the system was heated from 300 K to 2200 K under a zero pressure so that we could obtain a completely melted configuration and then isothermally relaxed for a supplementary 100 ps at 2200 K to make sure the melt was in equilibrium as determined by the energy change of the system. Subsequently, the melt in equilibrium was cooled down to 300 K under different cooling rates R varied from 9.5 × 1010 K s−1 to 9.5 × 1012 K s−1. The instantaneous spatial coordinates and the thermodynamic data were recorded every 4 ps and 1 ps for post-processing, respectively.

The curve of the heating process and snapshots of the atomic configuration under three characteristic temperatures are shown in Fig. 1. As the temperature increases, the volume of solid phase increases simultaneously due to the formation of vacancies and interstitial pairs. With further increasing the temperature up to 2000 K, the volume changes abruptly, indicating that the solid-to-liquid transition occurs as a consequence of executing large-amplitude vibrations about their equilibrium positions and from the subsequent instability of the atoms. Finally, the complete loss of crystalline order is obtained. Obviously, there exists superheating effects due to the application of 3D periodic boundary conditions to the perfect bcc configuration without various defects, such as free surface, defect clusters, voids, grain boundaries and dislocations, all of which could increase the energy of a solid and lead to a lower free-energy barrier of the melting process. Therefore, the melting temperature of iron in our simulation was higher than the experimental data of 1811 K. This mechanism is reminiscent of Lindemann's criterion and is in agreement with recent experimental observations of melting in superheated aluminum.33–35


image file: c8cp05654d-f1.tif
Fig. 1 The volume (V) variation with temperature (T) during the heating process and snapshots under three characteristic temperature: (a) VT curve, (b) 300 K, (c) 2018 K, (d) 2200 K.

2.2. Methodology of structural analysis

A. Radial distribution function. The radial distribution function (RDF), i.e., g(r), states the probability density to find an atom at a randomly chosen point at distance r from a central atom. It is one of the most popular and effective methods to characterize the structural features of liquids, and amorphous and crystalline solids. For a monatomic system, it can be written as:
 
image file: c8cp05654d-t3.tif(3)
where V is the volume of the simulated system, N is the total number of atoms, and n(r) is the number of atoms around a central atom within the interval between r and r + dr.
B. Voronoi tessellation analysis. The Voronoi tessellation method (VT) can serve as a geometric method to determine the nearest neighbors of a central atom (i.e., its coordination number, CN) by considering the faces of the Voronoi polyhedron (VP) enclosing the atom. Furthermore, the geometric configuration of the VP reflects the characteristic arrangement of near neighbors. Consequently, the VT method is widely used in the structural analysis of various properties of local atomic structures in melts and glasses.

In this method, the Voronoi index 〈n3, n4, n5, n6…〉 was adopted to differentiate different types of polyhedron, where ni is the number of i-edged faces of the VP and image file: c8cp05654d-t4.tif represents the number of nearest-neighbor atoms, i.e., its CN.

C. Bond order analysis. To characterize the local order of any atom and to distinguish between SRO, MRO, and long-range order (LRO) atomic structures during nucleation, we measured the averaged local bond order parameters (averaged-BOP) [q with combining macron]l(i) based on the bond orientational order parameters,11,36,37
 
image file: c8cp05654d-t5.tif(4)
with,
 
image file: c8cp05654d-t6.tif(5)
where, Nb(i) is the number of nearest neighbors of particle i, Ylm(θi,j,ϕi,j) is the set of spherical harmonics of degree l with m ∈ [−l,l], θi,j and ϕi,j are the polar and azimuthal angles of rij = rirj, and ri is the position vector of particle i.

Accordingly, the averaged-BOP can be formulated as follows,

 
image file: c8cp05654d-t7.tif(6)
where,
 
image file: c8cp05654d-t8.tif(7)
Here, l is a free integer parameter and determines the symmetry to which [q with combining macron]l(i) is sensitive, e.g., [q with combining macron]6(i) is sensitive to the local sixfold symmetry in a three-dimensional space. In this work, the parameter [q with combining macron]6(i) was used to identify atoms with SRO, MRO, and LRO. LRO are identified when [q with combining macron]6 ≥ 0.4, SRO are identified when [q with combining macron]6 ≤ 0.28, and other atoms with the value of 0.28 < [q with combining macron]6 < 0.4 are considered as MRO.

D. Adaptive common neighbor analysis. Adaptive common neighbor analysis (a-CNA),38 an extension of the conventional common neighbor analysis (CNA) method,39 is one of the most frequently used structural identification methods to study the solidification of metallic systems.40–42

Given a central atom to be analyzed, the list of Nmax nearest neighbors is first generated by a-CNA method and is sorted by the distance away from the central atom. Nmax is the maximum required number of neighbors for all considered reference structures; for example, Nmax(fcc) = 12 and Nmax(bcc) = 14. Then, a local cutoff radius rcutoff of each atom can be defined according to the average distance of these Nmax nearest neighbors, where those common neighbours which are separated less than the cutoff distance rcutoff are considered as geometrically bonded. Based on this, a-CNA computes a characteristic signature from the topology of bonds that connect the surrounding neighbor atoms. The result of such analysis is then summarized in a triplet ijk, where i is the number of common neighbours, j is the number of geometrical bonds found between the common neighbours, and k is the number of bonds in the longest chain of bonds. By comparing with a set of reference signatures, the structural type of the central atom can ultimately be determined.

Due to the variable cutoff distance employed, the a-CNA method is more suitable to handle a multi-phase system and can distinguish atomic configurations more precisely, such as FCC, BCC, HCP, ICO, and others (unknown coordination structures, such as amorphous liquid and solid, grain boundary, and so on). The obtained atomic configurations were visualized using the open visualization tool (OVITO).43 Within the framework of OVITO, the a-CNA and Voronoi analyses in our work were also achieved by this powerful post-processing software.

3. Results and discussion

3.1. Simulated solidification microstructures under different cooling rates

According to the above-described method, after preparing the initial configuration with high temperature iron melts, we cooled it down to 300 K under different cooling rates R. Fig. 2 demonstrates the variation of total energy (Etotal) of the system vs. temperature (T) during the cooling process. Generally, Etotal decreases with the decrease in temperature for all cooling rates. Nevertheless, there is a critical value of cooling rate, i.e., R0 = 3.2 K ps−1, below which Etotal changes sharply, indicating crystallization of the system. Furthermore, the onset temperature of crystallization becomes lower with the increase in cooling rate, implying the degree of undercooling in melts becomes larger. In contrast, when the cooling rates are above the critical value (R > R0), the energy decreases continuously and the curves resemble each other, indicating that the amorphous transition occurs, while crystallization is suppressed. Furthermore, the faster the cooling rate, the higher the total energy of the final state.
image file: c8cp05654d-f2.tif
Fig. 2 Temperature dependencies of total energy during solidification of iron under different cooling rates.

In order to further testify the effect of the cooling rate on the final microstructures, we utilized the adaptive-CNA (a-CNA) method to visualize simulated results for direct observations of the atomic configurations. Fig. 3 shows snapshots of the obtained microstructures after solidification under different cooling rates. When the cooling rate was 0.095 K ps−1, the obtained microstructure was a bcc single crystal, containing many amorphous atomic clusters acting as point defects (Fig. 3(a)). As the cooling rate increased to 0.95 K ps−1, microstructure of the sample consisted of several grains separated by others configuration atoms, which could be considered as a grain boundary, and thus the polycrystalline structure was achieved (Fig. 3(b)). However, when the cooling rates were above 0.95 K ps−1, the number of grains decreased with the increase in cooling rate (Fig. 3(c)). As the cooling rate further increased, the crystallization process was suppressed and the glassy state was observed (Fig. 3(d–f)). However, it was noticeable that the bcc atoms still existed in the glassy structures even when the cooling rate was as large as 9.5 K ps−1. This was mainly because of the fast nucleation and crystal growth kinetics, resulting in the extremely low glass-forming ability of monatomic metallic liquids and the existence of bcc clusters in the obtained MGs. Thus, the production of pure monatomic MGs in the experiments required an extremely high critical cooling rate of 1013–1014 K s−1.44


image file: c8cp05654d-f3.tif
Fig. 3 Snapshots of the simulated solidification structures under different cooling rates of: (a) 0.095 K ps−1, (b) 0.95 K ps−1, (c) 1.9 K ps−1, (d) 2.5 K ps−1, (e) 4 K ps−1, and (f) 9 K ps−1. The blue, red, green, yellow, and grey balls denote bcc, fcc, hcp, ICO, and the other configurations (unknown coordination structures, such as amorphous liquid and solid, grain boundary, and so on), respectively.

By means of the a-CNA method, the statistics in terms of the number of atoms with different configurations during evolution, including bcc, fcc, hcp, ICO, and the others configurations, were investigated. As shown in Fig. 4, it was found that the number of hcp atoms was slightly larger than those of bcc, fcc, and ICO configurations at the initial stage of the cooling process. With the time elapsing, the population of atoms with an ICO configuration increased at a higher rate compared to other configurations under all cooling rates (Fig. 4(a–e)). This implies that the ICO configuration is much preferred in undercooled iron melts. Under the condition of slower cooling rates (Fig. 4(a–c)), with the further decrease of temperature, the population of ICO atoms kept growing for a while, and then the number of bcc atoms was drastically increased, accompanied by a sharp decrease in ICO atoms at a certain temperature, indicating the formation of bcc nuclei with a critical size at this moment and a catastrophe crystallization. This also shows the strong competition between ICO and bcc configurations during crystallization. In contrast, the number of ICO atoms did not stop increasing until the completion of the liquid–glass transition when the cooling rate was large enough. Interestingly, this clearly shows that the number of fcc and bcc atoms grows suddenly under lower cooling rates (Fig. 4(a and b)); however, this trend gradually disappeared with the increase in cooling rate. Moreover, we observed that the content of bcc atoms in the final structures gradually decreased with the increase in cooling rate, accompanied by the increase in amorphous atoms, indicating the phase transition undergone by the undercooled melts changed from crystallization to amorphization.


image file: c8cp05654d-f4.tif
Fig. 4 Evolution in number of atoms with various configurations under different cooling rates. (a) 0.095 K ps−1, (b) 0.95 K ps−1, (c) 1.9 K ps−1, (d) 2.5 K ps−1, (e) 9.5 K ps−1 (the insets in (a–e) are enlarged ones around zero of the vertical axis).

Despite the a-CNA method being able to provide us information about the local SRO around the given atoms, it failed to indicate whether complete clusters existed. In addition, detailed information on the SRO and MRO configurations could not be unveiled during the cooling process. Considering the limitation of the a-CNA method, it was thus necessary to utilize some other characterization techniques (such as RDF, Voronoi tessellation, and bond order parameters) for performing a more comprehensive and detailed investigation on the microstructural evolution, so that the connections between SRO, MRO, and crystalline nucleation could be constructed.

3.2. Analysis of microstructural evolution by radial distribution function

In order to analyze the structural evolution during the cooling process under different cooling rates, we further adopted the RDF for structural characterization. RDF curves of the final structures under different cooling rates are given in Fig. 5(a). It can be observed that the maximum of the peaks in the RDF curves increase with the decrease in cooling rate, implying that the degree of ordering of the obtained structure becomes more and more intensive, which is in accordance with the critical cooling rate mentioned above. When the cooling rate is smaller than R0 (R < R0), all the peaks in the RDF curves are clear and sharp, indicating crystalline phases are obtained and the crystallization degree of the solidification structures increases with the decrease in cooling rate. Whereas, when the cooling rate is larger than R0 (R > R0), the peaks of RDFs become lower with the increase in cooling rate. In particular, the second peak of the RDFs is split into two sub-peaks. The left sub-peak becomes slightly higher than the right one when the cooling rate is raised to 4 K ps−1, which can be considered as the signature of amorphous transition. Therefore, we can conclude that the critical cooling rate for glass formation of iron is about 4 K ps−1.
image file: c8cp05654d-f5.tif
Fig. 5 RDF curves of the obtained structures at 300 K under different cooling rates (a). RDFs of the crystallization process under the cooling rate of 0.95 K ps−1 (b) and RDFs of the amorphization process under the cooling rate of 9.5 K ps−1 (c).

Furthermore, the RDF curves as a function of temperature during the crystallization and amorphization process are shown in Fig. 5(b and c), respectively, which could provide us with a more intuitive understanding of the structural variation during cooling as well as the distinctions between crystallization and amorphization. We can observe that there are only two noticeable peaks in the RDF curves for liquids at 2200 K, implying the existence of SRO due to the topological dense packing of atoms and the deficiency of MRO and LRO in iron melts at high temperatures. As the temperature decreases down to 1000 K, all the peaks in the RDF curves become increasingly sharp and distinct under both circumstances, indicating the enhancement of SRO and the formation of MRO during this stage. However, the discrepancy clearly appears at the temperature of 800 K. For the cooling rate of 0.95 K ps−1, the peak intensities become much stronger due to the formation of an ordered crystalline phase, while the second peak starts to decompose into two separated peaks when the cooling rate is 9.5 K ps−1 (see the inset of Fig. 5(c)). Also, the positions for the splitting of the second peaks were in the ratios of 1.78 and 1.98 with respect to the first peak position. This implied that the glass transition temperature Tg of iron was about 800 K, which is slightly lower than the results of Liu45 mainly due to the smaller cooling rate used in our work. Additionally, when the system was cooled down to 300 K, the peak intensities in the glassy state were stronger than those in the liquid state, indicating that the glassy state possessed a higher degree of ordering with respect to the liquid state. Also, the second peak of RDF was split into two sub-peaks at the second shell around 4.55 Å, where both sub-peaks were around 4.25 Å and 4.85 Å.

Besides the dependence of the peak amplitudes on temperature, we found that the positions of the peaks were also influenced by the variation of temperature (Fig. 5(b and c)). With the decrease in temperature, the position of the higher order peak gradually moves towards a smaller r. According to the definition of RDF, the ith peak position Ri represents the average distance from the ith nearest-neighbor shell to the centered atom. This implies that the inter-distance between the atoms gets much closer during cooling, resulting in a dense packing and ordering arrangement of atoms. Actually, the absolute value of Ri varies with the chemical composition and atomic size as well. Inspired by the splitting of the second peak in g(r) of MGs and the spherical-periodic order (SPO),46 Liu et al.47 carried out a statistical analysis of the RDF curves and revealed a general feature for different systems. By normalizing all the peak positions Ri (i = 1, 2…5) to R1, they suggested that SRO and MRO are a consequence of the spherical-periodic order (SPO) as well as the local translational symmetry (LTS) of atomic arrangements in MGs. Accordingly, we scrutinized the RDF curves of the undercooled liquid (1500 K) and the glassy state (300 K) in the iron system, as Fig. 6 shows. We calculated the ratios of Ri (i = 1, 2, 3, 4 for undercooled melt and i = 1, 2, …, 5 for glassy state) to the first peak position R1 and the results are shown in Table 1. We found that the values of undercooled melts agreed well with the constants calculated with SPO (1.0, 1.8, 2.6, and 3.4). For iron MGs, the calculated results were basically identical to the statistical results reported by Liu. However, the third peak at R3/R1 with the value of about 1.98 was unable to be predicted by SPO. This comparison suggests that the atomic arrangement in iron undercooled melts can be described by the SPO theory very well, and that the significant difference in atomic structures between the undercooled melts and the obtained glassy structures was the formation of a microstructure with MRO, which was related to the local translational symmetry. Thus, the global atomic packing in iron MG can essentially be described by combining the SPO with LTS. Moreover, this also verifies the accuracy and rationality of our modeling.


image file: c8cp05654d-f6.tif
Fig. 6 RDF curves of the supercooled liquid, crystalline state, and obtained glassy state.
Table 1 Ratios Ri/R1 calculated for undercooled liquid and the obtained glassy state. Values from SPO theory and the statistical values by Liu are also given for comparison
R 1/R1 R 2/R1 R 3/R1 R 4/R1 R 5/R1
Undercooled melt 1.0 1.86 2.71 3.57
MGs (9.5 K ps−1) 1.0 1.73 1.98 2.59 3.53
SPO theory 1.0 1.8 2.6 3.4
Statistical values image file: c8cp05654d-t9.tif image file: c8cp05654d-t10.tif image file: c8cp05654d-t11.tif image file: c8cp05654d-t12.tif image file: c8cp05654d-t13.tif


3.3. Analysis of SRO atomic clusters by Voronoi tessellation method

Due to the limited three-dimensional structural information provided by statistically averaged RDF, we further utilized the VT method to obtain the more detailed local environment of each atom. The calculated populations of the top ten VP types during different cooling processes are shown in Fig. 7. One can see that although there are numerous VP types, certain polyhedrons are always dominant upon cooling. For superheating liquids, various types of VP are well-distributed. The predominant Voronoi polyhedrons (VPs) are 〈0, 3, 6, 4〉, 〈0, 2, 8, 4〉, and 〈0, 3, 6, 5〉, with the percentages mainly ranging from 2% to 4%. This indicates the structural homogeneity of iron liquid above the melting temperature. With the temperature decreasing, the population of these dominant VPs gradually increases. At 1500 K, the fractions of four types of VPs are beyond 4%. It is noticeable that the number of VP 〈0, 0, 12, 0〉 is doubled during this stage. For the iron MG, the predominant VPs are 〈0, 1, 10, 2〉, 〈0, 2, 8, 4〉, 〈0, 3, 6, 4〉, with the percentages of 11.54%, 11.50%, and 11.16%, respectively. While the proportions of another three popular VPs, such as 〈0, 1, 10, 3〉, 〈0, 3, 6, 5〉, 〈0, 0, 12, 0〉, are all between 5% and 6%. We can obviously find that the proportions of 〈0, 1, 10, 2〉, 〈0, 2, 8, 4〉, 〈0, 1, 10, 3〉, and 〈0, 0, 12, 0〉 increase remarkably during the process of glass formation, which is realized through exhausting various low-population VP types with low symmetry. This implies that the atomic structure in undercooled liquid and MG becomes heterogeneous compared to in the superheated liquid state, and the increasing population of atomic clusters with these dominant VP types leads to the first peak in RDF curves getting sharper with the decrease in temperature.
image file: c8cp05654d-f7.tif
Fig. 7 Fractions of the predominant Voronoi polyhedral types during cooling process under the cooling rate of 9.5 K ps−1.

For the sake of a detailed analysis of SRO evolution upon cooling under different cooling rates, we classified these dominant VPs in three categories according to Cheng et al.:48 (a) fcc-like SRO (including 〈0, 3, 6, 4〉, 〈0, 3, 6, 5〉, 〈0, 4, 4, 6〉, 〈0, 4, 4, 7〉), (b) bcc-like SRO with 〈0, 6, 0, 8〉 Voronoi index, and (c) ICO-like SRO (ISRO), including 〈0, 0, 12, 0〉, 〈0, 1, 10, x〉, and 〈0, 2, 8, x〉, x = 1, 2…4. The type of 〈0, 0, 12, 0〉 represents a perfect ICO cluster with a fivefold symmetry, while clusters with a Voronoi index 〈0, 1, 10, x〉, 〈0, 2, 8, x〉 can be considered as distorted ICO clusters. Based on this, we plotted the fractions of these dominant SRO as a function of temperature under cooling rates of 0.95 K ps−1 (crystallization process) and 9.5 K ps−1 (amorphization process) in Fig. 8(a and b), respectively. Initially, the ICO-like and fcc-like clusters had been already populated in the systems due to the pre-existing population of 〈0, 3, 6, 4〉, 〈0, 1, 10, 2〉, and, 〈0, 2, 8, 4〉 in both liquid and glassy states. However, the fraction of bcc-like clusters is negligible (Fig. 8(a1 and b1)). With the decrease in temperature, the latter two types of ICO-like clusters take over and become the most popular under both cooling rates. The increasing population of various atomic clusters with SRO could be responsible for the energy reduction during the early stage of solidification seen in Fig. 2. However, for solidification under different cooling rates, a discrepancy appears when the temperature is close to the crystallization point Tc or amorphization point Tg. When crystallization takes place, the population of fcc-like clusters is always larger than that of ICO-like clusters upon cooling, which means that fcc-like configurations are much more favored by undercooled liquids under small cooling rates (Fig. 8(a2)). With temperature decreasing, all types of atomic clusters reach the maximum in population around Tc. Meanwhile, the clusters with a Voronoi index of 〈0, 6, 0, 8〉, representing a bcc configuration, begin to form catastrophically at the cost of other atomic clusters. Moreover, we can observe that the number of atomic structures with a Voronoi index of 〈0, 4, 4, 6〉 increases sharply before the onset of crystallization, implying its important role in nucleation. It is noticeable that there are many fcc-like clusters still maintained in the sample, probably acting as defects when crystallization is finished, mainly including 〈0, 4, 4, 6〉 and 〈0, 3, 6, 4〉 (Fig. 8(a1)). In contrast, the population of ICO-like clusters grows much more dramatically than that of fcc-like clusters, and the former ones dominate as soon as the temperature is cooled down to 1812 K under the cooling rate of 9.5 K ps−1. After the occurrence of glass transition at Tg, both ICO-like and fcc-like clusters keep increasing and ultimately are retained in the glassy structure, acting as its primary compositions. However, the crystallization of bcc crystals is markedly restrained because of the prevalence of ICO-like configurations and the deficiency of a relaxation time (Fig. 8(b1 and b2)).


image file: c8cp05654d-f8.tif
Fig. 8 Temperature dependence of dominant VPs under the cooling rates of (a) 0.95 K ps−1 and (b) 9.5 K ps−1.

Furthermore, the variation of ICO-like clusters with the decrease in temperature under different cooling rates is shown in Fig. 9. We clearly found that the population of ICO clusters was strongly dependent on temperature, and the maximum fraction grows remarkably with the increase in cooling rates. Meanwhile, there existed a critical cooling rate at about 3.2 K ps−1 as well. When the cooling rate is smaller than 3.2 K ps−1, all the curves start to jump steeply at a certain temperature, indicating the occurrence of crystallization and the transformation of ICO clusters into bcc structures suddenly. Moreover, the slope becomes increasingly steeper with the decrease of cooling rate, implying that the crystallization process happens much faster and easier. On the one hand, the onset temperature of crystallization increases with the decrease in cooling rate, whereby the central atoms of ICO-like clusters have enough activation energy and sufficient relaxation to escape from their sites surrounded by neighboring atoms, i.e., breaking the “cage-trapping” effect.49–51 On the other hand, the population of ICO clusters becomes smaller with the decrease in cooling rate, which could lead to a decrease of the nucleation free-energy barrier of bcc crystal from the iron melts. These two aspects make the bcc structure more favored than the ICO configurations under small cooling rates. Conversely, the fraction of ICO-clusters increases monotonously and the fraction curves tend to overlap each other with the further decrease of temperature when the cooling rate is higher than the critical one. These findings agree well with the results of the RDF analysis and the changing rule of the total energy under different cooling rates, proving that the ICO-like clusters play a negative role in liquid–crystal transitions but a positive role in liquid–glass transformations.


image file: c8cp05654d-f9.tif
Fig. 9 Evolution of icosahedral-like polyhedron during the cooling process under different cooling rates.

Based on the above analysis, we can conclude that the transition from liquid to a bcc crystal and/or glassy state is a result of competition between crystal-like and ICO-like (including perfect and distorted ones) clusters. During the process of solidification under different cooling rates, the dominant VPs are distorted ICO clusters with a Voronoi index of 〈0, 1, 10, 2〉 and 〈0, 2, 8, 4〉 rather than the perfect icosahedra 〈0, 0, 12, 0〉 in many alloy systems. Moreover, this shows that the ICO-like clusters play a role in suppressing crystallization, due to the undercooling of crystallization gradually increasing with the growth in number of ICO-like clusters. However, these SRO clusters could connect together and evolve into MRO clusters upon cooling, which are more popular in undercooled liquid and glassy states.

3.4. Mechanism of crystal nucleation from undercooled iron melts

In order to further investigate the process of nucleation from iron melts in detail, snapshots of atomic configuration during nucleation and subsequent grain growth are also shown in Fig. 10, with the enlarged graphs (b1), (c1), (d1) corresponding to the area inside the red box in (b), (c), (d), respectively. Clearly shown in the enlarged graphs is the existence of atoms with ICO, fcc, and bcc configurations in undercooled melts. With the evolution of time, the population of ICO-like atoms gradually increases, as shown in Fig. 10(b and c). When the temperature decreases down to 1071 K, several bcc nuclei begin to appear in the undercooled melts. This agrees well with the result of Fig. 8(a2), implying the density of ICO-like atoms reaches the maximum value at this moment. Since then, those nuclei with a critical size begin to grow drastically, accompanied by the formation of newborn nuclei from the undercooled melts, ultimately leading to the formation of polycrystalline structures (Fig. 10(c–f)).
image file: c8cp05654d-f10.tif
Fig. 10 Snapshots of the nucleation and subsequent solidification process under the cooling rate of 0.95 K ps−1. (a) 100 ps, (b) 736 ps, (c) 1172 ps, (d) 1236 ps, (e) 1352 ps, (f) 2000 ps. Blue, yellow, red, green and grey atoms represent atoms with bcc, ICO, hcp, fcc and others configurations. (b1), (c1), (d1) are enlarge graphs corresponding to the area enclosed by red pane in (b), (c), (d) respectively, in which those atoms in pink are with the atomic configuration of 〈0, 4, 4, 6〉. Atoms defined as “others” configurations are not shown in (b)–(d) for clarity.

For an in-depth understanding of the mechanisms of crystal nucleation, we utilized the averaged-BOP to scrutinize the nucleation process in a small region enclosed by the red pane in Fig. 10, in which more than 100[thin space (1/6-em)]000 atoms are contained. Fig. 11 shows the formation of bcc nuclei and its subsequent temporal evolution. Initially, the system is full of SRO clusters with a low degree of ordering (Fig. 11(a), blue regions). With the decrease of temperature, the degree of ordering in the system gradually increases. It was found that many atoms with MRO begin to appear and were distributed homogeneously in the whole simulation box (Fig. 11(b), atoms in green and light blue). When the temperature decreases down to 1100 K, atoms with MRO aggregate together, leading to the formation of MRO clusters and local structural heterogeneity in undercooled melts of iron. However, structural fluctuations can be obviously observed upon cooling. Many MRO clusters are not steady thermodynamically and could decompose into atoms with SRO configurations (i.e., atoms enclosed by a yellow circle in Fig. 11(b and c)). With time elapsing, some MRO clusters continually grow, as shown by the temporal evolution of the regions enclosed by red circles in Fig. 11(c–e). When the temperature decreases down to 1071 K, atoms locating at the center of MRO clusters start to get a high degree of ordering via a slight adjustment of their spatial locations, leading to the nucleation of bcc crystals from the inner of MRO clusters (Fig. 11(d and d1)). This adequately indicates that the preordered MRO clusters in undercooled iron melts act as a precursor of bcc crystalline nucleation, which could reduce the crystal–liquid interfacial energy and thus promote crystal nucleation.52 With the temperature further decreasing, SRO atoms in the melts continually attach onto the surface of MRO clusters. Meanwhile, the bcc nuclei keep growing via the consumption of MRO atoms. As the bcc nuclei reach the critical size, they can dramatically grow. When the temperature decreases from 1071 K to 1026 K, the number of atoms belonging to bcc nuclei grows from 57 to 3931 simultaneously (Fig. 11(d1 and e1)). Eventually, the whole simulation region is completely transformed into a bcc crystalline structure, leaving many atoms with MRO or SRO configurations acting as defects. Meanwhile, the positional ordering of atoms with bcc configurations takes place, leading to a higher degree of ordering and a lower free energy of the whole system (Fig. 11(h and i)). This finding in our work agrees well with the study of Ten Wolde in the Lennard-Jones system, further proving that the nucleation of bcc iron from undercooled melts follows the Ostwald's step rule.


image file: c8cp05654d-f11.tif
Fig. 11 Crystal nucleation and subsequent growth processes upon cooling under the cooling rate of 0.95 K ps−1. (a) 100 ps, (b) 1052 ps, (c) 1156 ps, (d) 1172 ps, (e) 1236 ps, (f) 1256 ps, (g) 1284 ps, (h) 1472 ps, (i) 2000 ps. (c1), (d1), (e1) are enlarge graphs corresponding to the area enclosed by red circle in (c), (d), (e) respectively. (c1) shows the MRO cluster, (d1) shows the formation of bcc crystal from the inner of MRO cluster and (e1) shows the growth of bcc nuclei. The color of atoms is determined by the averaged bond orientational order parameter [q with combining macron]6. Atoms with the value of [q with combining macron]6 ≤ 0.28 are defined as SRO (not shown in (b)–(e) for clarity), atoms with the value of 0.28 < [q with combining macron]6 < 0.4 are MRO atoms. The bcc atoms with LRO are identified as those atoms with the value of [q with combining macron]6 ≥ 0.4.

Additionally, we carried out an analysis of the structural composition in the MRO precursor by means of the VT method. We found that the MRO precursor formed at 1100 K (as shown in Fig. 11(c1)) was composed of crystal-like atomic structures, while the ICO-like atomic clusters could hardly be found. These atoms contained in the MRO precursor possessed a similar or identical CN compared to the bcc configuration, such as 〈0, 4, 4, 6〉, 〈0, 3, 6, 4〉, and 〈0, 3, 6, 5〉. Subsequently, we reversely tracked these atoms belonging to the MRO clusters formed at the temperatures of 1100 K and 1071 K (corresponding to the clusters in Fig. 11(d1 and e1) respectively), with the onset temperature selected as 1500 K. By doing this, the compositional origin of the MRO cluster and of the bcc nuclei could be revealed. The results showed that, besides the pre-existing ones, the crystal-like atomic structures contained in the MRO cluster were primarily transformed from distorted ICO-like clusters 〈0, 2, 8, 4〉 and 〈0, 1, 10, 2〉. Fig. 12 shows the temporal evolution of the number of atoms with different atomic structures contained in the MRO cluster within the above-mentioned temperature ranges. Fig. 12(a) clearly shows the competitive relation between ICO-like and crystal-like configurations during the forming process of the MRO precursor. At high temperatures, the number of crystal-like and ICO-like atomic structures fluctuates with temporal evolution, which means that they could transform into each other during this stage. With the temperature decreasing, the ICO-like clusters are more preferred by the undercooled iron melts (from 1000 ps to 1120 ps). However, when the temperature decreases down to 1136 K at 1120 ps, the crystal-like clusters become more favorable. The number of crystal-like clusters dramatically increases while that of ICO-like clusters decreases sharply. This indicates the transformation of atomic configurations from ICO-like to crystal-like, which could contribute to part of the nucleation free-energy barrier and lead to the undercooling of iron melts. Therefore, the formation of the MRO precursor is a result of the competition between ICO-like and crystal-like configurations. The ICO-like clusters play a negative role in the formation of MRO clusters, leading to the effect of suppressing crystal nucleation. Accordingly, we can conclude that the formation of MRO clusters is the precondition of bcc crystalline nucleation other than the accumulation of ICO-like atoms, leading to the nucleation of bcc nuclei taking place in a more ICO-like cluster atmosphere at a low temperature (as shown in Fig. 8(a2)). This fact can also explain the phenomenon described in previous sections that the atomic configuration with a Voronoi index of 〈0, 4, 4, 6〉 increases sharply before the onset of crystallization (the curve displayed in dark brown in Fig. 8(a1)). We can observe directly from the enlarged graphs of Fig. 10 that the number of pink atoms significantly increases before the occurrence of bcc crystal nucleation (as shown in Fig. 10(b1), (c1), and (d1)), which could be considered as a signature for the formation of MRO clusters. Therefore, the crystal-like clusters with SRO, especially for the atomic structure 〈0, 4, 4, 6〉, act as the primary composition of the MRO precursor and play an important role in connecting the structures between undercooled melts and bcc crystals.


image file: c8cp05654d-f12.tif
Fig. 12 Temporal evolution of the number of different atomic structures within the MRO cluster, the temperature ranges from (a) 1500 K to 1100 K and (b) 1500 K to 1071 K.

Moreover, as Fig. 12(b) shows, with the decrease in temperature, the MRO cluster keeps growing and reaches a momentarily steady state containing about 70 iron atoms. Meanwhile, the number of atoms with bcc configurations dramatically increases, indicating that the nucleation of bcc nuclei takes place at 1168 ps (corresponding to the temperature of 1091 K). From the above analysis, we can find that the process of crystal nucleation from iron melts can be divided into three stages: first, the competition between crystal-like and ICO-like clusters; second, the formation of MRO clusters via the transformation of atomic configurations from ICO-like to crystal-like; finally, the nucleation of bcc nuclei from the core of steady MRO clusters. The pathway of crystal nucleation revealed in our work agrees well with the experimental results in colloidal systems. Zhang et al.53 experimentally identified that the initial structure of the nucleus is liquid-like, which can be considered as a precursor, where the core first becomes ordered but the outer layer remains liquid-like during subsequent crystal growth. Lu et al.8 suggested that relatively ordered structures behave as the precursors for crystal nucleation, providing strong evidence for relatively high-ordered structures foreshadowing crystal nucleation. However, the direct experimental evidence in iron is still absent as is the role of ICO configurations in crystal nucleation, which has rarely been revealed by previous atomic simulations.

4. Conclusions

In this study, we carried out MD simulations in the iron system with 3[thin space (1/6-em)]456[thin space (1/6-em)]000 atoms to investigate the locally atomic structures and their effects on nucleation during rapid solidification. Based on several different kinds of structural analysis method, we carried out in-depth investigations of the evolution of microstructures during solidification on different scales and constructed the connections between SRO, MRO, and the final ordered crystal phase.

We found that there existed abundant types of atomic configurations with SRO in undercooled iron melts, among which the ICO-like and fcc-like configurations were prevalent. With the decrease in temperature, the number of predominant atomic configurations continually grows under various cooling rates. However, there existed a critical cooling rate about 4 K ps−1. Below this critical value, the fcc-like configurations take the advantage upon cooling and crystallization could take place. On the contrary, the ICO-like configurations are favored by undercooled melts, leading to the glass phase finally being obtained. The ICO-like atomic configurations play a role in suppressing nucleation due to the competition with crystal-like clusters. With the increase in cooling rate, the obtained structures transformed from single-grain to multi-grain microstructures, and ultimately to amorphous structures when the cooling rate was beyond 4 K ps−1.

We also proved that the nucleation of bcc crystal from undercooled melts follows the Ostwald's step rule. The process of nucleation from iron melts can be divided into three stages: first, the fluctuation and competition between crystal-like and ICO-like clusters in undercooled melts; second, the formation of MRO clusters via the transformation of atomic configurations from ICO-like to crystal-like; finally, the nucleation of bcc nuclei from the core of steady MRO clusters. Furthermore, the MRO clusters were mainly composed of crystal-like SRO configurations, such as 〈0, 4, 4, 6〉, 〈0, 3, 6, 4〉, and 〈0, 3, 6, 5〉, indicating that MRO clusters act as precursors of crystal nucleation from undercooled melts. Our findings may provide new perspectives on the structural origin and atomistic mechanism of the crystalline nucleation and growth from metallic melts.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

The work was supported by the National Natural Science foundation of China (Grants No. 51871183, 51571165). We also thank the Center for High Performance Computing of Northwestern Polytechnical University, China for computer time and facilities.

References

  1. B. J. Alder and T. E. Wainwright, J. Chem. Phys., 1957, 27, 1208–1209 CrossRef CAS.
  2. S. Auer and D. Frenkel, Nature, 2001, 409, 1020–1023 CrossRef CAS PubMed.
  3. J. Bokeloh, R. E. Rozas, J. Horbach and G. Wilde, Phys. Rev. Lett., 2011, 107, 145701 CrossRef CAS.
  4. S. R. Wilson, K. G. S. H. Gunawardana and M. I. Mendelev, J. Chem. Phys., 2015, 142, 134705 CrossRef CAS.
  5. R. L. Davidchack, J. Chem. Phys., 2010, 133, 234701 CrossRef PubMed.
  6. S. Alexander and J. McTague, Phys. Rev. Lett., 1978, 41, 702–705 CrossRef CAS.
  7. W. Klein and F. Leyvraz, Phys. Rev. Lett., 1986, 57, 2845–2848 CrossRef CAS.
  8. Y. Lu, X. Lu, Z. Qin and J. Shen, Solid State Commun., 2015, 217, 13–16 CrossRef CAS.
  9. T. Schilling, H. J. Schöpe, M. Oettel, G. Opletal and I. Snook, Phys. Rev. Lett., 2010, 105, 025701 CrossRef CAS.
  10. H. J. Schöpe, G. Bryant and W. van Megen, Phys. Rev. Lett., 2006, 96, 175701 CrossRef PubMed.
  11. P. Tan, N. Xu and L. Xu, Nat. Phys., 2014, 10, 73–79 Search PubMed.
  12. G. I. Tóth, T. Pusztai, G. Tegze, G. Tóth and L. Gránásy, Phys. Rev. Lett., 2011, 107, 175702 CrossRef.
  13. C. Guo, J. Wang, J. Li, Z. Wang and S. Tang, J. Phys. Lett., 2016, 5008–5014,  DOI:10.1021/acs.jpclett.6b02276,.
  14. G. Díaz Leines, R. Drautz and J. Rogal, J. Chem. Phys., 2017, 146, 154702 CrossRef.
  15. P. R. Ten Wolde, M. J. Ruiz-Montero and D. Frenkel, Phys. Rev. Lett., 1995, 75, 2714–2717 CrossRef CAS.
  16. C. Desgranges and J. Delhommelle, J. Am. Chem. Soc., 2006, 128, 10368–10369 CrossRef CAS.
  17. H. Wang, H. Gould and W. Klein, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2007, 76, 031604 CrossRef PubMed.
  18. W. C. Swope and H. C. Andersen, Phys. Rev. B: Condens. Matter Mater. Phys., 1990, 41, 7042–7054 CrossRef.
  19. S. Tang, J. C. Wang, B. Svendsen and D. Raabe, Acta Mater., 2017, 139, 196–204 CrossRef CAS.
  20. Y. Shibuta, S. Sakane, E. Miyoshi, S. Okita, T. Takaki and M. Ohno, Nat. Commun., 2017, 8, 017 CrossRef.
  21. G. Kurtuldu, A. Sicco and M. Rappaz, Acta Mater., 2014, 70, 240–248 CrossRef CAS.
  22. G. Kurtuldu, P. Jarry and M. Rappaz, Acta Mater., 2013, 61, 7098–7108 CrossRef CAS.
  23. F. C. Frank, Proc. R. Soc. London, Ser. A, 1952, 215, 43–46 CrossRef CAS.
  24. M. Leocmach and H. Tanaka, Nat. Commun., 2012, 3, 974 CrossRef.
  25. J. Taffs and C. Patrick Royall, Nat. Commun., 2016, 7, 13225 CrossRef CAS.
  26. F. Spaepen, Nature, 2000, 408, 781–782 CrossRef CAS.
  27. T. Schenk, D. Holland-Moritz, V. Simonet, R. Bellissent and D. M. Herlach, Phys. Rev. Lett., 2002, 89, 075507 CrossRef CAS.
  28. H. Reichert, O. Klein, H. Dosch, M. Denk, V. Honkimaki, T. Lippmann and G. Reiter, Nature, 2000, 408, 839–841 CrossRef CAS.
  29. M. I. Mendelev, S. Han, D. J. Srolovitz, G. J. Ackland, D. Y. Sun and M. Asta, Philos. Mag., 2003, 83, 3977–3994 CrossRef CAS.
  30. S. J. Plimpton, Fast Parallel Algorithms for Short-Range Molecular Dynamics, 1995.
  31. S. Nosé, J. Chem. Phys., 1984, 81, 511–519 CrossRef.
  32. W. G. Hoover, Phys. Rev. A: At., Mol., Opt. Phys., 1985, 31, 1695–1697 CrossRef.
  33. B. J. Siwick, J. R. Dwyer, R. E. Jordan and R. J. Miller, Science, 2003, 302, 1382–1385 CrossRef CAS.
  34. M. Forsblom and G. Grimvall, Nat. Mater., 2005, 4, 388–390 CrossRef CAS.
  35. A. Samanta, M. E. Tuckerman, T.-Q. Yu and W. E, Science, 2014, 346, 729–732 CrossRef CAS.
  36. P. J. Steinhardt, D. R. Nelson and M. Ronchetti, Phys. Rev. B: Condens. Matter Mater. Phys., 1983, 28, 784 CrossRef CAS.
  37. W. Lechner and C. Dellago, J. Chem. Phys., 2008, 129, 114707 CrossRef.
  38. S. Alexander, Modell. Simul. Mater. Sci. Eng., 2012, 20, 045021 CrossRef.
  39. J. D. Honeycutt and H. C. Andersen, J. Phys. Chem., 1987, 91, 4950–4963 CrossRef CAS.
  40. Y. Shibuta, S. Sakane, T. Takaki and M. Ohno, Acta Mater., 2016, 105, 328–337 CrossRef CAS.
  41. Y. Shibuta and T. Suzuki, J. Chem. Phys., 2008, 129, 144102 CrossRef.
  42. M. Avik, Z. Mohsen Asle and I. B. Michael, Modell. Simul. Mater. Sci. Eng., 2018, 26, 025007 CrossRef.
  43. S. Alexander, Modell. Simul. Mater. Sci. Eng., 2010, 18, 015012 CrossRef.
  44. L. Zhong, J. Wang, H. Sheng, Z. Zhang and S. X. Mao, Nature, 2014, 512, 177–180 CrossRef CAS.
  45. X. J. Liu, Y. Xu, Z. P. Lu, X. Hui, G. L. Chen, G. P. Zheng and C. T. Liu, Acta Mater., 2011, 59, 6480–6488 CrossRef CAS.
  46. P. Häussler, Phys. Rep., 1992, 222, 65–143 CrossRef.
  47. X. J. Liu, Y. Xu, X. Hui, Z. P. Lu, F. Li, G. L. Chen, J. Lu and C. T. Liu, Phys. Rev. Lett., 2010, 105, 155501 CrossRef CAS.
  48. Y. Q. Cheng and E. Ma, Prog. Mater. Sci., 2011, 56, 379–473 CrossRef CAS.
  49. B. Doliwa and A. Heuer, Phys. Rev. Lett., 1998, 80, 4915–4918 CrossRef CAS.
  50. C. Patrick Royall, S. R. Williams, T. Ohtsuka and H. Tanaka, Nat. Mater., 2008, 7, 556 CrossRef CAS.
  51. Z.-Y. Hou, L.-X. Liu, R.-S. Liu, Z.-A. Tian and J.-G. Wang, Chem. Phys. Lett., 2010, 491, 172–176 CrossRef CAS.
  52. T. Kawasaki and H. Tanaka, Proc. Natl. Acad. Sci. U. S. A., 2010, 107, 14036–14041 CrossRef CAS.
  53. T. H. Zhang and X. Y. Liu, Angew. Chem., Int. Ed. Engl., 2009, 48, 1308–1312 CrossRef CAS.

This journal is © the Owner Societies 2019