Home | About | Journals | Submit | Contact Us | Français |

**|**Scientific Reports**|**PMC5543070

Formats

Article sections

Authors

Related links

Sci Rep. 2017; 7: 7194.

Published online 2017 August 3. doi: 10.1038/s41598-017-06482-8

PMCID: PMC5543070

Masato Wakeda, Email: pj.ca.u-akaso.se.em@adekaw.

Received 2017 April 3; Accepted 2017 June 13.

Copyright © The Author(s) 2017

Theoretical prediction of glass forming ability (GFA) of metallic alloys is a key process in exploring metallic alloy compositions with excellent GFA and thus with the ability to form a large-sized bulk metallic glass. Molecular dynamics (MD) simulation is a promising tool to achieve a theoretical prediction. However, direct MD prediction continues to be challenging due to the time-scale limitation of MD. With respect to practical bulk metallic glass alloys, the time necessary for quenching at a typical cooling rate is five or more orders of magnitude higher than that at the MD time-scale. To overcome the time-scale issue, this study proposes a combined method of classical nucleation theory and MD simulations. The method actually allows to depict the time-temperature-transformation (TTT) diagram of the bulk metallic glass alloys. The TTT directly provides a prediction of the critical cooling rate and GFA. Although the method assumes conventional classical nucleation theory, all the material parameters appearing in the theory were determined by MD simulations using realistic interatomic potentials. The method is used to compute the TTT diagrams and critical cooling rates of two Cu-Zr alloy compositions (Cu_{50}Zr_{50} and Cu_{20}Zr_{80}). The results indicate that the proposed method reasonably predicts the critical cooling rate based on the computed TTT.

Metallic glasses possess brilliant properties as structural materials including high elastic limit^{1}, high toughness^{2}, and high corrosion resistance^{3}. However, since reachable size of bulk metallic glass has been limited within centimeters even for highly selected alloy compositions, the application of metallic glasses is substantially restricted as structural materials. Alloys with a lower glass forming ability (GFA) require a higher cooling rate in the melt-quenching process to realize a glass state eventually. In other words, the quenching must be finished prior to the spontaneous nucleation of a critical-sized crystal nucleus in the molten alloy. The critical size is the minimum size required for thermodynamic downhill crystal growth. A critical cooling rate ${\dot{T}}_{\mathrm{c}}$ is defined as the minimum cooling rate that results in a glass state. When we quench a plate-like shaped molten alloy with thickness *L* from very high temperature *T*
_{i} to ambient temperature *T*
_{s}, the spatial distribution of time dependent cooling rate in the molten alloy $\dot{T}(x,t)$ could be derived by solving the 1D heat conduction equation:

$$\frac{\partial T}{\partial t}=\alpha \frac{{\partial}^{2}T}{\partial {x}^{2}},$$

1

under the conditions *T*(*x*, 0)=*T*
_{i} and *T*(0, *t*)=*T*(*L*, *t*)=*T*
_{s}, where *t* denotes time, *x* (0≤*x*≤*L*) denotes coordinate of out-of-plane direction, and *α* denotes thermal diffusivity. Assuming constant thermal diffusivity, the solution is (see Supplemental Information for the details)

$$\dot{T}(x,t)=\frac{4\pi \alpha ({T}_{\mathrm{i}}-{T}_{\mathrm{s}})}{{L}^{2}}\sum _{m=0}^{\mathrm{\infty}}(2m+1)\mathrm{exp}[\phantom{\rule{-.20em}{0ex}}-\phantom{\rule{-.20em}{0ex}}\alpha {\left\{\frac{(2m+1)\pi}{L}\right\}}^{2}t]\mathrm{sin}\left\{\frac{(2m+1)\pi}{L}x\right\}.$$

2

