Generalized stacking fault energy of carbon-alloyed paramagnetic γ-Fe

Generalized stacking fault energy (GSFE) is an important parameter for understanding the underlying physics governing the deformation mechanisms in face-centred cubic (fcc) materials. In the present work, we study the long-standing question regarding the influence of C on the GSFE in austenitic steels at paramagnetic state. We calculate the GSFE in both -Fe and Fe–C alloys using the exact muffin-tin orbitals method and the Vienna Ab initio Simulation Package. Our results show that the GSFE is increased by the presence of interstitial C, and the universal scaling law is used to verify the accuracy of the obtained stacking fault energies. The C-driven change of the GSFE is discussed considering the magnetic contributions. The effective energy barriers for stacking fault, twinning and slip formation are employed to disclose the C effect on the deformation modes, and we also demonstrate that the magnetic structures as a function of volume explain the effect of paramagnetism on the C-driven changes of the stacking fault energies as compared to the hypothetical non-magnetic case.

(Some figures may appear in colour only in the online journal) Original content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. theoretical tools is not a straightforward task, the oversimplified models such as the non-magnetic (NM) state are often adopted. Calculations performed for the NM fcc Fe show that C increases the SFE by more than 300 mJ m −2 per wt.% C [6], whereas the C effect is significantly suppressed at the PM state according to our previous work (96 mJ m −2 per wt.% C [7]).
Noticeably, the reported C effects on SFE are not always consistent in the literature [8,9]. For example, the reported C effect on the SFE in austenites by Schramm and Reed [8] ( ∼ 410 mJ m −2 per wt.% C) is one magnitude larger than that given by Brofman et al ( ∼ 26-40 mJ m −2 per wt.% C) [9]. An accurate determination of SFE is limited by the complexity of the atomic-scale defects [10] using the existing experimental methods, e.g. x-ray diffraction (XRD), neutron diffraction (ND), and transmission electron microscopy (TEM) [11][12][13].
Recently, Hickel et al [14] demonstrated that the electron beam may locally heat the samples and causes C atoms to diffuse out of the stacking fault during the examination period in TEM observation. Such a mechanism reduces the measured SFE by around 30% and may cause an artificial interpretation of the C effect on the SFE [14]. Theoretical calculations, such as the ab initio methods, provide alternative ways to tackle the problem and give a consistent description of the composition and magnetic effects on the SFE. Furthermore, the ab initio calculations are capable of describing generalized stacking fault energy (GSFE or γ-surface), which describes the energy change in the process of shearing a crystal along given crystallographic directions on the slip plane, e.g. in fcc materials commonly the 1 1 2 direction on the {1 1 1} plane. The GSFE includes several stacking fault energies, such as intrinsic stacking fault energy (γ isf ), unstable stacking fault energy (γ usf ), unstable twinning fault energy (γ utw ) and extrinsic stacking fault energy (γ esf ), which are significant materials parameters for quantitatively predicting the critical twinning stress [15] and evaluating twinnability [16]. Recently, the most advanced theory for plasticity in fcc materials shows that the competition between plastic deformation mechanisms can be rationalized from the intrinsic energy barriers given by GSFE [17,18]. Furthermore, the GSFE can be used as input in the elasto-plastic models, e.g. the Peierls-Nabarro and phasefield models [19,20].
The GSFE in the NM C-alloyed γ-Fe has been studied by Gholizadeh et al [21] and Medvedeva et al [22]. Their work reveals that C significantly increases the entire γ-surface, and the C effect rapidly decays with increasing distance between C atoms and the stacking fault. This rapid energy drop is believed to provide a strong thermodynamic driving force for the out-diffusion of C from a stacking fault/twin at elevated temperatures (anti-Suzuki effect) [14,23]. In our previous work [7], the PM state in γ-Fe and Fe-Cr-Ni austenitic steels was described by the disordered local magnetic moment (DLM) approximation [24]; the C effect on SFE (γ isf ) was investigated at the PM state using the axial interaction model (AIM), which is however not suitable for calculating the GSFE. The significantly reduced C effect on SFE at the PM state draws further interest in the paramagnetism-altered C effect on the GSFE as compared to the NM behaviours.
Hitherto, no complete description of the C effect on the GSFE in PM γ-Fe or austenitic steels has been reported in the literature. This paper aims to fill the gap and uncover the C effect on the GSFE of the PM γ-Fe.
The rest of this paper is organized as follows. In section 2, the computational models, tools and numerical details are presented. In section 3, we first compare the obtained GSFE at the NM state by the exact muffin-tin orbitals (EMTO) method [25,26] and the Vienna Ab initio Simulation Package (VASP) [27,28], in order to evaluate the accuracy of the results obtained using EMTO with current settings. Then we present the GSFE of the C-alloyed γ-Fe at the PM state. In section 4, the results are discussed in terms of the universal scaling law [29] and the effects of paramagnetism are further investigated. The conclusions are given in section 5.

