First Principles Investigation on Thermodynamic Properties and Stacking Fault Energy of Paramagnetic Nickel at High Temperatures

: Reliable data on the temperature dependence of thermodynamic properties of alloy phases are very useful for modeling the behavior of high ‐ temperature materials such as nickel ‐ based superalloys. Moreover, for predicting the mechanical properties of such alloys, additional information on the energy of lattice defects (e.g., stacking faults) at high temperatures is highly desirable, but difficult to obtain experimentally. In this study, we use first ‐ principles calculations, in conjunction with a quasi ‐ harmonic Debye model, to evaluate the Helmholtz free energy of paramagnetic nickel as a function of temperature and volume, taking into account the electronic, magnetic, and vibrational contributions. The thermodynamic properties of Ni, such as the equilibrium lattice parameter and elastic moduli, are derived from the free energy in the temperature range from 800 to 1600 K and compared with available experimental data. The derived temperature dependence of the lattice parameter is then used for calculating the energies of intrinsic and extrinsic stacking faults in paramagnetic Ni. The stacking fault energies have been evaluated according to three different methodologies, the axial ‐ next ‐ nearest ‐ neighbor Ising (ANNNI) model, the tilted supercell approach, and the slab supercell approach. The results show that the elastic moduli and stacking fault energies of Ni decrease with increasing temperature. This “softening” effect of temperature on the mechanical properties of nickel is mainly due to thermal expansion, and partly due to magnetic free energy contribution. of chemically and magnetically disordered metallic systems.


Introduction
Nickel-based superalloys are an important group of high-temperature materials. They are widely used for components servicing at high temperatures, such as the blades of power station turbines and aircraft engines, for their excellent corrosion resistance, oxidation resistance, hightemperature strength, and creep strength [1][2][3]. Since the relationships between the microstructure and high-temperature properties are complex, research into nickel-based superalloys is still a challenge that attracts many scientists in industry and academia. Understanding the physical properties of pure nickel at elevated temperatures is of pivotal importance for further research progress because nickel is the main component of the matrix in nickel-based superalloys.
Temperature-dependent elastic properties, such as the shear modulus and lattice constant, are important parameters in modeling the high-temperature strength (in particular, creep strength) of alloys; for example, for evaluating solid solution strengthening (SSS) [4][5][6][7]. Despite the importance of such parameters, few experimental data can be found. During the past decades, the first-principles method based on density functional theory (DFT) was widely used to investigate the elastic properties of metal alloys [7][8][9][10]. Although the elastic moduli of nickel alloys were calculated in some studies, most of them reported results for 0 K [7,11,12]. Shang et al. [13][14][15] calculated elastic moduli at high temperatures by combining first-principles calculations and the Debye model. However, in most of their calculations, nickel was treated in the magnetically ordered ferromagnetic state even at high temperatures. Actually, at the temperatures above the Curie point, nickel is paramagnetic. In many previous studies, the paramagnetic state of Ni was described as non-magnetic (i.e., using spinunpolarized DFT calculations), and the calculations of elastic constants at high temperatures were often performed assuming the lattice constant of Ni at room temperature. Clearly, all these inconsistencies (ignoring the effects of magnetic disorder and thermal expansion) may introduce uncontrollable errors in the description of the paramagnetic state of nickel and its alloys.
Thermodynamic properties, such as entropy and heat capacity, are useful in the research and development of metallic materials. These properties are also of direct relevance to Calphad modeling, which is in the core of the Materials Genome. Due to the difficulties of experimental measurements at high temperatures, such properties are not easily available for some materials. In such cases, firstprinciples-based simulations are becoming an indispensable tool for predicting the thermodynamic properties at finite temperatures. The work by Shang et al. [13] is one of the first attempts to predict the entropy, enthalpy, and heat capacity of nickel based on first-principles calculations. However, nickel has been treated in the ordered ferromagnetic state at a high temperature as well.
First-principles calculations of stacking fault energy (SFE) is another domain of materials research where computational modeling is becoming important as a source of data. One reason is the lack of experimental data due to difficulties of measurements, especially at elevated temperatures. At the same time, stacking fault energy is one of the factors determining the mechanism of plastic deformation, and it, therefore, plays an important role in modeling the plastic behavior of metals. Stacking fault energy has been regarded as a critical parameter in some studies of creep properties [4,16,17]. As has been mentioned above, most of the computational studies of nickel correspond to T = 0 K. Although Zhao et al. [18] considered the effect of temperature, nickel was treated as ferromagnetic and non-magnetic at high temperatures.
In present work, a quasi-harmonic Debye model has been incorporated into the framework of first-principles calculations to describe the vibrational contribution to the Helmholtz free energy. Nickel has been considered at a finite electronic temperature (above the Curie point) to include thermal electronic excitations, and in the paramagnetic state represented by the disordered local moment (DLM) model to describe the magnetic moment fluctuations at the same temperature. The temperature-dependent lattice parameter, thermal expansion coefficient, elastic properties and other thermodynamic properties were derived from the Helmholtz free energy, calculated as a function of temperature and volume. The theoretical results had been compared with experimental data. The energies of intrinsic and extrinsic stacking faults were also calculated as a function of temperature. The effects of thermal expansion and magnetic fluctuations on the stacking fault energy are singled out and discussed.

