Quantum-mechanical analysis of effect of alloying elements on ε-martensite start temperature of steels

With regard to the transformation mechanism of austenitic high manganese steel, the prediction of the ε-martensite start temperature is a critical consideration in alloy design. Evaluation of the ε-martensite start temperature makes it possible to predict the microstructure and to understand the phase transformation occurring during deformation. Here we use the quantum mechanical calculation of random alloys to understand the physics for ε-martensitic transformation in steels. We could find the linear relationship between the measured ε-martensite start temperatures and the crystal structure stability for various compositions. We also could estimate the effect of several alloying elements. It is expected that the effect of decreasing the temperatures for the same amount of alloying elements addition will be larger moving farther from Group VIII. By creating a free-energy model that reflects the temperature effect, we were able to calculate the average driving force required for the ε-martensitic transformations.


Results and Discussions
Multiple linear regression model approach. The multiple linear regression model for ε M S temperature was constructed using the least square fit method with 13 input parameters corresponding to the amount of each element added based on 322 data sets. The statistics of the collected data are summarized in Table 1, and the results of the multiple linear regression analysis is as follows,  where w M is the amount of alloying element added, with wt% for M = C, Mn, Ni, Cr, Al, Si, Mo, Co, Cu, Nb, Ti, V, and W. The correlation-coefficient (R) of the model was 0.93, and the squared correlation-coefficient (R 2 ), which indicates the ratio of the variability of the response explained by the model, was 0.87. Figure 1 shows the results, including predicted values and coefficients for each input variables, according to the model. This model can be correlated well with the ε M S temperature depending on the alloying elements, but there are some limitations. (1) It cannot exhibit a non-linear relationship and cannot justify the dependencies between the variables. (2) This model cannot separate the effects of the alloying elements on the ε M S temperature and other effects, such as the initial austenite grain size and errors due to the measurement methods. (3) The effects of other alloying elements (e. g. phosphorus, sulfur, silver and gold) that were not included in model construction cannot be explored.
Thermodynamic driving forces using the CALPHAD approach. The stability of austenite and ε-martensite during cooling is dependent on ΔG γε (T), the Gibbs free energy difference between FCC and HCP structure at a given temperature (T). Here, thermodynamic quantities were evaluated based on the TCFE7.0 database 32 , and they were summarized in Fig. 2. Figure 2a and b show the calculated free energy differences at ε M S temperature and enthalpy differences at 300 K between the FCC and HCP structures, respectively. We calculated the enthalpy difference at 300 K, which is the lower limit of supported range of the database. The results of the free energy difference showed a wide range (−0.5 to +2.0 kJ mol −1 ) for both the part supported by the database and the part unsupported over the composition range. A lot of data show very low free energy changes near zero, and some data even have positive values. This is unreasonable considering the driving force required by the accompanying strain energy during martensitic transformation, thus limiting the prediction of the ε M S temperature 12 . There are 53 data corresponding to Δ > . γε G 0 5 kJ mol −1 , which means that austenite is more stable compared to ε-martensite. Included were corresponding data on 4.6-11.6 wt% Cr and 4.2-7.03 wt% Si, with C, Mn, and Ni. Adding large amounts of Cr and Si result in a large ΔG γε value, so the thermodynamic data corresponding to this concentration is relatively unreliable. When comparing the stability of the two structures through enthalpy difference shown in Fig. 2b, we could not find any relationship between the calculated data and ε M S temperature. It is believed that the thermodynamic database might be less accurate than the experimental error in measuring the ε M S temperature, since the distribution of thermodynamic quantities is not normal.
Quantum mechanical calculation for lattice stability. Using the calculated total energies from quantum-mechanical calculations, we could determine the effect of alloying elements on the lattice stability based on the following equation.
where ε E tot and γ E tot represent the equilibrium total energy per atomic site for the paramagnetic HCP and anti-ferromagnetic FCC structures, respectively. To calculate the total energy of the FCC and HCP structures, the CPA method was applied to reflect the disordering of the alloys. It is important to note that the CPA method has shortcomings that it cannot include the effect of the inhomogeneous distributions and short-range ordering of alloying elements because it uses single site approximation to simulate random alloys. Since the ε-martensitic transformation is a diffusionless transformation from austenite at high temperature, it is believed that the effect of inhomogeneous distribution of alloying elements on the crystal structure stability is not significant.
The characteristics of ε-martensitic transformations of high Mn steels are largely associated to their magnetic properties, and the role of anti-ferromagneitsm has been carefully analyzed in several studies [33][34][35] . The austenite and ε-martensite are stable at high temperatures in a paramagnetic state, but a magnetic transition occurs in an anti-ferromagnetic state as they fall below their respective Néel temperatures (T N (γ) and T N (ε)). Experimental measurements show that ε-martensite is anti-ferromagnetic in some Fe-X systems, but its T N (ε) is approximately 230 K, which is below the region of interest in the ε-martensitic transformations 34,35 . In contrast, the parent austenite is strongly stabilized by anti-ferromagnetic ordering below 500 K in a Fe-Mn binary system 35 . Most of the ε M S temperature data collected in this study are in the range of 230-460 K, and anti-ferromagnetic FCC and paramagnetic HCP are regarded as stable in this temperature range. It is therefore reasonable to compare the ε M S temperature with the crystal structure stability based on these magnetic state. One aspect to keep in mind is that other elements besides Fe and Mn affect T N (γ), at which the magnetic structure of the parent austenite is changed during the ε-martensitic transformation 36 . This is likely to be one of the main causes of the prediction error of the following model.
The ΔH γε (0 K) increased with the addition of Mn from −1.72 kJ mol −1 at 15 at% to 1.48 kJ mol −1 at 40 at% as shown in Fig. 3a. At a composition of more than 28 at%, the lattice stabilities of the FCC and HCP structures reversed, which is similar to the results predicted for 26 at% by the existing CALPHAD calculations 15 . The compositional dependence of ΔH γε (0 K) in the Fe-20Mn-X at% system including Al, Si, Ni, and Cr are shown in Fig. 3b. When the same amount of each element was added, Al had the greatest effect on increasing the austenite stability, and was effective in the order of Ni, Cr, Mn, and Si. Although, the effect of increasing the stability was mostly linear, non-linear relationship between Si content and lattice stability was obtained; the value of Si was highest at around 7 at% and further addition of Si decreased the stability of austenite. This is consistent with the effect of Si on SFE, which shows the highest effect at about 3-4 wt% 37 .
The ΔH γε (0 K) values for alloys against measured ε M S temperature are represented in Fig. 4. Overall, the lower ΔH γε (0 K) was, the higher the stability of HCP, was obtained, which corresponded to inverse correlationship, the higher ε M S temperature. The following relation was obtained based on linear regression models This model has a close relationship with R(0.87) and R 2 (0.77), taking into account the errors in the experimental measurements and the numerical errors. The effect of carbon (C) on the ε M S temperature is well predicted by the assumption that the C-C interaction is insignificant in most csompositions. In addition, the structural stability considering the effect of the interstitial alloying element, C, was not significantly different from that of the pure Fe alloy, or from that after the addition of other substitute alloying elements.
A detailed analysis of the data can provide insight into measurement errors. As can be seen in Fig. 4b, there is a closer linear correlation for the alloy containing C, compared to the C-free alloys. It is believed that the initial austenite grain size difference is large and the ε M S temperature error is also large in C-free steel, because the C-free steel has higher grain growth rate of austenite 38 . Yang and Bhadeshia proposed a model that shows the influence of the parent austenite grain size on the M s temperature and compared it with experimental measurements 39 . The ε M S temperature is also affected by the austenite grain size, similar to α′-martensite, and shows a difference at about 30 K with grain size differences of 10 μm and 200 μm. The relationship derived here reflects the effect of alloying elements alone, but has limitations in that it cannot reflect data scattering regarding the difference in the austenite grain size. In order to reduce the scattering of the ε M S temperature for the same compositions, it is necessary to modify the ε M S temperature reflecting the austenite grain size effect. In the case of the Fe-0.0017C-24.03 Mn wt% alloy, ΔH γε (0 K) was calculated to be −0.505 kJ mol −1 and the ε M S temperature measurement result showed a 54 K difference (421 K 40 and 367 K 41 ) depending on the measurement method. Given the error, it can be seen that there is a fairly precise linear relationship for C containing alloys. Figure 4b shows the Fe-Mn binary system, and the higher the amount of Mn, the lower the obtained ε M S temperature was. The data for approximately the same composition, 25.0 to 25.5 at% Mn, showed a large scatter of measurements (291-390 K), and this is corresponding to the magnetic ordering of austenite. If the parent austenite is in a paramagnetic state assuming completely disordered localized moments, then the magnetic entropy contribution S mag (μ i ) for the magnetic moment μ i of atom i can be estimated S mag (μ i ) = k B log(μ i + 1) with Boltzmann's constant k B 42 . Below T N (γ), this magnetic entropy contribution is reduced by anti-ferromagnetic ordering, making the parent austenite more stable. The T N (γ) of austenite is raised by the addition of Mn, and this has an enhanced effect on suppressing the ε-martensitic transformation if the ε M S temperature is below T N (γ) 34 . Thus, in the ε-martensitic transformation, there is non-linear behavior around T N (γ), where the slope of the ε M S temperature changes depending on the amount of Mn. In addition, overlapping of the magnetic transition and structure transformation makes it difficult to detect the ε M S temperature, which causes large data scattering. This is an essential factor underlying the prediction error.   Figure 5a shows the change in lattice stability ΔH γε (0 K) per 1 at% addition of element X in the Fe-20Mn-X at% system. The ε M S temperature was estimated to be as 397.5 K by Eq. (3) with ΔH γε (0 K) = −1.085 kJ mol −1 for Fe-20Mn at%. The change of ΔH γε (0 K) per 1 at% of element X was investigated for the interstitial element C; for the elements Al, Si, P, and S belonging to 3rd period; for the elements Sc, Ti, V, Cr, Mn, Co, Ni, Cu, and Zn belonging to 4th period; for the elements Y, Zr, Nb, Mo, Tc, Ru, Rh, Pd, Ag, and Cd belonging to 5th period; and for the elements W, Re, Os, Ir, Pt, Au, and Hg belonging to 6th period. The change in ε M S temperature as determined by Eq. (3) is plotted on the right axis. Figure 5b shows the effect of each element predicted through multiple linear regression, which is summarized in Fig. 1b, according to the atomic number with a factor corresponding to at% instead of wt%. The same ε M S temperature change scale is shown for comparison with Fig. 5a.
The In the case of Si addition, the ΔH γε (0 K) was estimated to be increased by 0.090 kJ mol −1 and have a correspondingly lower ε M S temperature by 6.6 K, while the coefficient of the multiple linear regression was negative. As shown above, the multiple linear regression model does not reflect non-linearity, while the Si effect increases ΔH γε (0 K) from 0 to 7 at%, and decreases it in the region above 7 at%. It is expected that the addition of Ti, V, Nb, Mo, and W at 1 at% will increase ΔH γε (0 K) to of 0.501, 0.253, 0.799, 0.558, and 0.628 kJ mol −1 , respectively, and lower the ε M S temperature to 36.9, 18.6, 58.8, 41.0, and 46.2 K, respectively. This is predicted to be small overall compared with the coefficients of multiple linear regression models 74.9, 31.7, 87.3, 34.3, and 43.8 K, but it shows qualitatively the same trend as ΔH γε (0 K). The large error in the Ti, V, Nb, Mo, and W coefficients is probably due to insufficient number of data used to predict the coefficients. In addition, the influence due to other independent elements such as Mn is not reflected.
We can also deduce the effect of other elements that have not been subjected the actual experiment. For example, addition of 1 at% P will lower ΔH γε (0 K), which will increase the ε M S temperature. As all the elements used in this search except P were added, its addition is expected to increase ΔH γε (0 K) and thereby lower the ε M S temperature. These effects show a certain regularity according to the atomic number. Elements in the Group VIII (Fe, Ru, Os) showed the lowest effect on the lattice stability, and the wider the distance from Group VIII, the greater the effect on the lattice stability per 1 at%. This shows that the FCC and HCP lattice stability tended to be similar to that predicted by Pettifor's theory 43 , which is determined by the electron occupancy of the d-orbital. Therefore, the effect on the ε M S temperature is mainly dominated by the occupancy of the d-orbitals. In order to confirm the effect of alloying elements in the proposed model, two samples having compositions with different P content were prepared and ε M S temperature was measured. A sample was prepared by adding a certain amount of P based on Fe-Mn binary. P is an alloying element whose influence on ε M S temperature has not been clarified by previous experiments. The chemical compositions of the prepared samples were Fe-20.8Mn at% and Fe-20.7Mn-0.27P at%, respectively. For the Fe-20.8Mn at% sample, the measured ε M S temperature, 387 K, is about 3 K lower than the value of 390 K predicted using Eq. (3). The measured ε M S temperature for the Fe-20.7Mn-0.27P at% sample is 391 K, which is consistent with the result that the ε-martensitic transformation is enhanced by adding P. Also, the effect of ε M S temperature increment of 8 K per 1 at% is expected, which corresponds to 2.2 K per 0.27 at%. The predicted effect of P appears to be little underestimated compared to a 4 K temperature increase with experimental measurements. This is thought to be caused by the difference in Mn content. Experimental verification confirmed the reliability of the ε M S temperature prediction model based on the quantum-mechanical calculations. Driving force evaluation for ε-martensitic transformation. The original formulation of density functional theory was for the T = 0 K ground state 44 . Although it can be shown that the concept can be extended to finite electronic temperatures, historically, most calculations were restricted to zero temperature properties due to computational costs 45 . The molar Gibbs free energy of the crystal at a given temperature T and pressure P can be approximately decomposed into several contributing terms [45][46][47][48] .
The leading term E tot is the total energy which is directly obtained from the electronic structure calculation allowed to be independent of temperatures. The term F vib , F conf , F el , and F mag account for the vibrational, configurational, electronic, and magnetic contributions, respectively. The last one is the pressure-volume dependent term. The free energy change ΔG γε (T) accompanying the ε-martensitic transformation involves ΔH γε (0 K), and can be described as follows.
where ΔH ext includes all temperature independent contributions to free energy except for electronic total energy. This includes differences due to atomic volume differences, magnetic transitions, and errors due to the calculation approximations and parameters commented above. Here, ΔS is an entropy term including electronic, magnetic, vibrational, and configuration excitations. A linear relationship between ε M S temperature and ΔH γε (0 K) appears in the free energy configuration of Eq. (6). Assuming ΔH ext is independent of T, the free energy change at ε M S temperature is expressed as follows, The driving force of the ε-martensitic transformation, Δ γ ε ε G M ( ) S , varies slightly depending on the alloy composition, but the difference is less than 0.1 kJ mol −112, 49,50 . The linear relationship between ΔH γε (0 K) and ε M S temperature shows that ΔH ext and ΔS in Eq. (7) can be approximated as a constant due to the small degree of variation depending on the alloy composition. In the case of the ΔS value, it can be approximated as −0.0136 kJ mol −1 K from the comparison of Eqs (4) and (7). The ΔH ext value can be determined from the T 0 temperature, which corresponds to ΔG γε (T 0 ) = 0 in the Fe-Mn binary system. For the Fe-28.45Mn at% system, ΔH γε (0 K) = −0.025 kJ mol −1 was calculated and ΔH ext = −4.671 kJ mol −1 was determined from the experimentally measured result with T 0 = 345.36 K in the system 49 .
The T 0 temperatures and free energy changes during martensitic transformation in the Fe-Mn binary system are evaluated based on Eq. (7). The prediction of T 0 temperatures is well correlated with the existing measurement data as shown in Fig. 6a. As mentioned above, the Fe-Mn system exhibits a non-linear tendency in which the slope of ε M S temperature as a function of the amount of Mn changes near T N (γ). The T 0 temperature is mainly determined by the ε M S temperature, and there is also a non-linear tendency due to the magnetic ordering. Below 18 at% Mn, the mangetic entropy contribution increases, and the free energy model derived from linear relationship. Nevertheless, the error is within the range of data scattering by measurement within 30 K.
The free energy changes during martensitic transformation, Δ γ ε ε G M ( ) S , show nearly normal distribution in Fig. 6b with mean of −0.281 kJ mol −1 and standard deviation of 0.39 kJ mol −1 , unlike in thermodynamic calculations. This standard deviation corresponds to ±28.7 K at ε M S temperature. It is considered that the distribution of Δ γ ε ε G M ( ) S is reasonable when the measurement error considered. The average driving force, 281 J mol −1 , is within the existing results of 68-120 J/mol 49 , 200-300 J/mol 50 , and 500-900 J/mol 51 . Figure 6b also shows the free energy changes at ε M S temperature in the ternary Fe-17Mn-Al, Fe-17Mn-Cr, Fe-17Mn-Mo, and Fe-17Mn-W in the unit of at%. The Fe-Mn binary and Fe-17Mn-Al systems show underestimated ε M S temperature with a lower Δ γ ε ε G M ( ) S value than the average driving force. In contrast, the values of the Fe-17Mn-Mo and Fe-17Mn-W systems are located above the average driving force line and show a higher measured ε M S temperature than the predicted value. However, the overall free energy change is close to the average driving force and the ε M S is within the prediction error. This method allows a temperature-dependent free energy model to be established based on the quantum-mechanical calculations.