Atomic model for GSFE calculation
In the present work, the GSFE calculations follow the methodology developed in [21]. A hexagonal supercell with 6 (1 1 1) layers along c direction is employed. The lattice vectors a, b and c are (a/2)[1 1 0] fcc , (a/2)[0 1 1] fcc and (2a)[1 1 1] fcc , respectively, in which a is the lattice parameter of a conventional fcc cell. To obtain the energy variation (GSFE), the c axis is tilted along the [11 2] direction by a vector u, i.e. c changes to c + u. When u = b = (a/6)[112] fcc , the supercell contains one ISF. The extrinsic stacking fault (ESF), which is also seen as a two-layer twin nucleus, is obtained by creating a new stacking fault on an adjacent (1 1 1) plane to the exiting stacking fault plane. The GSFE is obtained by monitoring the energy variation when the structure changes gradually from the ideal fcc to the ISF, and then to the ESF, via where E is the total energy of the corresponding structure and A is the area of the stacking fault. We illustrate the main temperature effect on the GSFE through the thermal lattice expansion and neglect the electronic, magnetic and vibrational terms. We notice that at temperatures less than 300 K, contributions from the electronic, magnetic and vibrational entropies are estimated to be small [30,31]. For the GSFE calculation in the C-alloyed γ-Fe, we double the size of the above supercell along the a and b directions. One C atom is placed at the octahedral interstitial site. This system is denoted as Fe24C1. The corresponding C concentration is 4.0 at.% ( ∼ 0.9 wt.%). Depending on the distance between the C atom and the stacking fault plane, we denote the obtained GSFE as γ d0 , γ d1 , γ d2 , γ d3 , ... , where d x ( x = 0, 1, 2, 3, ...) indicates that the C atom is placed at x interlayer spacing away from the slip plane. When C is in the slip plane, its position is properly relaxed at the constraint of a generalized stacking fault. To achieve the GSFE for the Fe-C alloy with homogeneously distributed C, the energy surface corresponding to the shear process from fcc to ISF can be evaluated by taking an average of the γ d0 , γ d1 , γ d2 and γ d3 values, i.e.

Ab initio methods
We adopt both the VASP and EMTO methods for NM calculations. In the VASP calculations, the Projector Augmented Wave (PAW) scheme [32,33] is used. The exchange-correlation functional is the Perdew-Burke-Ernzerhof (PBE) [34] version of the generalized-gradient approximation (GGA). A plane-wave basis set with energy cut-off 500 eV is adopted. The k-point meshes of 13 × 13 × 5 and 9 × 9 × 5 are generated using the Monkhorst and Pack scheme [35] for the γ-Fe and Fe24C1 supercells, respectively. The interatomic forces are converged to less than 0.02 eV Å −1 in the cases where the atomic relaxation is performed.
The PM state with DLM approximation [24] is approximated using the coherent potential approximation (CPA) [25,36] as implemented in the EMTO package. In the EMTO calcul ations, the one-electron equations are solved within the scalar-relativistic approximation and the soft-core scheme. The calculation of the Green function for valence states is implemented with 16 complex energy points. The s, p, d, f orbitals are included in our muffin-tin basis set. The PBE approx imation [34] is employed for the exchange-correlation functional. The potential sphere radius of interstitial C (w C ) is optimized to be 0.77 w 0 C , where w 0 C is the atomic sphere radius of the corresponding Voronoi polyhedron around the C site. With such settings, the obtained C effect on the SFE of γ-Fe is shown to be well reproduced by the AIM model [7]. The detailed assessment of the EMTO method for the GSFE calculations using the supercell technique is presented in section 3.1. We select the k-point meshes of 9 × 17 × 1 and 7 × 13 × 1 for the NM and PM calculations in the EMTO method, respectively, in order to ensure the same conv ergence criterion for γ isf at both magnetic states (the difference between the γ isf with two proximate k-point samplings is less than 5.0%). In the EMTO calculations, the symmetry of the computational cell is set to be triclinic.
Previous studies have shown that the C-induced atomic relaxation is most significant in the first coordination shell and decreases fast in the further shells [37,38]. Therefore, in the present work, only the atomic relaxation in the the first coordination shell is considered for both EMTO and VASP methods. Figure 1 shows the relative total energies of the Fe24C1 configurations with different expansions for the first coordination shell of C, which is calculated by the EMTO method for the NM and PM states. We find that EMTO gives similar relaxations for the first coordination shell in the NM and PM Fe24C1: both relaxations corresponding to the lowest total energy are approximately 7.0%, which are consistent with the value given by VASP for the NM Fe24C1 (6.8%) [39,40]. Considering the weak magnetic state dependence of the local relaxation and the good parallelism between the EMTO and VASP results, in the following we employ the optimized atomic configurations obtained by VASP in both the NM and PM EMTO calculations for Fe24C1.