Helmholtz Free Energy and Lattice Parameter
All first-principles calculations in this work were performed by the exact muffin-tin orbital (EMTO) method [19][20][21] combined with the coherent potential approximation (CPA) [22,23]. In the self-consistent of all calculations, the exchange-correlation energy was taken in the form suggested by Perdew, Burke, and Ernzerhof (PBE) for the generalized gradient approximation (GGA) [24]. A 33 × 33 × 33 Monkhorst-Pack [25] grid of k-points was used in the integration over the Brillouin zone; the full-charge density (FCD) technique [26,27] was employed in the total energy calculations. The magnetic contribution of paramagnetic nickel was evaluated using the self-consistent disordered local moment (DLM) approach [28] which has been implemented into the EMTO calculation scheme.
The Helmholtz free energy of a system with volume V and at temperature T can be expressed as: where V is the volume of primitive cell, is the electronic free energy (including the magnetic fluctuation part) and is the contribution that comes from lattice vibrations. The equilibrium lattice volume (and, therefore, the lattice parameter) is obtained by minimizing the Helmholtz free energy represented as a function of volume by a four-parameter Murnaghan equation of state. More details will be given in Section 2.3.

Free Energy Contributions Due to Electrons, Magnetism and Lattice Vibrations
In the electronic free energy, the Fermi-Dirac distribution function was used to model the excitation of valence electrons from the energy states below to the energy states above the Fermi level [29,30]: where ϵ is the energy of an electron state, is the Fermi level, is the Boltzmann constant and T the electronic temperature. Then, the contributions of electronic energy and entropy were evaluated by a contour integration in the complex energy plane. A piecewise-elliptical contour was used for Green's function integration, with a total number of 56 complex energy points. The bottom of the contour was chosen some 1Ry below the Fermi level to envelope all populated valence states and four Matsubara poles [29].
Above the Curie temperature (about 631 K for nickel [31]), the ferromagnetic order is destroyed and Ni transforms to a paramagnetic state. In the paramagnetic state, fluctuating local magnetic moments are assumed to exist and to be distributed randomly when treated in the self-consistent DLM model [28]. The corresponding magnetic entropy can be estimated as follows, where ms is the calculated local spin moment on a nickel atom and is the Boltzmann constant. Together with the terms due to thermal excitations of electrons across the Fermi level, the -TSmag term and the corresponding total energy contribution were included into the electronic part of the Helmholtz free energy, Equation (1).
In previous studies, several versions of the Debye model were used to model contribution of phonons. The Debye-Grüneisen model is the most popular one which employs the so-called Grüneisen parameter to describe the volume dependence of the Debye temperature, as well as some empirical scaling factors [13,32,33]. In this work, the vibrational free energy was estimated in the quasi-harmonic Debye model [30,34] as follows, where kB is the Boltzmann constant and Θ is the Debye temperature. The Debye temperature Θ was evaluated at every volume V and temperature T from the mean acoustic velocity as suggested by Anderson [35], where , and are mean acoustic velocity, acoustic velocities for longitudinal and transversal sound waves, respectively. The longitudinal and transversal velocities can be determined from volume-and temperature-dependent bulk modulus and shear modulus. ,