Conclusions
On the basis of the quantum-mechanical calculations, we discuss the relationship between the stability of FCC and HCP lattice structures and the ε M S temperature of various alloy systems. We first note that there is a close linear relationship between the measured ε M S temperature and the crystal structure stability of the anti-ferromagnetic FCC as well as the paramagnetic HCP at the 0 K. Thus, the effect of each element on the ε M S temperature could be determined. The effect of each element on ε M S temperature is related to the occupancy of electrons filling in the d-orbital. As the distance between the element and the Group VIII widens, it was found that the degree of lowering the ε M S temperature by the addition of the same amount of element, was increased. Based on the measured T 0 and ε M S temperature, the force driving the transformation to ε-martensite was calculated to be −0.281 kJ mol −1 on average, which is within the predicted range.
Our study is the first to show the relationships between the experimentally measured data of various compositions of steels, and quantum-mechanical calculations. Applying the alloy composition and temperature dependency to the results based on the density-functional theory was difficult problem. The coherent-potential approximation method was used to simulate the effect of alloys; so that the effect of alloying elements on the ε M S temperature could be deduced. It is believed that the quantum-mechanical calculations can be used to predict the Scientific RepoRTS | (2017) 7:17860 | DOI:10.1038/s41598-017-18230-z characteristics of steel alloys, which are typical structural materials, and can now be used to design new alloys by predicting their martensitic transformation start temperature.