NM and PM γ-surfaces of fcc Fe at 0 K
According to recent studies [15,16], γ utw is a critical parameter to build the relationship between GSFE and plasticity. Unfortunately, γ utw is missing in most of previous γ-surface calculations for γ-Fe [14,19,21]. In this section, the GSFE of pure γ-Fe is calculated at both NM and PM states with their corresponding equilibrium volumes at 0 K and the γ utw is provided. Further, the magnetic effect on the GSFE is demonstrated for pure γ-Fe.
The calculated GSFE for the NM γ-Fe at 0 K is featured in figure 2. The numerical values of the stacking fault energies  are summarized in table 1, together with some previous results. Pure γ-Fe is known to be stable at high temper atures, which is also indicated by the negative γ isf at 0 K. Previous theoretical studies [6,14,19,21] with different ab initio methods have reported consistent stacking fault energies, with γ isf and γ usf in the ranges of approximately −460 to −390 and 500 to 520 mJ m −2 , respectively. The deviations are due to the calculation details, e.g. different exchange-correlation functionals, relaxation schemes or Brillouin zone samplings. In figure 2, we compare the results obtained with VASP and EMTO for the unrelaxed lattices (rigidly shifting the upper half of the crystal). The γ isf and γ esf are in perfect agreement, while the γ usf and γ utw differ by 6.6% and 10.3%, respectively. The discrepancy may partly be ascribed to the muffin-tin approximation that leads to relatively large errors in the less close-packed structures (e.g. for the unstable stacking faults) [14]. Based on this observation, we introduce a parameter ∆γ to bring the EMTO results to the level of the VASP results. The ∆γ is defined as the difference between the stacking fault energies obtained with VASP and EMTO in the same configuration (u), ∆γ(u) = γ VASP (u) − γ EMTO (u) [14]. We assume that ∆γ is independent of the magnetic state, and use the ∆γ values obtained from the NM results to correct the stacking fault energies calculated by EMTO at the PM state. For reference, wherever is possible we also present results obtained without this correction. The VASP results with a fullrelaxation and the EMTO results with an interlayer relaxation in the vicinity of a slip plane are also plotted in figure 2 for comparison in order to probe the effect of atomic relaxation on the GSFE of γ-Fe. Whereas both methods show weak influence of the relaxation on γ isf (<6.0%), the atomic relaxation causes a much stronger effect on γ usf and γ utw , which is consistent with previous studies [14].
In order to demonstrate the magnetic effect on the GSFE of γ-Fe, we list the stacking fault energies obtained at the PM state in the lower part of table 1. The stacking fault energies are given with the correction ∆γ. The equilibrium lattice parameter of γ-Fe at the PM state (3.498 Å ) increases as compared to that for the NM state (3.457 Å ). The optimized lattice parameter of the PM Fe agrees better with the experimental lattice parameters obtained by extrapolating from Fe-Cu alloys (3.562 Å ) [41] and from the high-temperature lattice parameter of γ-Fe (3.569 Å ) [42]. Here, we should notice that the difference between the theoretical and extrapolated experimental data is mainly due to the employed exchangecorrelation approximation. Recently, Dong et al [43] used the quasi-non-uniform approximation (QNA) [44,45] and found nearly perfect agreement between experimental values and theoretical predictions. Regarding the stacking fault energies, the γ isf is significantly increased (from −393 to −171 mJ m −2 ) by the introduction of local spin polarization. γ isf was previously calculated at various magnetic states for γ-Fe [14], including NM, low-spin and high-spin ferromagnetic (LS-FM and HS-FM), PM, single-and double-layer antiferromagnetic (AFMS and AFMD) states. Different magnetic states lead to different equilibrium lattice parameters, which significantly affect the obtained γ isf . The NM and LS-FM calculations give similar equilibrium lattice parameters ( ∼ 3.460 Å for NM and ∼ 3.470 Å for LS-FM) and the magnetic moment at LS-FM state is very small, thus resulting in similar γ isf of γ-Fe ( ∼ 3.0% difference) [6,22]. In contrast, the HS-FM simulation shows a much larger lattice parameter ( ∼ 3.640 Å ), which in turn leads to much larger γ isf [22]. For PM γ-Fe, the equilibrium lattice parameter ( ∼ 3.498 Å ) is only slightly larger than that at NM state, but the paramagnetism itself (local magnetic moments) significantly influences the γ isf . The present DLM approximation for the PM calculations gives similar conclusion as obtained previously by the special quasi-random structures (SQS) approximation for the PM state [46]. Whereas the γ isf is increased with the PM coupling, the γ usf is decreased from 562 to 337 mJ m −2 . However, the magnetic effect on the γ utw is relatively small. The reason is that an unstable twinning fault may geometrically be approximated as an ISF plus an unstable stacking fault (USF) thus the two magnetic effects cancel each other to a large extent.