Elastic Constants
The flow chart of finite-temperature free energy modelling used in the present work is presented in Figure 1. In order to take into account all contributions due to electronic, magnetic, and vibrational degrees of freedom, three steps were required. In the first step of the calculations, contribution of lattice vibration was not included, the EMTOcomputed partial free energy Fel(V,T) only contained the electronic and magnetic contributions. A four-parameter Murnaghan equation of state [36], in the form of Equation (10), was fitted to the calculated partial free energies at temperature T as a function of primitive cell volume V. Thereby, an "adiabatic" bulk modulus (that is, not containing the phonon contribution) was obtained from the equation of state: where is the volume corresponding to the minimum of partial free energy at temperature T, is the bulk modulus determined at temperature T and volume , and is the pressure derivative of the bulk modulus.
For a cubic structure, there are, three independent single-crystal elastic constants in total, i.e., C11 C12 and C44. At every volume V and temperature T, the three elastic constants can be derived from the bulk modulus B, Equations (10) and (11), and two elastic shear constants corresponding to two volume-conserving distortions. In the elastic constant calculations for cubic systems, orthorhombic and monoclinic distortions are commonly applied to the primitive cell, given, respectively, by Equations (12) and (13), For each considered volume and temperature, the and C44 elastic shear constants were derived from the strain-dependent energy changes as follows [8]: Five strain values were used in the calculations for each type of distortion, i.e., , 0.00, 0.01, 0.02, 0.03, 0.04.
In the second step, indicated in Figure 1 using boxes with dashed boundaries, the results of the adiabatic volume-dependent bulk modulus, volume-dependent C11, C12, and C44 were used for calculating the temperature-and volume-dependent average shear modulus that appears in Equations (8) and (9) through the Voigt-Reuss-Hill averaging method [37] as follows, where and are the shear modulus estimated in the Voigt [38] and Reuss [39] limits, respectively, and is the Hill average of the shear modulus. Then the adiabatic volume-dependent bulk modulus calculated according to Equation (11) and the shear modulus calculated according to Equation (20) were plugged into Equations (8) and (9) to get the averaged acoustic wave velocities and the Debye temperature, Equation (6).
In the third step, the temperature-and volume-dependent Debye temperature was used as the parameter to calculate the vibrational free energy in the quasi-harmonic Debye model according to Equations (4) and (5). After the total Helmholtz free energy with the contribution of phonons was obtained, its volume dependence was fitted with the four-parameter Murnaghan equation of state to obtain the equilibrium volume V0 and isothermal bulk modulus B0 corresponding to temperature T.

Stacking Fault Energy
In this work, three approaches were used to calculate the intrinsic and extrinsic stacking fault energy. In the first approach, the stacking fault energy was calculated within the axial-next-nearestneighbor Ising (ANNNI) model which is a widely used method in the stacking fault energy calculation of semiconductors, metallic materials, etc. [40][41][42][43]. In this model, the energies of intrinsic stacking fault (ISF) or extrinsic stacking fault (ESF) are approximated by certain combinations of the Helmholtz free energies of different short-period stacking sequences of close-packed (111) planes. The ISF and ESF energies can be estimated within the first-and second-order ANNNI models. In the first-order ANNNI model, the energies of extrinsic and intrinsic stacking faults are equal, i.e., The second-order ANNNI formulae read as follows: where , and are the Helmholtz free energies (per atom) of the perfect hexagonalclose-packed (HCP) structure, double hexagonal-close-packed (DHCP) structure and face-centeredcubic (FCC) structure, respectively. The stacking sequence along the [111] direction is ABAB for HCP structure, ABAC for DHCP structure and ABC for FCC structure. The equilibrium volume of FCC structure at a given temperature T was used to calculate all free energies proposed above. For HCP and DHCP structures, the c/a ratio was set to the ideal value of approximately 1.633. The interface area can be derived from the equilibrium lattice constant as follows, where is the equilibrium lattice constant of FCC crystal structure at a given temperature. The second approach, i.e., the supercell method, is the most popular method used in SFE calculations. In the present work, the supercell containing nine (111)   The third approach, i.e., the slab supercell method, is the most comprehensive method to calculate SFE. As Figure 3 shows, The intrinsic and extrinsic stacking faults can be generated by, respectively, removing from or inserting into the supercell one layer of atoms. Although the method of Anderson [35] allows one, in principle, to derive the Debye temperature from the elastic constants for any structure, in practice, the Debye model cannot accurately represent the differences in harmonic vibrational free energy between different close-packed structures. For large supercells, the numerical errors of elastic constant calculations (translated into Debye temperature uncertainty) may overshadow the sought-after free energy differences with the FCC structure. In this work, to describe all close-packed structures on an equal footing, we approximate their harmonic vibrational free energies (per atom) by those calculated for the FCC structure at the same temperature and atomic volume. Within this approximation, the vibrational part of stacking fault energy is reduced to the anharmonic contribution due to thermal expansion. The harmonic contribution to SFE of Ni has been studied in detail in [42]. This work is focused on the anharmonic and magnetic contributions to the stacking fault energy.