Since $\dot{T}>{\dot{T}}_{\mathrm{c}}$ should be always satisfied throughout the molten alloy sample to obtain a glass state, the critical cooling rate ${\dot{T}}_{\mathrm{c}}$ actually limits the possible sample size *L* as it can be understood from Eq. (^{2}). Thus, ${\dot{T}}_{\mathrm{c}}$ is potentially a good measure of GFA^{4}. The key to obtaining a larger sample size involves searching for alloy composites with a slower critical cooling rate. The critical cooling rate can be estimated after depicting the time-temperature-transformation (TTT) diagram of an alloy^{5}. Therefore, it is necessary to establish a method to predicting the TTT diagram as it can lead to a computational high-throughput screening of alloy composition to obtain a larger bulk metallic glass sample. Furthermore, it is essential to compute the incubation time for the nucleation of the critical crystal nucleus to illustrate the TTT^{5}. Currently, molecular dynamics (MD) simulation is the best tool available for the computation of the incubation time because the critical crystal nucleus typically corresponds to the nanometer range^{6}, and thus it is necessary for the crystal nucleation process to consist of atomic scale events. Additionally, the quantitative analyses of the event properties are possible given that reliable interatomic potentials are provided. Recently, a study demonstrated solidification from melting by using direct MD quenching simulations^{7, 8}. This is followed by estimating the TTT and critical cooling rate from direct MD observation of the incubation time. However, these studies only focused on highly simplified model materials, such as the Lennard-Jones system, which possess a very high critical cooling rate because of the time scale limitations of MD simulations. Conversely, the composition of target alloys must possess the potential to generate a centimeter sized metallic glass with a very long incubation time and a very slow critical cooling rate when compared to simple system models. Thus, it is necessary for the incubation time to considerably exceed typical MD time scales such as microseconds. Thus, a direct MD simulation does not work well due to the fore-mentioned limitation. Therefore, in this study, classical nucleation theory is employed to compute the incubation time as opposed to the direct MD simulations. However, the intrinsic parameters of the materials that appear in classical nucleation theory are determined using MD simulations such as the free energy difference between a melt and crystal and the interfacial free energy between a melt and crystal. Finally, the TTT diagram is depicted, and the critical cooling rate is evaluated based on the TTT diagram.

Two alloy compositions possessing different GFA, such as Cu_{50}Zr_{50} and Cu_{20}Zr_{80}, are examined in the study to demonstrate the proposed method.

Firstly, the critical nucleation radius *r** is determined as follows. A cubic simulation cell is set with a periodic boundary condition (PBC) containing *N* atoms, which consist of supercooled liquid (melt) and a spherical crystal nucleus of a radius *r* (see Fig. ^{S1} and Tables ^{S1} and ^{S2} for the details). The Finnis-Sinclair (FS) potential^{9} is used to describe the interatomic interaction for both alloy systems. The B2 crystal structure is assumed as the crystal structure of Cu_{50}Zr_{50}
^{10}, while a distorted BCC like structure determined by MD quenching simulation is used as the crystal structure of Cu_{20}Zr_{80}. The melts of both alloys are prepared by annealing at 2,500K for 100ps under a zero pressure condition. The crystal nucleus is then embedded into the melt, while omitting all atoms in the melt that overlap with the embedded crystal nucleus. The systems are subsequently quenched at 10K for 2ps to fill the space gap at the interface between the melt and the crystal nucleus. The MD time step of 2fs is used throughout the study for all the MD simulations. The constructed melt-crystal models are employed as the initial configuration to determine the temperature dependent critical nucleus radius *r**(*T*). The NPT ensemble MD simulations are performed with respect to the models with different nucleus radius at different temperatures under a zero pressure condition to determine the critical temperature as the middle point between the highest nucleus growth and the lowest nucleus shrinking temperatures (see Fig. ^{S2}). Ten MD simulations are performed for each condition to reduce the statistical error. Actual growing and shrinking behaviors of the Cu_{50}Zr_{50} crystal nucleus of *r*=2.0nm immediately above and below the critical temperature *T*=1,175K are shown in Fig. 1. The atoms in the melt and crystal are identified and colored using bond-order analysis^{11}. The relationship between temperature and inverse of radius is summarized in Fig. ^{S2}, which clearly shows a trend in which the inverse of critical radius approximately linearly increases with increases in the temperature. Figure 2 represents the relationship between supercooled degree Δ*T*=*T*
_{m}−*T* and the inverse of *r**. The results indicate that Δ*T* is reasonably assumed as proportional to the inverse of *r**:

$$\mathrm{\Delta}T={T}_{\mathrm{m}}-T\simeq \frac{k}{{r}^{\u204e}},$$

3