Data.
A search of the literature revealed a large number of data sources some 322 combinations on chemical composition and ε M S temperatures 52 . Included were data on Mn corresponding to 13-35.9 wt%, and on twelve other elements (C, Ni, Cr, Al, Si, Mo, Co, Cu, Nb, Ti, V and W). Of the total data, 232 do not contain carbon. In addition, 68 data are distributed from 0 to 0.08 wt% or less carbon content, 10 data from 0.09 to 0.16 wt% or less, and 12 data from 0.16 to 0.35 wt%. For Mn, there were 59 data for Fe-Mn binary alloys, with Mn proportion ranging from 13.6 to 29.3 wt%. In the case of Mo, Ti, Nb, V, and W, there were 4, 4, 3, 3, and 3 data, respectively, with a specified amount added; so it is necessary to consider that there is an insufficient number of data to effectively grasp the tendency. In contrast, in cases of C, Mn, Ni, Cr, Al, Si, Co and Cu, a relatively wide variety of alloy systems are included, thereby ensuring reliability in understanding the tendency.
DFT calculations. The lattice stability of random alloys was calculated using the EMTO-CPA method 53 . In the present calculations, the one-electron equations were solved using the full charge density (FCD) 31 method and frozen-core approximation, i.e., the core states were fixed to the initial atomic states. All self-consistent calculations were performed using the local density approximation (LDA) of the effective exchange-correlation potential, and the total energies were obtained using the Perdew-Burke-Ernzerhof (PBE) realization of the generalized gradient approximation (GGA) 54 . Poisson's equation was solved within the spherical cell approximation(SCA) 53 . The Green's function for the valence states was calculated for 16 complex energy points distributed exponentially on a semi-circular contour. The EMTO basis set included s, p, d, and f orbitals. We used more than 2000 inequivalent k-points in the irreducible wedge of the FCC, HCP, and BCT Brillouin zone, which ensured relative accuracy 1.0 × 10 −3 eV/atom in the total energy differences. Self-consistency was assumed when the distances between the input and output total energies and fermi level were less than 1.3 × 10 −7 eV and 1.3 × 10 −6 eV, respectively. The convergence of these computational parameters was carefully checked.
Austenite with a FCC structure was simulated using a primitive unit cell containing one lattice site and a BCT structure including two lattice sites. We include different magnetic structures for austenite, including nonmagnetic (NM), paramagnetic (PM), and anti-ferromagnetic (AFM) state. For austenite, the AFM state is known to be stable at low temperature 55,56 . This was realized using a BCT unit cell for the structure in which the spin up-down state is crossed at every (001)-plane. According to the results of previous studies, the structure in which the spin direction changes every two layers was calculated to be the most stable, but the energy difference was found to be insignificant 57,58 . The ε-martensite, which has the HCP structure, was simulated assuming a para-magnetic state based on a primitive unit cell containing two lattice sites. The paramagnetic configuration of the FCC and HCP alloys was simulated by means of the disordered local moments (DLM) model 59 . The DLM model was shown to describe the random distribution of the local magnetic moments of the paramagnetic state of metals, above the magnetic transition temperature.
The crystal structures of different alloys and their magnetic states were optimized. For FCC, the energy for nine points was calculated in 0.5% increments from −2% to 2% based on the lattice parameter of a γ = 3.6 Å, which is the experimental value of austenite. In the case of ε-martensite with HCP structure, seven c/a ratios  580, 1.585, 1.590, 1.595, 1.600, 1.605, and = . 8/3 1 632, the ideal ratio for FCC structure) were included. The energy for nine points was calculated in 0.5% increments from −2% to 2% based on = = .
γ a a / 2 2 74 0 Å which corresponds to the austenite lattice parameter. The optimum lattice parameter, atomic volume, and c/a ratio were determined based on calculations for different 63 combinations of a 0 and c 0 . Equilibrium volume, equilibrium total energy, and bulk modulus were calculated using a third-order Birch-Murnaghan fit after calculating energy values per atomic volume in each structure 60 . The optimum c/a ratio was determined by fitting the equilibrium total energy function to the c/a ratio with a fourth order polynomial. The element concentrations of each alloy were converted to mole fractions and rounded to an even number at the third decimal place to achieve the same amount of alloying in NM, AFM, and DLM simulations. The lattice stability effect of carbon as an interstitial element, was applied using the results calculated based on the FCD-EMTO method in a previous study 14 . Based on the value of ΔH γε (0 K) = 3.66 kJ/mol calculated for the 32Fe-1C structure, which corresponds to 0.667 wt% C, the result was 5.49 kJ/mol per 1 wt%. The carbon range of interest is less than 0.35 wt%, and the dilute solution linearly approximated the effect of C, assuming that the effect of C-C interaction is small. We neglected all thermal contributions in the present calculations. The datasets generated during the current study are available in Supplementary Table 1 and the ref. 61 .

Measurements of M S
ε temperature. A set of alloys was prepared as 400 g melts, which were cast and vacuum sealed. They were then homogenized for two days at 1473 K, after which their chemical compositions were measured. Cylindrical samples with 3 mm diameter and 10 mm length were machined and studied using a dilatometer. The samples were heated at 2 K s −1 under a vacuum to 873 K for 1 min, and then quenched at 10 K s −1 using argon gas to room temperature. The offset method was used to determine the ε M S temperature, and the critical strain value was set at 2% of the strain amount at the completion of the transformation. More than two experiments were conducted on the same composition, and the measured values were distributed within 1 K.