Lattice Misfit and Thermal Expansion Coefficient
The lattice parameter and thermal expansion coefficient (TEC), calculated as a function of temperature, are shown in Figure 4 and Figure 5 where the calculated results are compared with available experimental information. It can be seen in Figure 4 that the theoretical values of the lattice parameter, as well as its temperature dependence, are similar to the experimental data. Although the theoretical values of TEC are smaller than experimental values at about 800 K, the theoretical and experimental TEC values are close to each other at temperatures above 1000 K. The temperature dependencies of the lattice parameter and TEC agree with the experimentally observed trends at high temperatures.
The good agreement between theoretical results and experimental observations indicates that, with contributions due to magnetic fluctuations, electronic excitations, and lattice vibrations, the EMTO-CPA based calculations can predict the equilibrium atomic volume of paramagnetic nickel at high temperatures reasonably well. The accuracy of the quasi-harmonic Debye model turns out to be sufficient for representing the free energy contribution of lattice vibrations (phonons). In combination with the EMTO-CPA methodology, it provides an efficient technique for evaluating the free energy (and related thermodynamic properties) of chemically and magnetically disordered metallic systems.

Temperature-Dependent Elastic Properties
Summing together the EMTO-calculated partial free energy , and the vibrational free energy , given by the quasi-harmonic Debye model, one obtains the total Helmholtz free energy, Equation (1), from which the equilibrium volume and isothermal bulk modulus can be determined as specified at the end of Section 2.3 (the third step). At the equilibrium volume, the shear elastic constants C' and C44 can be evaluated (either by interpolating the previously calculated values or by making additional EMTO calculations at the equilibrium volume and temperature T). Thus, one obtains the set of temperature-dependent single-crystal elastic constants C11, C12, and C44, from which the isothermal shear modulus G0 of paramagnetic nickel can be determined as the Hill average. With these results at hand, the polycrystalline Young's modulus Y can be calculated as follows: The theoretical results are presented and compared with experimental data in Figure 6 and Figure 7. The calculated values and temperature dependence of the isothermal bulk modulus are similar to those of the experimental data. Essentially, the same conclusion holds for the values of single crystal elastic constants C11 and C12; however, the slopes of theoretical and experimental temperature dependencies of these elastic constants are slightly different. Moreover, EMTO calculations seem to substantially overestimate the value of C44, which also explains why the calculated Young's modulus and shear modulus deviate from the corresponding experimental data. Except for these problems, the finite-temperature EMTO calculations, combined with quasi-harmonic Debye model, correctly convey the tendency of elastic moduli to decrease with increasing temperature.

Thermodynamic and Magnetic Properties
In Figure 8, the calculated magnitudes of the paramagnetic spin moment on a nickel atom in FCC and HCP crystal structures are plotted as a function of the calculated equilibrium lattice parameter for the FCC structure (HCP lattice spacing is taken to be the same) at temperatures in the range shown in Figure 4. Therefore, the monotonous increase in magnetic moments presented by Figure 8 is not solely an effect of volume; the magnitude of magnetic fluctuations tends to increase also as a function of temperature, as it yields a higher magnetic entropy according to Equation (3). This magnetic entropy contribution is plotted as a function of temperature in Figure 9 (red squares). According to the comparison in Figure 8, the paramagnetic moment of HCP Ni is obviously larger than that of FCC Ni. This information is in contrast with the assumptions in Calphad modeling, where the magnetic moment of HCP nickel is assumed to be equal to that of FCC Ni [53] because experimental information is lacking. In a previous ab inito study, [42], the contribution to SFE due to magnetic disorder was neglected because the magnetic moments of HCP and FCC Ni were calculated to be quite similar. However, in this work, we show that the magnetic contribution is quite sizeable and, together with the vibrational and electronic contribution, makes the SFE decrease as a function of temperature.
In the present work, the electronic and magnetic entropy contributions have been obtained selfconsistently from the EMTO calculations, while the vibrational entropy was estimated in the Debye model [29] as 4 3 ln 1 (26) where (V,T) is the Debye temperature, is the Boltzmann constant. In the present work, the total entropy . Knowing the thermal expansion coefficient, bulk modulus and entropy, the isobaric heat capacity can be estimated as follows [13], where is the isobaric heat capacity, is the isochoric heat capacity, is the linear thermal expansion coefficient, B(T) is the isothermal bulk modulus calculated at the equilibrium lattice parameter corresponding to temperature T, and S(T) is the entropy.
In order to assess the relative importance of the electronic, magnetic, and phonon contributions to the total entropy and isobaric heat capacity, the corresponding parts and the total values have been plotted in Figure 9 and Figure 10, where they are also compared with experimental values for paramagnetic nickel. Although the values of both total entropy and heat capacity are slightly smaller than the experimental values, the temperature dependences are very similar.
The largest contribution to both entropy and heat capacity comes from the lattice vibrations. Compared to the influence of phonons, the contributions due to electronic excitations and magnetic fluctuations are smaller. This means that a good model description of the vibrational free energy contribution is a key point in the modeling of temperature-dependent properties. In the present work, we have demonstrated that the EMTO-CPA method calculations, including the Fermi-Dirac distribution and self-consistent DLM formalism and combined with a quasi-harmonic Debye model, can reasonably well describe the contributions from electronic excitations, magnetic disorder, and lattice vibrations. This makes the EMTO-based approach a powerful and efficient tool for predicting the thermodynamic properties of materials.