The impact of C on the γ-surface in NM Fe
Several attempts have been made to study the γ-surface of the γ-Fe with interstitials at the NM state [6,19,21]. Before extending the GSFE calculation to PM state, the performance of EMTO on describing the interstitial structure is assessed by comparing the GSFE obtained by EMTO and VASP for the NM Fe24C1. The accuracy of the EMTO method is similar to the full-potential density functional theory (DFT) methods for close-packed structures but might be somewhat lower for open or low symmetry structures, which is actually the case of interstitial impurities.
The dependences of the γ isf and γ usf on the distance of the C atom from the slip plane in the NM Fe24C1 are illustrated in figure 3. Both the EMTO and VASP results show that only the C atom close to the slip plane has a strong effect. When the distance between C and the slip plane is larger than 2-3 interlayer spacing, the chemical effect diminishes and the remained difference may be ascribed to the lattice expansion effect [6]. In addition, a Friedel-like oscillation is found for γ isf and γ usf , explaining why γ d3 is larger than γ d2 [6,21]. The γ isf values obtained by EMTO and VASP are in good agreement for different C sites. However, the quantitative agreement for the γ usf is less impressive, especially for the C sites d2 and d3. The larger difference between the γ usf for the d3 configuration and the corresponding bulk Fe obtained in EMTO compared to that in VASP also points in this direction. To look into the origin of this discrepancy, we consider the accuracy of the muffin-tin approximation when the USF is formed. In the configuration with both USF and interstitial C, the muffin-tin zero discontinuity, which represents the jump in the muffintin potential at the potential sphere boundary [26], is found enlarged as compared to the fcc configuration. The overlap errors are expected to decrease the accuracy when muffintin zero discontinuity is increased [48]. We perform some additional calculations where the C atom and its first-nearest neighbouring Fe atoms are excluded from establishing the muffin-tin zero potential. The new γ d2 and γ d3 values decrease to 509 and 521 m −2 from 568 and 623 mJ m −2 , respectively (see figure 3(b)), which agree better with the VASP results. On this ground, we conclude that at least part of the deviations observed between the EMTO and VASP results for the NM Fe24C1 (see figure 3(b)) are due to the increased overlap error around the unstable stacking fault. In the following, we do not consider separately the overlap error since that is explicitly included in the ∆γ correction.
The γ isf , γ usf , γ utw and γ esf of the NM Fe24C1 calculated according to equations (2) and (3) are listed in the upper part of table 2, together with available previous theoretical results [6,21]. Despite the discrepancies shown in figure 3 (solid lines), the averaged stacking fault energies obtained by EMTO and VASP are consistent with each other. Similarly to the case of the pure γ-Fe, the energy differences between the VASP and EMTO results obtained for the NM Fe24C1 (∆γ) are used to correct the GSFE of the PM Fe24C1. The  uncorrected data for the PM Fe24C1 (without ∆γ) are shown for reference in parentheses. Here, we should stress that the differences between the EMTO and VASP values for γ isf , γ usf and γ utw are very similar for the NM γ-Fe and NM Fe24C1. Therefore, we believe that the EMTO method is suitable for describing the interstitial impurity in γ-Fe.