where the values of *k* are estimated as −3.29×10^{2}nm/K (Cu_{50}Zr_{50}) and −4.66×10^{2}nm/K (Cu_{20}Zr_{80}), respectively. The melting temperatures *T*
_{m} are determined from the ordinate intercept of the linear fitting line for critical temperatures in Fig. ^{S2}. Since Eq. (^{3}), *T*
_{m} is assumed as the critical temperature at *r**→∞. The obtained melting temperature is shown in Table 1. The melting temperature of Cu_{50}Zr_{50} is lower than that of Cu_{20}Zr_{80} by 100K. The difference is consistent with an experiment involving Cu-Zr alloys^{12} while the actual melting temperature of Cu_{50}Zr_{50}, *T*
_{m}=1,208K^{13} is slightly lower (by 134K) than the computed value. The discrepancy can be mainly attributed to an underestimation of melting temperature subject to FS potential energetics. With respect to Cu_{50}Zr_{50}, the melting temperature is also estimated using a melt-crystal biphase model with PBC and a flat interface that is perpendicular to the (100) crystal plane (see Fig. ^{S3}). Determination of temperature involves an immobile interface under NPH ensemble MD simulation at a zero pressure condition^{14} to provide the melting temperature. Although the melting temperature determined using biphase model can depend on interfacial crystal orientation, both are in reasonable agreement with each other as shown in Fig. ^{S2(a)}.