Stacking Fault Energy
Solid lines in Figure 11 and Figure 12 show the free energies of formation for an intrinsic stacking fault (ISF) and an extrinsic stacking fault (ESF) that have been calculated as a function of temperature at the theoretical equilibrium volumes shown in Figure 4. The calculations were done using the different approaches described in Section 2.4. The difference between the results of tilted supercell and slab supercell approaches is found to be negligible, so their dependencies in Figure 11 and Figure  12 are practically on top of each other. In some previous studies, a vacuum layer was added to the slab supercell geometry to yield stress-free boundary conditions [42,44,55]. However, in the present work, the slab supercell model can reasonably well predict the stacking fault energies even without the vacuum layer. At 1000 K, the ISF energy is calculated to be 113 mJ/m 3 which is close to the value 118 mJ/m 3 in the work of Zhao et al. [18], where the paramagnetic state of nickel had been treated as a nonmagnetic state; the smallness of the magnetic contribution should be the reason for the small difference between the ISF energy value obtained in the present work and Zhao's result.
In the ANNNI model, the ISF and ESF energies can be calculated by the first-order and secondorder formulae [40,42,43]. Often the same expression for the ESF energy is used in the first and second order ANNNI models [56], but in this work we use a more accurate second-order expression in Equation (23). ANNNI models yield lower values of stacking fault energy compared to the results of supercell calculations, both for ISF and ESF. Additionally, as can be observed in Figure 11, compared to the first-order ANNNI model, the ISF energies calculated by the second-order ANNNI model are closer to the results of presumably more rigorous supercell calculations, which indicates a rapid decay of the effective inter-layer interactions in the ANNNI model for paramagnetic Ni. Therefore, the second-order ANNNI model is a computationally efficient way to estimate the ISF energy with a sufficient accuracy.
In order to separate the effect of thermal expansion from that of electronic/magnetic excitations on the temperature dependence of stacking fault energy, in Figure 11 and Figure 12 we also plot (with dashed lines) the results of stacking fault energy calculations performed at a fixed atomic volume corresponding to the calculated equilibrium lattice parameter for paramagnetic Ni at 1000 K. In these fixed-volume calculations, the effect of thermal expansion had been excluded. One can see from the difference in slopes between the dashed and solid lines in Figure 11 and Figure 12 that the contribution of thermal expansion into the temperature dependence of stacking fault energy is quite substantial because, according to Equations (26) and (27), the stacking fault energy is inversely proportional to the square of the equilibrium lattice parameter. Since the thermal expansion is mainly a result of lattice vibrations, the dominant contribution to the decrease in SFE with temperature in paramagnetic Ni comes from phonons. This finding coincides with the conclusion made in [57]. In both studies, the harmonic phonon contribution (due to the differences of Ni phonon spectra in different crystal structures, such as FCC and HCP) is neglected because it cannot be fully captured by the Debye model which assumes a model phonon spectrum. For close-packed lattices, that contribution is often assumed to be small, although recent ab initio studies show that it is quite substantial for Cu [58] and Ni [42] metals.  The ISF and ESF energies of paramagnetic Ni decrease with temperature even if evaluated at a fixed volume. This phenomenon indicates that effect of electronic and magnetic fluctuations also results in the lowering of the SFE. The ANNNI model allows us to trace this effect back to the difference in entropy (electronic and magnetic) of Ni in the FCC and HCP structure. Indeed, Figure  8 shows that the paramagnetic moment of HCP Ni is higher than that of FCC Ni, and this translates into the higher magnetic entropy of Ni atoms at the stacking fault, where they effectively have an HCP local environment. The electronic contribution is roughly proportional to the density of states near the Fermi level [57]. Similar to magnetic entropy, because of the higher density of states at the Fermi level, the electronic entropy of paramagnetic HCP Ni is higher than that of FCC Ni, as Figure  13 shows. This entropy effect lowers the free energy of stacking faults in paramagnetic Ni. The decreasing tendency of SFE found in the present work agrees with the DFT result obtained in Zhang's work for Ni [42]. It worth mentioning that the effect of magnetic disorder on the stacking fault energy in Ni is opposite to that in -Fe and austenitic stainless steel [59]. As discussed above, the higher magnetic entropy of HCP Ni contributes to the decrease in SFEs. On the contrary, in the paramagnetic systems based on FCC Fe, the magnetic entropy contribution makes the stacking fault energy increase with increasing temperature, due to the fact that the magnitudes of paramagnetic moments of FCC Fe (and eventually those of Mn and Co solutes in austenitic steel) are larger than the corresponding values for the HCP structure. This phenomenon has been recently reported for the stacking fault energy of an Fe-Mn alloy [60].
In Calphad modeling, the magnetic moment of HCP Ni was treated empirically or by extrapolations from other systems, due to the absence of direct experimental data. For instance, in [61], the magnetic entropy of HCP Ni was said to be equal to that of FCC Ni. The incorrect treatment of electronic and magnetic entropy contributions for HCP Ni may result in the wrong temperature dependence of stacking fault energy being obtained [42]. This has been indirectly confirmed in a study of a NiCo system in [57], where the presence of Co caused a small difference between the densities of states for the HCP and FCC crystal structures which increased the structural free energy difference and, thereby, the SFE. Therefore, the influence of HCP Ni cannot be neglected and should be treated appropriately.