The impact of C on the γ-surface in PM Fe
Based on the work shown in sections 3.1 and 3.2, the EMTO method is believed to be able to describe the GSFE of a Fe-C system, and the PM coupling is a more accurate approximation for γ-Fe and Fe-C alloys. The GSFE calculations for Fe24C1 are then extended to the PM state, which is the main focus in the present work.
The equilibrium lattice parameters for the PM Fe24C1 and γ-Fe are 3.536 and 3.498 Å , respectively. The calculated volume expansion is approximately 0.043 Å per wt.% C at the PM state, which is similar to the C effect on the volume expansion at the NM state (0.045 Å per wt.% C). The stacking fault energies of the PM Fe24C1 corresponding to the equilibrium (0 K) volume are listed in the lower part of table 2. The γ isf increases dramatically from −171 to 57 mJ m −2 at the PM state when 4.0 at.% of C is added to γ-Fe. The γ usf and γ utw increase with the C addition from 337 to 450 mJ m −2 and from 281 to 471 mJ m −2 , respectively. The above figures correspond to the results including the ∆γ correction obtained from the NM calculations.
Considering the strongly volume-dependent C effect on the γ-surface [18], here we also present the volume dependence of the GSFE in the PM γ-Fe and Fe24C1 (see figure 4). The γ usf and γ utw decrease while γ isf and γ esf increase with volume in Fe24C1 (see figure 4(b)) and γ-Fe (see figure 4(c)). In particular, a dramatic growth of the γ isf (δγ isf /δa, see figure 4(a)) for the PM γ-Fe appears around its equilibrium lattice parameter (a = 3.498 Å , see table 1), which is attributed to the magnetic transition from low spin to high spin state. Details of the magnetic structures are presented in the next section. In figure 4(a), we also compare the volume dependence of the γ isf obtained using the supercell technique with previous results using the supercell [18] and the AIM model [7]. It is found that the supercell and the AIM techniques display quantitative agreement. The somewhat weaker agreement at lower volumes is ascribed to the deviating magnetic structure descriptions around the stacking fault in these two methods [18].
We should emphasize here that the C effect on the γsurface is strongly volume-dependent because of the magnetic-volume coupling effect in PM γ-Fe [7,18]. In practice, the C effect on the SFE (γ isf ) at room temperature is usually discussed in stable/metastable austenitic steels. The austenitic phase of γ-Fe can be stabilized at room temperature by alloying with e.g. C, N, Ni. The room temperature lattice parameter of γ-Fe can be obtained by extrapolating from Fe-Mn alloys with high Mn contents [49,50], Fe-C alloys [41], Fe-Cr-Ni austenitic steels [51] or high-temperature γ-Fe [41,42]. In the present work, the room temperature volumes of γ-Fe and Fe24C1 are taken according to an exper imental relation (a = 3.562 + 0.037x C ) given in [41], in which a denotes the room temperature volume of Fe-C alloys and 3.562 is the lattice parameter of pure γ-Fe at room temper ature with a unit of Å . x C is the weight percent of the C atom in γ-Fe. Correspondingly, the room temperature volume of Fe24C1 is estimated to be 3.595 Å .
The stacking fault energies corresponding to the room temper ature volumes for the PM γ-Fe (3.562 Å ) and PM Fe24C1 (3.595 Å ) are listed in the lower parts of tables 1 and 2, respectively. The differences in the stacking fault energies between the PM Fe24C1 and PM γ-Fe are calculated and denoted corre spondingly as ∆γ isf , ∆γ usf , ∆γ utw and ∆γ esf . Further, the dependences of the stacking fault energies on C concentration are calculated by ∆γ isf /∆x C , ∆γ usf /∆x C , ∆γ utw /∆x C and ∆γ esf /∆x C , respectively, in which ∆x C is equal to ∼ 0.9 wt.%. These slopes are used to represent the C effects on the stacking fault energies. Using the stacking fault energies obtained for the PM γ-Fe (see the lower part of table 1) and PM Fe24C1 (see the lower part of table 2) at their corresponding equilibrium volumes at 0 K (3.498 Å for the PM γ-Fe  3.536 Å for the PM Fe24C1), ∆γ isf /∆x C , ∆γ usf /∆x C , ∆γ utw /∆x C and ∆γ esf /∆x C are approximately 256, 127, 213, 188 mJ m −2 per wt.% C, respectively. Similarly, using the stacking fault energies corresponding to the room temperature volumes of the PM γ-Fe (3.562 Å , see the lower part of table 1) and PM Fe24C1 (3.595 Å , see the lower part of table 2), ∆γ isf /∆x C , ∆γ usf /∆x C , ∆γ utw /∆x C and ∆γ esf /∆x C are approximately 130, 86, 125 and 101 mJ m −2 per wt.% C, respectively. Both of the above sets of data take into account the correction ∆γ. Obviously, the theoretically predicted C effects on the stacking fault energies are dramatically reduced with volume in the PM γ-Fe.