Snapshots of crystal nucleus growth and shrink processes of Cu_{50}Zr_{50} model. (**a**) denotes the growth process (*T*=1,100K, *r*=2.0nm), and (**b**) denotes the shrink process (*T*=1,200K, **...**

The free energy barrier of crystal nucleation Δ*G** is computed using Eq. (^{16}). The interfacial free energy *σ* is computed using Eqs. (^{14}) and (^{15}). The latent heat per unit volume Δ*H*
_{m} in Eq. (^{15}) is evaluated from the enthalpy difference between the melt and crystal parts in the melt-crystal biphase model (Fig. ^{S3}). The melt and crystal models used in the biphase model are obtained by annealing the melt for 100ps at 2,500K and relaxing the crystal structure for 100ps at 0.1K under a zero pressure condition, respectively. The two structures are subsequently attached as shown in Fig. ^{S3}, and 600ps MD relaxation is performed at the melting temperature *T*
_{m} and zero pressure condition. The computed Δ*H*
_{m} values are shown in Table 1. The findings reveal that Δ*H*
_{m} of Cu_{50}Zr_{50} is almost twice that of Cu_{20}Zr_{80}. The temperature dependent enthalpy is obtained by averaging the total energy over a 1ns NPT ensemble MD simulation using a *N*=16,000 atoms simulation cell with PBC at different temperatures ranging from 50K to *T*
_{m} with 50K intervals and a zero pressure condition. The obtained temperature dependent enthalpies of the melt and crystal are fitted by the following quadratic polynomial function (see Fig. ^{S4}):

4

The volumetric specific heats *C*
_{p} of the melt and crystal are then obtained as the derivative of enthalpy (Eq. (^{4})) with respect to the temperature at constant pressure as follows:

$${C}_{\mathrm{p}}={\left(\frac{\partial H}{\partial T}\right)}_{\mathrm{p}}\simeq 2aT+b\mathrm{.}$$

5

Additionally, Δ*C*
_{p} in Eq. (^{14}) is estimated by $\mathrm{\Delta}{C}_{\mathrm{p}}={C}_{\mathrm{p}}^{\mathrm{melt}}-{C}_{\mathrm{p}}^{\mathrm{crystal}}$. The parameters determined in Eqs. (^{4}) and (^{5}) are summarized in Table ^{S3}. Subsequently, the temperature dependent interfacial free energy *σ* is computed and shown in Fig. 3. The interfacial free energy of Cu_{50}Zr_{50} is consistently higher than that of Cu_{20}Zr_{80}.

Eventually, the free energy barrier of crystal nucleation Δ*G** is computed as a function of temperature and shown in Fig. 4.

The temperature dependent nucleation rate *J* and the incubation time *t* are directly computed from Eq. (^{8}). The atomic number density of the melt *ρ*
_{melt} and that of the critical crystal nucleus ${\rho}_{\mathrm{crystal}}^{\u204e}$ are computed and summarized in Table 2 by taking the average volume of simulation cell *V* (*N*=16,000) with PBC over 100ps NPT ensemble MD simulation at the critical temperature and a zero pressure condition. The Zeldovich factor *Z* is estimated by using Eqs. (^{18}) and (^{19}). With the exception of the parameters obtained in the previous sections, the attachment rate *f*
^{+} that appears in the pre-exponential factor *J*
_{0} (Eq. (^{19})) is still unknown, and thus it is necessary to obtain *J*. In order to compute the attachment rate, 4ns NPT ensemble MD simulations on the same melt - spherical crystal model with *r*=*r**(*T*) used in the previous analysis are performed at different temperatures ranging from 1,060 to 1,160K under a zero pressure condition. The number of atoms belonging to the embedded crystal, Δ*n**(*t*)=*n**(*t*)−*n**(0), are recorded during the MD simulations. This is followed by obtaining the MSD Δ*n**^{2}(*t*) and attachment rate (Eq. (^{20})). It should be noted that ten independent MD simulations are performed for each temperature, and they are averaged to reduce the statistical error. Figure 5 shows time evolution of the MSD for Cu_{50}Zr_{50} at *T*=1,155K. The slope of the MSD-time plot corresponds to twice that of the attachment rate, 2*f*
^{+}(*n**). The computed attachment rates are plotted in Fig. 6 for Cu_{50}Zr_{50} and are summarized for both alloys in Table 2. The attachment rate increases with increases in the temperature. The attachment rate of Cu_{20}Zr_{80} is an order of magnitude higher than that of Cu_{50}Zr_{50}. The computed attachment rate is used to eventually compute the nucleation rate *J* of unit volume and summarize it in the Table 2. The nucleation rate rapidly increases with increases in the temperature. The rapid change in the nucleation rate is mainly attributed to the exponential term in Eq. (^{19}) as opposed to the pre-exponential factors. As previously stated, the direct MD analyses of the attachment rate *f*
^{+}(*n**) are only performed above *T*=1,060K because a large statistical error due to the lower number of samples is expected below this temperature. In order to predict the attachment rate at low temperatures, the following Arrhenius type exponential function fitted to high temperature MD data of the attachment rate is used instead of performing direct MD analysis at the low temperatures.

6

where *A* and *B* denote fitting parameters. The values of *A* and *B* are summarized in Table ^{S4}. It is worth noting that recent theoretical crystallization study for a Lennard-Jones system^{15} shows a super-Arrhenius behavior of the attachment rate. Actually our attachment rate data of Cu_{20}Zr_{80} also seem to exhibit super-Arrhenius like behavior, while data of Cu_{50}Zr_{50} do not clearly exhibit at least within the examined temperature range (see Fig. ^{S5}). However, since the deviation from Arrhenius behavior is not significant, here the Arrhenius behavior was assumed for simplicity.

Attachment rate - temperature relationship (Cu_{50}Zr_{50}). The dashed line represents the fitted curve (Eq. (^{6}) and Table ^{S4}).

The TTT diagrams for Cu_{50}Zr_{50} and Cu_{20}Zr_{80} are depicted with respect to the temperature dependent incubation time *t*=1/*J* using the parameters obtained in previous analysis and shown in Fig. 7. The critical cooling rate ${\dot{T}}_{\mathrm{c}}$ is approximately calculated from the following equation:

$${\dot{T}}_{\mathrm{c}}\simeq \frac{{T}_{\mathrm{m}}-{T}_{\mathrm{n}}}{{t}_{\mathrm{n}}},$$

7

where *T*
_{n} and *t*
_{n} denote the nose temperature of the TTT diagram and incubation time at the nose temperature, respectively. Table 3 shows the estimated critical cooling rate for a local volume of Ω=343nm^{3}. Fortunately, with respect to Cu_{20}Zr_{80}, the crystallization can be observed even in a direct MD quenching simulation due to its very low incubation time. Therefore, the obtained critical cooling rate using the proposed method is compared with that using direct MD simulation for the same average volume of the simulation cell. In order to directly estimate the critical cooling rate from MD, 1ns NPT ensemble MD simulations are performed using *N*=16,000 atoms simulation cell with PBC at different temperatures ranging from 800K to *T*
_{m} with 50K intervals. The incubation time for each temperature is estimated as an incubation time for a rapid reduction of system potential energy, which is a sign of crystallization. The nose temperature corresponds to *T*
_{n}=950K, and the incubation time corresponds to *t*
_{n}=4.5×10^{−10}s. The average volume of the simulation cell at *T*
_{n} corresponds to Ω=343nm^{3}. The critical cooling rates obtained in the direct MD and experiment are also shown in Table 3. The experimental critical cooling rate ${\dot{T}}_{\mathrm{c}}$ is approximately estimated using Eq. (^{2}) as the slowest cooling rate throughout the sample when the local temperature $T(x,{t}_{\mathrm{n}}(x))={T}_{\mathrm{n}}\simeq ({T}_{\mathrm{m}}+{T}_{\mathrm{g}})/2;\phantom{\rule{thinmathspace}{0ex}}\phantom{\rule{thinmathspace}{0ex}}{\dot{T}}_{\mathrm{c}}={min}_{x\in L}(\dot{T}(x,{t}_{\mathrm{n}}(x)))$, *t*
_{n}(*x*): time when *T*(*x*, *t*
_{n}(*x*))=*T*
_{n} at *x*, with *L*=2mm^{16}, *α*=2.5×10^{−2}cm^{2}/s^{17}, *T*
_{m}=1,208K^{13}, *T*
_{m}−*T*
_{g}=400K^{17}, *T*
_{i}=1,500K (>*T*
_{m}), *T*
_{s}=300K (room temperature) and assuming that the volume with the slowest cooling rate in the system is the same as the volume of the simulation cell (Ω=343nm^{3}) (see Fig. ^{S6}). Note that the order of the slowest cooling rate has a weak dependence on the temperature *T*. The results in Table 3 indicate that the proposed method can reasonably reproduce the critical cooling rate obtained by the direct MD simulation and experiment. The computed critical cooling rate of Cu_{50}Zr_{50} is significantly lower than that of Cu_{20}Zr_{80}. This corresponds to the fact that Cu_{50}Zr_{50} does actually form a bulk metallic glass at the usual experimental cooling rate^{16}. In contrast, Cu_{20}Zr_{80} does not form a bulk metallic glass even at highest possible cooling rate in the experiment. Hence, the proposed classical nucleation theory - MD combined method allows the prediction of the critical cooling rate.

It should be noted that a homogeneous nucleation is assumed in this study. However, in practice, the nucleation is mostly heterogeneous due to the existence of impurities and/or the wall surfaces of containers. Thus, the proposed method underestimates the critical cooling rate. Nevertheless, it is necessary to adequately reproduce the GFA ranking, which is the most important factor in the high-throughput alloy design. It is worth noting that the target crystal structure and interatomic potential function are usually unknown for less familiar alloys. Although determination of them is non-straightforward task, recent development of crystal structure prediction methods^{18, 19} and interatomic potential construction methods^{20, 21} from first-principles may permit the task.

In this study, an atomistically-informed method is proposed to predict a TTT diagram of alloys over a wide temperature range by combining classical nucleation theory and MD simulations. The proposed method is used to depict the TTT diagram and critical cooling rate of two Cu-Zr alloys, namely Cu_{50}Zr_{50} and Cu_{20}Zr_{80}. The results indicate that the method reasonably reproduces the GFA ranking and the critical cooling rate by comparing the direct MD simulation and experimental knowledge. Hence, it is expected that the method can open up a computational high-throughput screening of higher GFA alloys.

The incubation time *t* corresponds to the inverse of nucleation rate *J* for a certain piece of volume Ω^{22} of melt, which is subject to the Arrhenius equation^{23} as follows:

$$J=\mathrm{\Omega}{J}_{0}\mathrm{exp}\left(-\frac{\mathrm{\Delta}{G}^{\u204e}}{{k}_{\mathrm{B}}T}\right),$$

8

where *k*
_{B} denotes the Boltzmann constant, and *T* denotes temperature of the supercooled melt. Additionally, Δ*G** denotes the free energy barrier of the nucleation process of a critical sized crystal nucleus. Subject to the classical nucleation theory, Δ*G** is formulated as follows.

A change in the Gibbs free energy during the crystal growth process Δ*G* is expressed as a function of crystal nucleus radius *r* and interfacial free energy *σ* subject to the approximation of spherical nucleus shape^{24} as follows:

$$\mathrm{\Delta}G=-\phantom{\rule{-.25em}{0ex}}\frac{4\pi {r}^{3}}{3}\mathrm{\Delta}{G}_{\mathrm{v}}+4\pi {r}^{2}\sigma ,$$

9

where Δ*G*
_{v} denotes free energy difference per unit volume between the melt and crystal. Extant studies proposed several approximate expressions of Δ*G*
_{v}
^{25, 26}. In this study, it is directly computes by using atomistically computed isobaric heat capacity per unit volume (volumetric specific heat) of melt ${C}_{\mathrm{p}}^{\mathrm{melt}}$ and crystal ${C}_{\mathrm{p}}^{\mathrm{crystal}}$. The Δ*G*
_{v} is given as follows:

Δ*G*_{v} = Δ*H* − *T*Δ*S*,

10

where Δ*H* and Δ*S* denote enthalpy and entropy differences per unit volume between the melt and crystal, respectively. The Δ*H* and Δ*S* can be formulated using the isobaric volumetric specific heat difference between the melt and crystal, $\mathrm{\Delta}{C}_{\mathrm{p}}={C}_{\mathrm{p}}^{\mathrm{melt}}-{C}_{\mathrm{p}}^{\mathrm{crystal}}$ as follows:

$$\mathrm{\Delta}H=\mathrm{\Delta}{H}_{\mathrm{m}}-{\int}_{T}^{{T}_{\mathrm{m}}}\mathrm{\Delta}{C}_{\mathrm{p}}\mathrm{d}T,$$

11

$$\mathrm{\Delta}S=\mathrm{\Delta}{S}_{\mathrm{m}}-{\int}_{T}^{{T}_{\mathrm{m}}}\frac{\mathrm{\Delta}{C}_{\mathrm{p}}}{T}\mathrm{d}T,$$

12

where Δ*H*
_{m} denotes melting enthalpy (latent heat), and Δ*S*
_{m} denotes melting entropy. These are related to each other as follows:

$$\mathrm{\Delta}{S}_{\mathrm{m}}=\frac{\mathrm{\Delta}{H}_{\mathrm{m}}}{{T}_{\mathrm{m}}}\mathrm{.}$$

13

Substituting Eqs. (^{11}) and (^{13}) into Eq. (^{10}), the following expression is obtained:

$$\mathrm{\Delta}{G}_{\mathrm{v}}=\mathrm{\Delta}{H}_{\mathrm{m}}\left(1-\frac{T}{{T}_{\mathrm{m}}}\right)-{\int}_{T}^{{T}_{\mathrm{m}}}\mathrm{\Delta}{C}_{\mathrm{p}}\mathrm{d}T+T{\int}_{T}^{{T}_{\mathrm{m}}}\frac{\mathrm{\Delta}{C}_{\mathrm{p}}}{T}\mathrm{d}T\mathrm{.}$$

14

It should be noted that the crystal nucleus shape is not exactly spherical as described by like Wulff construction^{27} because the interfacial energy is typically anisotropic and additionally, the interface should not very smooth at the atomic level, and thereby it should not be possible to precisely describe the interface by a smooth function especially with respect to a very small crystal nucleus. Specifically, this fact is observed in atomic simulations^{28} and also observed in MD simulations (see Fig. ^{S1}). Nevertheless, the spherical approximation is sufficient for the incubation time estimation as demonstrated later in the study.

In the right hand side of Eq. (^{9}), the two terms compete with each other, and this leads to a crossover with respect to the crystal nucleus radius *r*. The crystal nucleus tends to shrink at a small *r* (dΔ*G*/d*r*>0) while the crystal nucleus tends to grow spontaneously at a large *r* (dΔ*G*/d*r*<0). The first term decreases in proportion to *r*
^{3} while the second increases in proportion to *r*
^{2} with increases in *r*. Therefore, Δ*G* is maximized at a critical nucleus radius *r*=*r**, which should satisfy the following:

$${\left(\frac{\mathrm{d}\mathrm{\Delta}G}{\mathrm{d}r}\right)}_{r={r}^{\ast}}=-\phantom{\rule{-.15em}{0ex}}4\pi {r}^{\ast 2}\mathrm{\Delta}{G}_{\mathrm{v}}+8\pi {r}^{\ast}\sigma =0.$$

15

The free energy barrier, which corresponds to energy at the critical nucleus radius, is given as follows:

$$\mathrm{\Delta}{G}^{\ast}=\mathrm{\Delta}{G}_{r={r}^{\ast}}=\frac{4}{3}\pi \sigma {r}^{\ast 2}.$$

16

It is assumed that the nucleus growth and shrink are achieved by attaching atoms to the nucleus, and thus the pre-exponential factor *J*
_{0} in Eq. (^{8}) is expressed as follows^{29}:

17

$$Z=\sqrt{-\frac{1}{2\pi {k}_{\mathrm{B}}T}{\left(\frac{{\partial}^{2}\mathrm{\Delta}G}{\partial {n}^{2}}\right)}_{n={n}^{\u204e}}},$$

18

$${\left(\frac{{\mathrm{\partial}}^{2}\mathrm{\Delta}G}{\mathrm{\partial}{n}^{2}}\right)}_{n={n}^{\ast}}=-\frac{8}{27}\pi \sigma {\left(\frac{3}{4\pi {\rho}_{\mathrm{c}}^{\ast}}\right)}^{\frac{2}{3}}{n}^{\ast -\frac{4}{3}},$$