Conclusions
The temperature-dependent properties, including the lattice parameter, thermal expansion coefficient, elastic moduli, and stacking fault energies of paramagnetic nickel have been derived from a quasi-harmonic free energy model parameterized using first-principles calculations at finite electronic temperature. The quasi-harmonic Debye model had been used to include the effect of lattice vibration and phonons. The contributions of electronic excitations, magnetic disorder, and lattice vibrations have been separated and compared. A good agreement is observed between theoretical values of thermal properties and experimental data for paramagnetic nickel in the temperature range from 800 to 1600 K. In particular, the total entropy and isobaric heat capacity, including the electronic, magnetic, and lattice vibration contributions, are very close to the corresponding experimental values. The analysis shows that the largest contribution comes from the lattice vibrations, while the contribution due to the disordered magnetic moments is the smallest.
The EMTO method combined with the quasi-harmonic Debye model can reproduce the experimentally observed tendency of elastic moduli to decrease with increasing temperature. Good agreement is found between the calculated and experimental values of the bulk modulus and C11 elastic constant, whereas the value of C44 and the temperature dependencies of C12 and C44 seem to be overestimated by the present calculations.
The developed methodology of free energy modeling has been used for calculating the temperature dependence of stacking fault energy in paramagnetic nickel. Several calculation schemes have been employed and their results compared against each other. Two versions of the supercell approach yield almost identical and accurate stacking fault energies. The second-order axial-nextnearest-neighbor Ising model is found to yield very similar values, but at a smaller computational cost. It can then be regarded as an efficient way to calculate stacking fault energies with sufficient accuracy. Thermal expansion is found to have a significant effect on the temperature dependence of stacking fault energy. Additionally, magnetic entropy is found to contribute to the calculated decrease in the stacking fault energy as a function of temperature, in contrast to the magnetic entropy effect on the stacking fault energy of -Fe and austenitic stainless steel.
In conclusion, the proposed quasi-harmonic Debye model, parameterized using first-principles calculations employing the EMTO method, is found to be a physically sound and efficient method to evaluate free energy contributions due to the magnetic, electronic and lattice vibration degrees of freedom, which opens the way to efficiently describing the temperature-dependent properties of metallic materials.