Simple assessment of the theoretical results and C effect on energy barriers
The stacking fault energies for γ-Fe and Fe24C1 calculated in this work are screened against the universal scaling law (see figure 5). The universal scaling law approximately holds for various intrinsic energy barriers in metals and alloys [29,52]. It expresses a relationship between unstable stacking fault energy, unstable twinning fault energy and intrinsic stacking fault energy. The physics behind the universal scaling law is the short range nature of the above planar faults and the fact that an unstable twinning fault is composed of a USF plus half of an ISF (since the intrinsic and unstable stacking faults share a common layer in an unstable twinning fault configuration). Figure 5 shows that all the present results fit well the universal scaling law, regardless of volume, magnetic state and configuration, thus verifying the accuracy of our theoretical calculations with respect to the universal relationship between the individual stacking fault energies.
We observe that the calculated γ isf of the PM Fe24C1 is 106 mJ m −2 at room temperature as shown in table 2, which is higher than the experimental SFE (γ isf ) shown in Fe-Mn-C TWIP steels (0-50 mJ m −2 ) [53][54][55][56]. This overestimation may be due to the current assumption that C is homogeneously distributed in γ-Fe. Instead, the interstitial C is found to be repelled from the stacking fault [57] in high-Mn steels and the C diffusion out of the slip plane would substantially decrease the SFE (γ isf ). We should also notice that the in-plane C atomic concentration of our supercell is 1/5 (20.0 at.%). However, in the pure Fe-C austenite with a high C content, the tendency towards Fe 8 C ordering is suggested [58], which shows a lower local C concentration ( ∼ 11.0 at.%) than our computational model. Therefore, the local C concentration in the present work is relatively high (nevertheless, the size of the supercell is sufficient to avoid the in-plane interaction of two adjacent C atoms). In order to illustrate this presumption, we may consider the present result obtained in the situation where C is located far away from the stacking fault. The corresponding γ isf is 35 mJ m −2 at the room temperature volume, which is in a much better agreement with the above experimental values [53][54][55][56]. In general, one should be cautious when comparing the present result to the experimental ones which are affected by many factors. For example, the TEM measurements may be strongly affected by the local temperature rise caused by electron beam, causing the local C depletion from stacking fault [57]. Moreover, the C distribution in reality is not perfectly random, and up to 25% fluctuation in concentration can occur on the scale of a few nanometres [57]. C may also segregate in grain boundary and causes the deviation of local C concentration around stacking fault from the nominal C content [23].
Our work for the first time presents the C effect on the stacking fault energy barriers in the PM γ-Fe. Using the theory for plasticity of fcc metals proposed by Jo et al [59], we transform the stacking fault energies of the PM γ-Fe and PM Fe24C1 to the effective energy barriers as a function of θ according to the equations given in [59]. Here, θ is the angle between the resolved shear stress direction and the stacking fault's easy direction. The effective energy barriers for stacking fault, twinning and full slip are respectively. The stacking fault energies of the PM γ-Fe and PM Fe24C1 obtained at the corresponding room temperature volumes are used to calculate the effective energy barriers. For the PM Fe24C1, γ sf , γ tw and γ sl are calculated using both the averaged γ usf , γ isf and γ utw and those for C located far away from the stacking fault (with the consideration that the C concentration may be small close to the stacking fault zone). The γ usf , γ isf and γ utw are 275, −8 and 275 mJ m −2 for the PM γ-Fe, respectively. The averaged γ usf , γ isf and γ utw (over layers) for the PM Fe24C1 are 351, 106 and 385 mJ m −2 , respectively. In the case with C located far away from the stacking fault, the γ usf , γ isf and γ utw are 348, 35 and 363 mJ m −2 , respectively. The effective energy barriers for stacking fault, twinning and full slip of the PM γ-Fe and PM Fe24C1 as a function of θ are then shown in figure 6. The effective energy barriers of the PM γ-Fe are increased with C addition. In the PM Fe24C1, the effective energy barrier for twinning (γ tw ) is lower than that for stacking fault (γ sf ), showing that twinning is the preferred deformation mechanism. However, the PM γ-Fe shows an opposite trend: namely the γ tw is larger than γ sf , indicating that pure Fe shows martensitic transformation rather than twinning and the presence of C stabilizes the austenitic phase. Moreover, in the PM γ-Fe, γ tw is equal to γ sl when θ is slightly above 30°. In the PM Fe24C1, the θ corresponding to γ tw = γ sl is smaller than 30°, indicating a smaller probability of activating twinning with C alloying theoretically. We also compare the effective energy barriers obtained with the averaged fault energies to those obtained in the case when C is located far away from the stacking fault for the PM Fe24C1 (see the inset in figure 6). Note that averaging the stacking fault energies over all the layers is an approximation of homogeneous C distribution in γ-Fe. On the one hand, the C segregation to the stacking fault increases the effective energy barriers for twinning and full slip but enhances the probability of activating twinning relative to the homogeneous C distribution. On the other hand, the C segregation case shows similar behaviours to pure γ-Fe, which makes sense since the stacking fault energies of Fe24C1 are close to those of pure γ-Fe when C is far away from the stacking fault. Nevertheless, the general C effect on the effective energy barriers is not altered by the C segregation as compared to the homogeneously distributed C. Hence, although C stabilizes the fcc phase it decreases the twinning probability in γ-Fe from a theoretical point of view regarding the GSFE calculations.