19

where *ρ*
_{m} denotes the number density of the melt, *n* denotes number of atoms in crystal nucleus, and ${n}^{\u204e}={n}_{r={r}^{\u204e}}=\frac{4}{3}\pi {r}^{\mathrm{\u204e3}}{\rho}_{\mathrm{c}}^{\u204e}$. Additionally, ${\rho}_{\mathrm{c}}^{\u204e}$ denotes atomic number density of the critical crystal nucleus. Furthermore, *f*
^{+}(*n**) denotes the attachment rate of atoms to the critical crystal nucleus, and *Z* denotes the Zeldovich factor^{30}. The attachment rate *f*
^{+}(*n**) can be expressed as^{31} follows:

$${f}^{+}({n}^{\u204e})=\frac{1}{2}\frac{\langle \mathrm{\Delta}{n}^{\mathrm{\u204e2}}(t)\rangle}{t},$$

20

where Δ*n**^{2}(*t*) denotes the mean square deviation (MSD) of the atoms attached to the critical crystal nucleus during a time interval *t*. In the study, the nucleation rates *J* and incubation time *t*=1/*J* are estimated using Eq. (^{8}) with Eq. (^{17}) using the atomistically determined Δ*G** and *f*
^{+}(*n**). The TTT diagram is then depicted using the obtained incubation time *t* at different temperatures below the melting temperature *T*
_{m}.

The authors thank Susumu Goto (Osaka University) for fruitful discussions. This work was supported by the Elements StrategyInitiative for Structural Materials (ESISM) and Grant-in-Aid for Young Scientists (A) (No. 17H04949) and Scientific Research (A) (No. 17H01238).

Author Contributions

Y.S. and C.N. conducted the numerical simulations. M.W. and S.O. designed this study. Y.S., M.W., and S.O. analyzed the results and drafted the manuscript. All authors contribute to discussion of the results.

The authors declare that they have no competing interests.

**Electronic supplementary material**

**Supplementary information** accompanies this paper at doi:10.1038/s41598-017-06482-8

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Masato Wakeda, Email: pj.ca.u-akaso.se.em@adekaw.

Shigenobu Ogata, Email: pj.ca.u-akaso.se.em@atago.

1. Masumoto T, Maddin R. The mechanical properties of palladium 20% silicon alloy quenched from the liquid state. Acta metallurgica. 1971;19(7):725–741. doi: 10.1016/0001-6160(71)90028-9. [Cross Ref]

2. Kimura H, Masumoto T. Fracture toughness of amorphous metals. Scripta Metallurgica. 1975;9:211–221. doi: 10.1016/0036-9748(75)90196-9. [Cross Ref]

3. Asami K, Kawashima A, Hashimoto K. Chemical properties and applications of some amorphous alloys. Material Science and Engineering. 1988;99:475–481. doi: 10.1016/0025-5416(88)90380-1. [Cross Ref]

4. Inoue A. High Strength Bulk Amorphous Alloys with Low Critical Cooling Rates (Overview) Materials Transactions. 1995;36:866–875. doi: 10.2320/matertrans1989.36.866. [Cross Ref]

5. Shibuta Y, Oguchi K, Takaki T, Ohno M. Homogeneous nucleation and microstructure evolution in million-atom molecular dynamics simulation. Scientific Reports. 2015;5:13534. doi: 10.1038/srep13534. [PMC free article] [PubMed] [Cross Ref]

6. Qi Y, Cagin T, Johnson WL, Goddard WA., III Melting and crystallization in Ni nanoclusters: The mesoscale regime. Journal of Chemical Physics. 2001;115(1):385–394. doi: 10.1063/1.1373664. [Cross Ref]

7. Peng Y, et al. Two-step nucleation mechanism in solid-solid phase transitions. Nature Materials. 2015;14:101–108. doi: 10.1038/nmat4083. [PubMed] [Cross Ref]

8. Shibuta Y, Sakane S, Takaki T, Ohno M. Submicrometer-scale molecular dynamics simulation of nucleation and solidification from undercooled melt: Linkage between empirical interpretation and atomistic nature. Acta Materialia. 2016;105:328–337. doi: 10.1016/j.actamat.2015.12.033. [Cross Ref]

9. Mendelev MI, et al. Development of suitable interatomic potentials for simulation of liquid and amorphous Cu-Zr alloys. Philosophical Magazine. 2009;89:967–987. doi: 10.1080/14786430902832773. [Cross Ref]

10. Arias D, Abriata JP. Cu-Zr (Copper-Zirconium) Journal of Phase Equilibria. 1990;11(5):452–459. doi: 10.1007/BF02898260. [Cross Ref]

11. Lechner W, Dellago C. Accurate determination of crystal structures based on averaged local bond order parameters. Journal of Chemical Physics. 2008;129:114707. doi: 10.1063/1.2977970. [PubMed] [Cross Ref]

12. Zhou SH, Napolitano RE. Phase stability for the Cu-Zr system: First-principles, experiments and solution-based modeling. Acta Materialia. 2010;58:2186–2196. doi: 10.1016/j.actamat.2009.12.004. [Cross Ref]

13. Lundin CE, McPherson DJ, Hansen M. System zirconium-copper. AIME Transactions. 1953;197:273–278.

14. Shibuta Y, Takamoto S, Suzuki T. A molecular dynamics study of the energy and structure of the symmetric tilt boundary of iron. ISIJ International. 2008;48(11):1582–1591. doi: 10.2355/isijinternational.48.1582. [Cross Ref]

15. Malek SMA, Morrow GP, Saika-Voivod I. Crystallization of Lennard-Jones nanodroplets: From near melting to deeply supercooled. The Journal of Chemical Physics. 2015;142:124506. doi: 10.1063/1.4915917. [PubMed] [Cross Ref]

16. Wang WH, Lewandowski JJ, Greer AL. Understanding the Glass-forming Ability of Cu_{50}Zr_{50} Alloys in Terms of a Metastable Eutectic. Journal of Materials Research. 2005;20:2307–2313. doi: 10.1557/jmr.2005.0302. [Cross Ref]

17. Lin XH, Johnson WL. Formation of Ti-Zr-Cu-Ni bulk metallic glasses. Journal of Applied Physics. 1995;78(11):6514–6519. doi: 10.1063/1.360537. [Cross Ref]

18. Wang Y, Lv J, Zhu L, Ma Y. CALYPSO: A method for crystal structure prediction. Computer Physics Communications. 2012;183(10):2063–2070. doi: 10.1016/j.cpc.2012.05.008. [Cross Ref]

19. Lyakhov AO, Oganov AR, Stokes HT, Zhu Q. New developments in evolutionary structure prediction algorithm USPEX. Computer Physics Communications. 2013;184(4):1172–1182. doi: 10.1016/j.cpc.2012.12.009. [Cross Ref]

20. Jelinek B, et al. Modified embedded atom method potential for Al, Si, Mg, Cu, and Fe alloys. Physical Review B. 2012;85:245102. doi: 10.1103/PhysRevB.85.245102. [Cross Ref]

21. Cui Z, Gao F, Cui Z, Qu J. A second nearest-neighbor embedded atom method interatomic potential for Li-Si alloys. Journal of Power Sources. 2012;207:150–159. doi: 10.1016/j.jpowsour.2012.01.145. [Cross Ref]

22. Dantzig, J. A. and Rappaz, M. Solidification. EPFL Press and CRC Press, Lausanne (2009).

23. Erdemir D, Lee AY, Myerson AS. Nucleation of Crystals from Solution: Classical and Two-Step Models. Accounts of Chemical Research. 2009;42(5):621–629. doi: 10.1021/ar800217x. [PubMed] [Cross Ref]

24. Kurz, W. and Fisher, D. J. Fundamental of Solidification 4th revised edition. Trans Tech Publication, Aedermannsdorf, 21 (1998).

25. Turnbull D. Formation of Crystal Nuclei in Liquid Metals. Journal of Applied Physics. 1950;21:1022–1028. doi: 10.1063/1.1699435. [Cross Ref]

26. Thompson CV, Spaepen F. On the approximation of the free energy change on crystallization. Acta Metallurgica. 1979;27:1855–1859. doi: 10.1016/0001-6160(79)90076-2. [Cross Ref]

27. Dobrushin, R. L., Kotecky, R. and Shlosman, S. Wulff construction: a global shape from local interaction. American Mathematical Society, Providence, Rhode Island (1992).

28. Lechner W, Dellago C, Bolhuis PG. Role of the Prestructured Surface Cloud in Crystal Nucleation. Physical Review Letters. 2011;106:085701. doi: 10.1103/PhysRevLett.106.085701. [PubMed] [Cross Ref]

29. Auer S, Frenkel D. Prediction of absolute crystal-nucleation rate in hard-sphere colloids. Nature. 2001;409(6823):1020–1023. doi: 10.1038/35059035. [PubMed] [Cross Ref]

30. Kelton KF. Crystal nucleation in liquids and glasses. Solid State Physics. 1991;45:75–177. doi: 10.1016/S0081-1947(08)60144-7. [Cross Ref]

Articles from Scientific Reports are provided here courtesy of **Nature Publishing Group**