Magnetic structure
As concluded in [7], the introduction of paramagnetism significantly weakens the C effect on γ isf with respect to the C effect at the NM state. This phenomenon is also illustrated in figure 7, where we show that the slope ∆γ isf /∆x C of 130 mJ m −2 /wt.% at room temperature for the PM state is consistent with [7]. Furthermore, the experimental measurement for Fe-18Cr-10Ni-C alloys with different C contents [60] gives a similar C impact like the one reported in this work. However, a four-dimensional regression analysis in Schramm's [8] and Brofman's [9] work yields different C-driven SFE (γ isf )  C-driven change of γ isf (SFE) for PM Fe24C1 calculated at room temperature volume as a function of C concentration, together with previous results reported by Schramm [8], Gholizadeh [21] and Abbasi [6], Petrov [60], Lu [7] and Brofman [9]. changes. One potential reason is that the regression analysis is adopted in the austenites with different C concentrations [7]. Note that the regression analysis is also based on the austenite stainless steels with different Cr and Ni contents. Hence, the regression analysis may lead to totally inconsistent results [30]. Theoretically, the C-driven SFE (γ isf ) change is computed using ab initio method by Gholizadeh et al [21] and Abbasi et al [6] at the NM state, which shows much larger SFE (γ isf ) change as compared to the PM calculation. Therefore, the neglect of the magnetic transition during the generation of stacking faults can be fatal for understanding appropriately the C effect on the GSFE.
Note that the ∆γ isf versus C concentration in figure 7 is obtained with the corresponding room temperature volumes of the PM γ-Fe and PM Fe24C1. Such a weakening impact of paramagnetism on the C effect can be also observed for the γ-Fe and Fe24C1 with their corresponding equilibrium volumes. The γ isf of the pure γ-Fe and Fe24C1 with equilibrium volumes are extracted from tables 1 and 2, respectively, and are illustrated in figure 8. Figure 8 illustrates that the slope of the line connecting the γ isf of γ-Fe and Fe24C1 at the PM state is smaller than that at the NM state, indicating that the C effect on enhancing γ isf is weakened by paramagnetism. Furthermore, the introduction of paramagnetism stabilizes fcc phase (see the upward arrow for γ isf ) and the paramagnetism-induced stabilization is more profound in γ-Fe with respect to Fe24C1. In contrast, the enhancing effect of C on γ usf is strengthened by paramagnetism. The downward arrows from the NM to the PM γ usf in figure 8 indicate that the introduction of paramagnetism stabilizes the USF configuration instead of fcc phase for both γ-Fe and Fe24C1, and this stabilization effect is stronger in γ-Fe than that in Fe24C1. Here, we need to stress that the volume expansion caused by C addition is almost the same at NM and PM states (0.045 and 0.043 per wt.% C, respectively) thus the volume effect on the γ isf and γ usf changes from γ-Fe to Fe24C1 can be neglected.
The above behaviour upon C alloying can be detected from the magnetic structures of γ-Fe and Fe24C1 (see figure 9). Figure 9(a) shows the local magnetic moments of the Fe atoms in the pure Fe with fcc, USF and ISF configurations as a function of lattice parameter. The top and the second panels are for Fe atoms at the stacking fault plane and in the next nearest neighbouring layer from the slip plane, respectively. For the Fe atoms farther away from the stacking fault, their local magnetic moments are almost the same as bulk Fe (not shown). In figure 9(a), from the equilibrium volume of the NM γ-Fe to that of the PM γ-Fe, the Fe atom in the γ-Fe containing the ISF gains local magnetic moment and its size is smaller than that of the bulk Fe, while the magnetic moment of the Fe atom in the Fe24C1 containing the ISF still remains zero within the corresponding volume range (see figure 9(c)). On the one hand, the smaller magnetic moment of the Fe atom in the ISF configuration for γ-Fe indicates the destabilized effect of paramagnetism on the ISF structure, and such an effect also stays for Fe24C1. However, the C addition to γ-Fe pushes this effect to a larger volume in Fe24C1 through suppressing the local magnetic coupling. Therefore, we conclude that paramagnetism stabilizes/destabilizes the fcc Fe/Fe with the ISF more strongly than the fcc Fe24C1/Fe24C1 with the ISF. Actually, this phenomenon is equivalent to the weakened C effect on γ isf with the consideration of paramagnetism (see figure 8). Similarly, the speculation that paramagnetism stabilizes/destabilizes the Fe with the USF/fcc Fe more strongly than the Fe24C1 with the USF/fcc Fe24C1 can be correlated with the magnetic behaviours shown in figures 9(a) and (d). The pure Fe with the USF shows relatively strong spin polarization at small volumes, while the first-nearest neighbouring Fe atoms of C are either non-magnetic or just starting to gain local magnetic moment within the volume range from the NM to the PM equilibrium in Fe24C1. Here, the stabilization effect on the USF configuration is constrained by the interstitial C due to its suppressing effect on the local magnetic moment. In general, the different paramagnetic effects appeared in γ-Fe and Fe24C1 rely on the C-induced magnetic suppression in the spin polarized environment.
A previous study shows that the local magnetic interactions are strongly affected by the interstitial C in Fe-Mn-C alloys; the ferromagnetic Fe-Fe and antiferromagnetic Fe-Mn interactions are in competition [61]. Such competition together with the C-driven frustrated magnetic coupling may result in the spin-glass-like behaviour observed in austenitic Fe-Mn-C alloys [61,62]. Moreover, the magnetic discrimination shown in fcc and hcp Fe-Cr-Ni alloys contributes to the different alloying effects on SFE (γ isf ) [30]. Therefore, an appropriate description on the magnetic state is significant for understanding the C effect on the GSFE. On the one hand, the introduction of paramagnetism causes the opposite changes of the C effects on γ isf and γ usf , i.e. the C effect is weakened for γ isf while it is enhanced for γ usf as compared to the NM case. On the other hand, the interstitial C, which suppresses the magnetic moments of its adjacent Fe atoms, results in the different stabilities of the pure Fe and Fe24C1 with stacking faults relative to their corresponding fcc phases.

Conclusions
The C effect on the GSFE of PM γ-Fe is investigated using two density functional theory solvers. Since VASP can account for the local lattice relaxation and EMTO for the disordered magnetic state, we employ both methods to provide a reasonable theoretical description of the complex problem of C alloying. The accuracy of the GSFE calculations for γ-Fe and Fe24C1 is tested against the universal scaling law. Based on the calculated stacking fault energies at room temperature volumes, the effective energy barriers are constructed for PM γ-Fe and Fe24C1. It is found that the C addition increases the effective energy barriers of pure γ-Fe and the twinning probability is slightly lower with C alloying. We also analyse the magnetic structures of the C-alloyed γ-Fe containing stacking faults in detail. The introduction of paramagnetism suppresses the C-driven change on γ isf while it enhances the C-driven change on γ usf . The different behaviours for γ isf and γ usf are correlated to the magnetic transitions within the NM and PM equilibrium volume ranges in the configurations with fcc phase and generalized stacking fault for pure Fe and Fe24C1. The suppression effect of C on the magnetic moments of its adjacent Fe atoms results in the different relative stabilities of pure Fe and Fe24C1 with stacking faults as compared to their corresponding fcc phases. In summary, the present theoretical results confirm the relatively weak C effect on γ isf at the PM state in line with recent experimental findings. The difference between theoretical and experimental γ isf at room temperature is mainly ascribed to the exchange-correlation functionals, the magnetic excitations and the anharmonic phonon contribution. Furthermore, the C effect on γ usf is found to be strengthened at the PM state, thus emphasizing the importance of investigating the GSFE with a proper magnetic state in Fe-C alloys.