1 Introduction

Metallic materials designed for high-temperature applications must resist both mechanical and environmental degradation. The environmental degradation of materials at high temperatures is an important issue in energy production, engines and industrial processes. It limits the useful life of installations, restricts the utilization of the fuel and hinders the development of more economical and environmentally friendly processes and systems.

In corrosion resistant Fe-based alloys, environmental degradation due to high-temperature corrosion is controlled by the properties of the oxide scales that develop on the material. Iron-based alloys are regularly used in high-temperature applications due to their good mechanical properties and corrosion resistance at a reasonable cost. Low-alloy steels are normally used at low/medium temperatures, stainless steels (iron-based alloys with \(>12\%\) Cr) are used at high-temperatures (\(T>500 ^\circ \)C) and/or in more corrosive environments while FeCrAl alloys are preferred at even higher temperatures. The protective layer of each of the mentioned alloys consists of iron oxide (steels), Cr-rich oxide of \(M_2O_3\) (stainless steels) and Al-rich \(M_2O_3\) (FeCrAl alloys), respectively.

Despite the large variety of alloy grades and compositions, the dominant oxides in the protective layer mentioned above are known to provide protection against fast oxidation and corrosion in different environments.[1] The protection is a result of the formation of a non-porous, well-adherent oxide on the metal surface and its resistance against cracking or spalling during heating and cooling. The growth rate of this protective scale is determined by its ability to transport ions through the oxide scale, which explains the different temperature regimes where they are considered protective.

\(Cr_2O_3\) and \(Al_2O_3\) feature lower intrinsic transport rates of ions than iron oxide, low vapor pressures and high melting points. To form such protective scales on the iron-based alloys several factors are important, including temperature, corrosive environment, alloy microstructure and the activity of the Al and Cr in the bulk material. There is a great number of publications concerning the oxidation properties of the protective scales in steels, stainless steels and FeCrAl alloys and they are summarized by Kofstad,[2] Birks and Meier[3] and Young.[4]

Calphad, Computer Calculation of Phase Diagrams, as a phenomenological method has played a vital role in the field of materials design for more than three decades.[5] In Calphad method thermodynamic properties of the materials are modelled based on atomistic simulations and experimental data. Stable phases under different conditions, e.g. temperature, pressure, composition etc., are calculated by help of modelling thermodynamic properties of the system. Thermodynamic and kinetic databases are built using a Calphad-based software. Such databases are the main input to predict materials microstructure and possible products in different circumstances. Calphad is a self-consistent method and well-developed databases can be used to study higher-order systems, i.e. from binary to ternary and so on.[6,7,8]

Calphad approach can be a valuable tool for the life assessment of materials used in different environments. However, despite the potential of this method, its application in the field of oxidation and corrosion of materials in an industrial level has been limited and a model to predict environmental degradation, i.e. high temperature corrosion, has up till now been lacking. This has been partly due to the complicated nature of these phenomena, but also because of the lack of computational efficiency to perform the demanding simulations required for high-temperature corrosion. The scientific developments in multiple fields, however, have over the last decade matured to such an extent that realistic modelling of the oxidation process is within reach.

Modelling of oxide scale growth traces its roots to the analytical work by Wagner.[9] An early example of a computer simulation using independently assessed diffusion coefficients is the work by Narita et al.[10] A number of research have used Calphad approach to investigate oxidation in different types of materials. For example Krupp et al.[11] used Chemapp software[12] to model diffusion-based oxidation in Ni-based superalloys and Cr-containing steels. Toffolon-Masclet et al.[13] and Glazoff et al.[14] also used Calphad-based programs to study high temperature corrosion in Zr alloys in nuclear industry. There have been some other attempts to use the Calphad approach for oxidation in metallic materials too.[15,16,17,18,19,20,21,22,23]

Detailed simulations of diffusion in the oxide scale itself are numerically difficult. However, for many cases of practical interest, lifetime predictions can be based on the rate of alloying element depletion and the microstructural evolution below the oxide scale. Several coupled models have thus been developed that to different extents used separate but coupled simulation models for the growth of the oxide scale and the sub-scale microstructural evolution.[15,16,17,18,24,25]

Calphad-based assessments of oxide diffusivities are very scarce,[19,20,21,22,26] though compilations of diffusivity data have been made by Freer[27] and Van Orman and Crispin[28]. The situation is complex due to the fact that in practice the overall diffusion rate in oxides is often dominated by grain boundary diffusion; this topic is discussed by Atkinson.[29] The general reference for grain boundary diffusion is the book by Kaur and Gust.[30] A classification scheme relating the contributions from bulk and grain boundary diffusion, originally due to Harrison,[31] proposes that there are three different regimes: type C (valid at short times), type B (intermediate times) and type A (long times). In the present work only type A kinetics will be considered. According to Atkinson,[29] type A kinetics applies when \(\sqrt{D_{total}t}\ge d\), i.e., when the average diffusion distance is larger than the grain size d; \(D_{total}\) is an average effective diffusivity incorporating the bulk and grain boundary diffusivities.

In the present work we investigate the nature of the protective layer during oxidation with the help of the equilibrium and kinetic-based calculations and develop hypotheses for the corrosion behaviour of stainless steels and FeCrAl alloys based on the results. We also develop some strategies to model high-temperature corrosion in Fe-based alloys using the Calphad approach.

We show that based on the experimental observations, the oxide scale formed on FeCr and FeCrAl alloys in air mostly consists of chromium oxide and aluminum oxide at 600 and 900 \(^\circ \)C, respectively. Thus, kinetic-based simulations for oxidation of Al and Cr accurately represent the lifetime of the protective layers in these alloys at different temperatures.

2 Method

2.1 Experimental Observations of High-Temperature Corrosion

In this section we summarise the experimental data on the nature of the oxide scale formed on the stainless steels, FeCr and FeCrAl steels. We show that based on detailed microscopic investigations found in the literature, the oxide layer formed on these steels is mainly formed of chromia and alumina at T = 600 \(^\circ \)C and \(T > 900\, ^\circ \)C, respectively.

There is a great number of publications covering the oxidation properties in the steels, stainless steels and FeCrAl alloys as mentioned in the introduction (literature summarized in the established books.[2,3,4] The main aim of the present work is to model the initial oxidation of the FeCr(Al) alloys. It is well known that small amounts of e.g. Si (Refs. 2,3,4, 32,33,34,35) or reactive elements (RE), e.g. Zr, Y, La, etc., (Refs. 36,37,38,39) may influence the oxidation but these elements will not be discussed in this work. The proposed modelling approach here is generic, but the main attention will be paid to two temperatures, i.e. 600 \(^\circ \)C and 900 \(^\circ \)C. The focus of this section is to give an overview, and especially examples, of the expected microstructure/composition of the scales formed in dry oxygen to form a foundation for the calculations and predictions.

The experimental investigations reviewed in this section show that FeCr and FeCrAl alloys can be categorized into two groups regarding the protective scale formed on the surface: Cr-rich oxide forming (chromia forming), and Al-rich oxide forming (alumina forming). Composition of the alloy and temperature determine the type of the material class. Literature related to these findings are reviewed in the following.

Detailed microstructural investigations of the very thin scales formed in dry oxygen are relatively scarce. A great amount of information about the oxide microstructure can be found in the literature where the microstructure of oxide in reference cases (dry \(O_2\)) is investigated. In many cases more corrosive environments have in addition been studied. The scale formed in \(O_2\) at 600 \(^\circ \)C, i.e. reference cases, is well documented in several studies in the high-temperature oxidation of the stainless steels[40,41,42,43,44] and high-Cr steels.[45,46,47,48,49,50]

Stainless steels are normally used in non-corrosive (mild) environments at temperatures up to about 800 \(^\circ \)C. The oxidation kinetics of stainless steels has been shown to follow approximately parabolic behaviour.[2,3,4,51,52] The formation of a Cr-rich oxide scale is known to provide this material class with its corrosion protection. In the upper range of the temperature span the oxide scale was shown to consist of almost pure chromia covered by a Mn-Cr spinel (see e.g. Ref. 53) while some iron was found in the Cr-rich spinel (about 30 at.\(\%\) cations) at lower temperatures, i.e. below 700 \(^\circ \)C (see e.g. Ref. 40). In high-Cr steels on the other hand, a thin layer of \((Fe,Cr)_2O_3\) is formed on the surface in exposure to the dry air (Fig. 1). This oxide layer grows slower than the pure \(Cr_2O_3\) formed on pure Cr.[45,46,47,48,49,50,54,55]

Fig. 1
figure 1

TEM bright field image of protective \((FeCr)_2O_3\) layer formed on the surface of a Fe-10Cr alloy after exposure to air at 600 \(^\circ \)C from Ref. 45

The oxidation of a model Fe10Cr alloy in Ref. 45 in oxygen may be used as a representative illustration. In that work the scale was investigated in detail using Transmission Electron Microscope (TEM). Convergent beam electron diffraction (CEBD) showed that the thin scale consists of corundum and Scanning Transmission Electron Microscope with Energy-Dispersive X-ray spectroscopy (STEM-EDX) showed that the scale consists of an inner Cr-rich part and an outer Fe-rich part. This is in agreement with several studies of stainless steels[42,43,44,51,52] and Fe-Cr steels[45,46,47,48,49,50] where only \((FeCr)_2O_3\) is observed at 500-900 \(^\circ \)C with a Cr-rich layer in the metal-oxide interface. These alloys thus are considered chromia forming.

FeCrAl alloys may form a very protective \(\alpha \)-alumina scale at higher temperatures (900–1300 \(^\circ \)C).[2,38,56] In this temperature range the continuous scale of \(\alpha \)-alumina acts as a strong diffusion barrier which protects the underlying metal from fast oxidation. The oxidation kinetics has been shown to follow approximately parabolic behaviour.[2,3,4] FeAl alloys may also form \(\alpha \)-alumina; addition of Cr has been shown to promote the formation of \(\alpha \)-alumina, enabling FeCrAl alloys to form this protective layer at lower Al-contents. The mechanism behind this effect is known as the third element effect, i.e. the addition of a third element of a tertiary alloy system which can form an oxide with a stability intermediate to those of the two other elements.[4]

The formation of \(\alpha \)-alumina occurs at such a low rate that transient alumina is usually observed instead at medium/high temperatures (500-700 \(^\circ \)C).[34,42,57,37,38,39,40] The transient alumina layer has been shown to be protective and limits the propagation of the oxidation process.[34,42,57,58] However, the oxide microstructure of the thin protective scales is not well characterized. Reference 59 showed presence of small amount of Al in the protective scale formed on the surface of the FeCrAl alloys at 600 \(^\circ \)C. The kinetics of the transient alumina or a Cr-rich oxide has properties similar to the scale found on stainless steels.[42,43,44] It should be mentioned that in FeCrAl alloys the oxide growth at lower temperatures is very slow and thus, there is little work related to this in the literature. This class of materials is alumina forming.

These experimental observations indicate that these two material classes may be divided into Cr-rich oxide forming (chromia forming), e.g. a Fe10Cr alloy at 600 \(^\circ \)C where the Cr oxidation is the rate-determining reaction in the oxide growth, and Al-rich oxide forming (alumina forming), e.g. a Fe10Cr6Al alloys at \(T > 900 ^\circ \)C in which the Al oxidation is the representative reaction. Addition of a third element to the bulk material may affect these reactions and decrease the rate of oxidation to be slower than the one for the “pure Cr/Al” oxidation. These results are used to develop a strategy for kinetics-based simulations for these alloys in which the simulations for oxidation of Cr and Al can accurately predict the oxidation rate of these alloys at T = 600 \(^\circ \)C and \(T > 900 ^\circ \)C, respectively.

In the next section, first a short summary about the Calphad method and different models available in the Thermo-Calc software[60] are reviewed and then our new strategy will be explained. In the results section this strategy is tested for FeCr and FeCrAl alloys while stainless steels will be investigated and published using similar methodology in a future work.

2.2 Modelling High Temperature Oxidation Using Calphad Approach

The basis of a Calphad-based software is finding the stable phases at equilibrium under different conditions by minimizing the Gibbs energy of the system. The aim of the equilibrium calculations in the present work using Thermo-Calc software (section 4.1) is therefore to find the stable phases formed during the oxidation process, where the results help to interpret the microstructure and the mechanisms of the process. Although such results can be very useful, it is very unlikely that in real-life processes equilibrium is reached, especially at low temperatures. In such cases kinetics of different reactions plays an important role. To predict the actual growth rate the diffusion through the protective scale must be simulated. The models available in Dictra, a module in Thermo-Calc software, focus on the diffusion-based models by assuming local equilibrium in the system by phase evolution with time. This enables us to model the oxidation process considering both kinetics and thermodynamics of the system.

In this section first a short summary of the diffusion theory is reviewed. Then different models available in Thermo-Calc software for modelling high-temperature oxidation are discussed. Finally, our approach for life prediction of these alloys using Calphad-type models is presented.

2.2.1 Theoretical Background: Theory of Diffusion

In general, diffusion is described mathematically by the Ficks-Onsager law in one dimension (along an x-axis)[61] as Eq 1:

$$\begin{aligned} J_k=-\sum D_{kj}\frac{\partial c_j}{\partial x} \end{aligned}$$
(1)

where \(J_k\), the flux of species k, shows diffusional transport of this species by the amount that crosses a unit area perpendicular to the direction of diffusion per time unit. \(D_{kj}\) is the diffusion coefficient of k with respect to the concentration gradient of j and is an indication of how mobile this species is. This coefficient is strongly temperature-dependent and usually varies with composition.

In the absence of chemical gradients the so-called tracer diffusivity \(D^*_k\) typically obeys an Arrhenius model, shown in Eq 2:

$$\begin{aligned} D^*_k=D_k^0exp\frac{-Q_k}{RT} \end{aligned}$$
(2)

where \(Q_k\) is the activation energy and \(D_k^0\) is the frequency factor with the same dimension as the diffusion coefficient. R is the gas constant and T is the temperature in kelvin. \(D_{kj}\) can be evaluated from Eq 3:

$$\begin{aligned} D_{kj}=c_kM_k\frac{\partial \mu _k}{\partial c_j} \end{aligned}$$
(3)

where \( M_k=\frac{D^*_k}{RT}\) is the mobility of the species k.

In the case of oxidation, as the contribution from grain boundary diffusion plays an important role in the growth of oxide layer,[29,30] a rather simple model can be used for oxidation in simple systems. It can be assumed that the width of grain boundaries \(\delta \) along with the grain size (d) yield an estimate of the area fraction \(f^{gb}\) of grain boundaries perpendicular to the direction of diffusion, shown in Eq 4.

$$\begin{aligned} f^{gb}=\frac{\delta }{d} \end{aligned}$$
(4)

Furthermore, the activation energy for diffusion along the grain boundaries is assumed to be a specific fraction \(\alpha \) of the bulk activation energy and same pre-exponential factor is also assumed for bulk and grain boundary diffusion (Eq 5 and 6). Using these assumptions, the overall effective diffusivities can be obtained as a linear combination of the bulk and grain boundary diffusivities weighted by the area fractions[30] as in Eq 7.

In Equations 56 and 7, \(M_{bulk}\), \(M_{gb}\) and \(M_{total}\) represent the mobility of species in the bulk and in the grain boundary and the overall effective diffusivity respectively.

$$\begin{aligned} M_{bulk}=M_0^{bulk}\times exp(\frac{-Q}{RT}) \end{aligned}$$
(5)
$$\begin{aligned} M_{gb}=M_0^{bulk}\times exp(\frac{-\alpha Q}{RT}) \end{aligned}$$
(6)
$$\begin{aligned} M_{total}=f^{gb} \times M^{gb}+(1-f^{gb})M_{bulk} \end{aligned}$$
(7)

where \(\alpha \) is a factor to take into account a lower activation energy for diffusion along grain boundaries compared to in the bulk.[23] This coefficient basically should have a value \(<1\) and it can differ for different oxides.

2.2.2 Available Models in Calphad-Based Software for High-Temperature Oxidation

Two different models are developed in the Dictra module in Thermo-Calc for modelling oxidation in multi-component systems based on the theory explained in the previous section: the Dictra moving phase boundary model and the so-called homogenization model.

The moving phase boundary model[60] is based on the assumption of local equilibrium at a sharp (zero-width) migrating interface. It is implemented in the Dictra module based on the Finite Element Method (FEM). A schematic of a concentration profile across a moving phase boundary is shown in Fig. 2(a). It may be applied to cases where multiple oxide phases form in a predefined sequence as successive layers. Nevertheless, this model cannot handle possible miscibility gaps in the diffusion path, nor cases where the microstructural evolution involves both internal oxidation (dispersed phases) and scale formation.

To overcome the limitations mentioned above, the homogenization model[62] was developed in which the calculation domain is divided into cells/slices and full local equilibration is assumed in each cell (Fig. 2b). This model obtains the local properties from equilibrium calculations in each cell, e.g. the local chemical potentials (the same in all phases present locally since local equilibration is assumed), phase fractions, phase compositions, etc. The local chemical potential gradients provide the driving forces for diffusion in the domain. This model obtains the local permeability (essentially the product of effective solubility and diffusivity) from a mathematical coarse graining procedure. There may thus be extreme differences in permeability between neighbouring grid points cells where phase transitions occur.

The model does not track explicitly the movement of migrating interfaces. To avoid numerical instability, it is essential to use a low order finite difference approximation of derivatives which is implemented using the Finite Volume Method (FVM). This model is computationally less efficient than the moving phase boundary model although the low order numerical method adopted makes the model robust. On the other hand, it can handle all the types of simulations of interest including different phase transformations, e.g., computational difficulties derived from the thermodynamics of the systems in question. A pre-knowledge about the sequence of different phase regions occurring along the diffusion path is not required using this model. The homogenization model, as the moving phase boundary model, allows for fully implicit time-stepping.

As Dictra simulations are computationally very demanding, in the work by Larsson et al.[23] a development version of Dictra was used in which the program was parallelized using the MPI protocol. Robustness of thermodynamic equilibrium calculations was also improved during the present work by modifying the routines that generate starting values for the optimization problem that constitute the equilibrium calculation.

Fig. 2
figure 2

Schematic of different models available in Dictra program for simulation of oxidation

2.2.3 A New Approach for Life Prediction of Fe-Based Alloys Using Calphad-Type Simulations

In the present work, the Dictra version from the work by Larsson et al.,[23] was used to simulate the growth of different oxide scales on the metal surface using both the moving phase boundary[60] and the homogenization model.[62] The thermodynamic and kinetic databases used in the simulations are development versions along with the kinetic database MOBFE4 and the thermodynamic database TCFE9 and kinetic data from Refs. 22, 63, and 64.

Although the work by Larsson et al.[23] was a great improvement, for the higher-order systems, e.g. ternaries, the simulations are not computationally efficient as the number of equilibrium calculations increases as more species are involved. Thus, we aim to use experimental observations and simplified theoretical considerations to develop a new strategy to overcome these difficulties for oxidation of alloys.

As mentioned in section 2.1, experimental observations show formation of a protective layer of \((Fe,Cr)_2O_3\) on the surface of the high-Cr steels[56,65,66,67] and stainless steels[41,42,43,44] at 500-900 \(^\circ \)C. This oxide contains a region that is rich in chromium close to the metal-oxide interface. Asokan et al.,[59] showed formation of an oxide with similar properties on the FeCrAl alloys at 600 \(^\circ \)C. However, Refs. 36,37,38 showed that protective scale formed on the FeCrAl alloys at 900 \(^\circ \)C and higher temperatures mainly consists of \(Al_2O_3\). They also showed that in these alloys the protective \(\alpha \)-alumina (with corundum structure) cannot form below approximately 900 \(^\circ \)C. Thus, one can assume that Cr oxidation provides the oxidation resistance properties in these alloys at 600 \(^\circ \)C. At temperatures higher than approximately 900 \(^\circ \)C Al oxidation may instead predict a reliable lifetime of FeCrAl alloys.

We used these findings to develop a tentative strategy to improve the efficiency of the calculations (mentioned above) and develop a simplified model for the effect of alloying elements on the oxidation rate. It is proposed that simulations for oxidation of the main element in the corundum, i.e. Al and Cr, may therefore be a reliable representation of the growth rate of the alloys at different temperatures for these three categories of steels, i.e. high-Cr, FeCrAl and stainless steels.[36,37,41,42,43,44,59] In the present work we use growth of Cr oxide at 600 \(^\circ \)C to model oxidation of chromia forming alloys (explained in section 3.1). In the alumina forming alloys on the other hand (section 3.1), Al oxide growth can model the oxidation of the alloys at T \(\ge \) 900 \(^\circ \)C. Although “pure Al” melts at 660 \(^\circ \)C, pure alumina scale formed on FeCrAl alloys does not melt up to 2000 \(^\circ \)C. Thus, it is reasonable to use the growth of alumina for these alloys at high temperatures. The simulations for the alloys at high temperatures therefore considers growth of a thin pre-formed \(Al_2O_3\) metastable oxide layer in contact with an infinite source of solid Al and a constant gas flux.

Addition of different amounts of a third element to the alloy may affect the oxidation properties compared to the “pure Cr” and “pure Al”[1, 9,68,69] in two ways: the grain size of the oxide and the activation energy needed for diffusion along the grain boundary in the oxide phase since the grain boundary diffusion is dominating in this temperature regime. The activation energy for the diffusion in the grain boundary of the oxide is assumed to relate to the bulk activation energy by a pre-factor (\(\alpha \) in Equ. 6). Thus, the variation in the bulk composition should in this case show a systematic trend for the values of \(\alpha \) for different compositions. To model the composition dependence of oxidation properties of the alloys we used the grain size of the oxide for each alloy extracted from the experimental results as an input in the simulations while \(\alpha \) was treated as a fitting parameter to reproduce the experimental oxide thickness of each alloy. It is shown in the Results section that different classes of steels, i.e. alumina and chromia forming alloys, show a systematic variation of the pre-factor \(\alpha \) with Cr+Al content. In the section 4 we use the results to suggest different strategies for a generic model based on the grain boundary diffusion for alloys.

As explained above, the grain size in the oxide was taken from the experimental results for different alloys[1,45]. In many cases the oxide is very thin, i.e. in the order of 10-100 nm, and it is impossible to measure the grain size experimentally. In such cases the grain size was estimated based on some assumptions drawn from the experimental observations for the alloys that form thicker oxide scales and oxidation of “pure Cr”, i.e. where the grain size can be measured.[1,20] A schematic representation of a cross section of the oxide scale microstructure and the oxide grains is shown in Fig. 3. The assumptions based on the investigated data are: i) the grains do not grow significantly with time, ii) a very thin oxide scale (100 nm) contains no more than 2-4 grains in the normal direction to the surface, iii) the grains’ shape can be simplified as a rectangle normal to the surface and iv) the length of the rectangular grain in the a direction is half of the length in b direction.

Fig. 3
figure 3

Schematic of the grains in the cross section of the oxide scale formed of \(Cr_2O_3\). The markings show the intersection with the grain boundaries in the normal direction to the surface, i.e. the number of the grains in that direction. It is also shown how the grain’s shape can be simplified to a rectangle

The assumptions explained above for the grain size estimation were based on the results of detailed microstructural investigations of the oxide scale formed on “pure Cr” at 600 and 700 \(^\circ \)C using TEM for exposure time of 24 h[20] and 480 h.[70]

In the case of oxide growth on “pure Cr” a thin layer containing fine grains can be observed in the bottom of the oxide,[70] i.e. the oxide-metal interface, which does not show any indication of grain growth at 600 or 700 \(^\circ \)C. The growth of this layer will result in parabolic growth kinetics in good agreement with the reported mass gains.[20] However, the work by Iyer et al.[70] showed some larger grains at 700 \(^\circ \)C in the outer part of the scale after longer exposures while the small grains remains intact in the inner part. This is expected to result in a small deviation from parabolic behaviour at higher temperature and longer exposure times, i.e. the model might overestimate the oxide growth under these conditions. While this might be a challenge for the higher temperatures and longer exposure times, it will not be a concern for the FeCr and FeCrAl alloys at 600 \(^\circ \)C as no such grain growth is observed in the chromia layer formed on these alloys[1] or on “pure Cr” for short time exposures and 600 \(^\circ \)C[20].

3 Results

Iron-based alloys consist of a very large number of grades and compositions. However, for most of the alloy classes studied in the present work, high-temperature corrosion applications and the corrosion protection the key elements are Fe, Cr and Al. In this section we first show the results from equilibrium calculations for the Fe-Cr-Al-O system. In the second part kinetics is taken into account to predict the oxide’s growth rates.

3.1 Equilibrium Calculations

Figure 4 shows equilibrium calculations for the quaternary Fe-Cr-Al-O system. In this figure the partial pressure of oxygen is fixed to \(P_{O_2}=5\times 10^{-20}\) atm, i.e. within the stability region of corundum. Similar results (with slight changes in the phase boundaries) are achieved at higher partial pressures of oxygen, i.e. \(P_{O_2}=5\times 10^{-5}\) atm. Since oxygen activity is constant, the system has the same topology as a ternary and an isothermal section can be calculated, shown in Fig. 4. Phase diagrams for the sub-systems at the same oxygen partial pressure, i.e. Fe-Cr-O, Fe-Al-O and Cr-Al-O, are used for the ternary construction in Fig. 4. The temperatures 600 and 900 \(^\circ \)C are marked in all sub-systems by red and blue dashed lines, respectively.

It is shown in Fig. 4 that the phase diagrams of all the three sub–systems, i.e. Fe-Al-O, Fe-Cr-O and Cr-Al-O, show a miscibility gap in the corundum phase. The calculations are performed for the \(M_2O_3\) (corundum) structure as this is known to provide the protection at high temperatures[1]. The range of decomposition temperatures is very different in different sub–systems. The decomposition can be observed in the Fe-Al-O and Cr-Al-O systems at 600 and 900 \(^\circ \)C (marked lines in Fig. 4) while lower temperatures are required for the decomposition in the Fe-Cr-O system. The marked temperatures are those of interest used for kinetic simulations.

All phase diagrams shown in Fig. 4 are calculated using database TCOX9 in Thermo-Calc software, except the Al-Fe-O system which is based on the work by Lindwall et al.[71]. In the case of Fe-Cr-O system, a T-term with opposite sign of TCOX9 description was added to the interaction parameters used. As TCOX is a commercial encrypted databases, the details of the description and parameters cannot be discussed here. Nevertheless, the aim of this adjustment was to improve the description of corundum at lower temperatures in this system, according to the experimental observation work by Murakami et al.,[72]. The small kink in the middle of the system corresponds to the magnetic transition.

Fig. 4
figure 4

Isothermal section for Fe-Cr-Al-O system at 600 \(^\circ \)C and \(P_{O_2}=5\times 10^{-20}\) atm showing a miscibility gap occurring in corundum phase. Same phenomenon is found in all sub-systems (shown as triangle construction). Compositions are in molar fraction and temperature 600 and 900 \(^\circ \)C are marked in all sub-systems by red and blue dashed lines respectively

Experimental studies by Somiya and Muan[73] showed the existence of a miscibility gap in the \(Fe_2O_3-Cr_2O_3-Al_2O_3\) system at the temperature range 1300-1500 \(^\circ \)C and high partial pressure of oxygen. There is no experimental equilibrium evidence for such phase transformation at lower temperatures, i.e. 600 \(^\circ \)C. However, in all the sub-systems the results indicate a miscibility gap at this temperature and low partial pressure of oxygen as shown in Fig. 4 based on the experiments in Refs. 71, 72, and 74. Thus, it is reasonable to assume that isothermal section of the quaternary extrapolated from the description of the sub-systems presented in Fig. 4 can be trusted. For the \(Fe_2O_3-Cr_2O_3\) system however, no decomposition occurs at 600 \(^\circ \)C (red dashed line does not intersect with the miscibility gap). This is in good agreement with the experimental observations for the Fe-Cr alloys exposed to air at 600 \(^\circ \)C.[56, 65,66,67]

The results presented in Fig. 4 suggest that oxidized FeCrAl alloys may form a non-stable corundum at high temperatures. If the system reaches equilibrium there will be a decomposition given that the oxide scale’s chemical composition is within the miscibility gap. In that case two types of corundum may form within the protective scale i.e. one Cr-rich and the other Al-rich. In the case when the formation of an Al-rich corundum is promoted, i.e. in alloys containing high amount of Al, it may lead to a more protective oxide layer.

Although the calculation results show that corundum will decompose at 600 \(^\circ \)C, the experimental observations suggest that the size of decomposed areas in the microstructure may not be large enough to play an effective role unless only the grain boundaries are affected. As the diffusion of Al ions is very slow at this temperature, the decomposed microstructure may take a very long time to develop. This is confirmed by the experimental results reported from Asokan et al.,[59] whose specimen (Fe10Cr4Al) had been exposed for 168 h at 600 \(^\circ \)C in air. They investigated the oxide in detail by TEM and only the very initial stage of decomposition could be observed. It may be concluded that the kinetics of the system does not act in favour of the formation of the Al-rich corundum at 600 \(^\circ \)C for alloys containing low/medium amount of Al. For the alloys containing high amount of Al (\(\ge 4 \%\)wt) on the other hand, even initiation of the decomposition can affect the corundum structure and form a very strong protective scale even though the equilibrium is not completely reached. Experimental investigations of the FeCrAl alloys containing 4-6 \(\%\)wt Al have shown this alloys resists breakaway of the oxide scale even in a very corrosive environment after a long time exposure.[1,59]

At 900 \(^\circ \)C and higher temperatures on the other hand, it is more likely that the equilibrium is reached as the kinetic of the system is much faster. References 36,37,38,39 showed clear regions of \(Al_2O_3\) in the corundum structure. It is shown in these works that the protective layer is formed of pure \(Al_2O_3\), thus oxidation of Al can model the oxidation properties accurately. While formation of alumina can be considered to relate to the presence of traces of reactive elements in Refs. 36,37,38,39, our experimental exposure for a FeCrAl model alloy at 900 \(^\circ \)C shows formation of a pure \(Al_2O_3\) scale too[63].

These results can be used as a basis for the life prediction model suggested in section 2, in which oxidation of Cr may be used to model oxidation of chromia forming alloys at 600 \(^\circ \)C and at 900 \(^\circ \)C and higher temperatures oxidation of Al can predict the oxide growth on alumina forming alloys. In the next section the diffusion-based simulation results are presented and discussed to test this new methodology.

3.2 Dictra Simulation Results

3.3 Oxidation at 600 \(^\circ \)C

In this section we aim to test and evaluate our simulation strategy described in the previous sections. We use the oxidation of Cr to predict the corrosion process of FeCr and FeCrAl alloys at \(T=600 ^\circ \)C. Simulation results are compared with oxide thickness values calculated from mass gain data in Ref. 1 for different model alloys.

It is well known that an oxide’s growth rate may change as a result of adding an element to the bulk. The effect is known as “the third element effect” explained by Wagner.[9] In addition, Kofstad[2] discussed the phenomenon in Fe-Cr-Co and Fe-Cr-Ni alloys. The additional element affects the grain boundary diffusion and consequently the oxide’s growth rate. It is shown in the present work that Al can play a similar role in the synergic effect with Cr in FeCrAl alloys.

Figure 5 shows Dictra simulation results for the Cr-O system after 24 h at 600 \(^\circ \)C in comparison with experimental data[1] for oxidation of FeCr and FeCrAl alloys. As mentioned in section 2.2.3, bulk composition of each alloy is assumed to affect only the grain size and grain boundary diffusion in the oxide. Grain size for each alloy was taken from the experiment[1] while the activation energy coefficients \(\alpha \) (Equ. 6) was treated as a fitting parameter in the simulations to find the best agreement with the experimental data.

Fig. 5
figure 5

Dictra simulation results for the Cr-O system after 24 h at 600 \(^\circ \)C fitted to the experimental results for FeCr and FeCrAl alloys.[1] Different curves show results for variation of the activation energy coefficients (\(\alpha \) in Eq 6), treated as a fitting parameter in the simulations. Values for the oxide grain size are taken from the experiment[1] for different model alloys

Table 1 shows the obtained values for fitting \(\alpha \) versus experimental oxide thickness for “pure Cr”, different FeCr and FeCrAl alloys in Fig. 5. These results are plotted in Fig. 6 with variation of Cr/Cr+Al content for each alloy as x axis.

Table 1 Obtained values for \(\alpha \) in the present work through fitting versus experimental results for oxide thickness of “pure Cr”, different FeCr and FeCrAl alloys, see Fig. 4

Figure 6 indicates that the fitted values for \(\alpha \) seem to indicate a systematic trend with the variation of Cr/Cr+Al content; while for FeCr alloys it has a larger value than for “pure Cr”, the FeCrAl alloys with low/medium Al content have \(\alpha \) closer to the “pure Cr”. The FeCrAl alloys with high Al content show significantly larger values for \(\alpha \) compared to “pure Cr”, suggesting the possibility of formation of Al-rich corundum and a decomposed microstructure, explained in the previous section.

As increasing the Cr/Cr+Al content in low/medium Al alloys increases the activation energy for grain boundary diffusion, it can be speculated that Al affects the grain boundary in the oxide in connection to the “third element effect” explained in the beginning of this section. The highest activation energy corresponds to the Fe25Cr3Al alloy which shows the lowest oxide thickness among FeCrAl alloys with low/medium Al content. This shows that we can accurately model dependence of the grain boundary diffusion in the oxide to the bulk composition, i.e. dependence of \(\alpha \) to Cr/Cr+Al, and use this model to predict the oxidation behaviour in these alloys. If more experimental data are available, a mathematical function can be fitted to the data points and a generic model can be developed for composition dependence of oxidation properties in the alloys based on the oxidation of Cr. This is summarized in the section 4.

Fig. 6
figure 6

Trend of obtained values for \(\alpha \) in the present work for “pure Cr”, different FeCr and FeCrAl alloys when Cr/Cr+Al varies

Fig. 7
figure 7

Dictra simulation results for oxidation of FeCr alloys for longer time at 650 \(^\circ \)C (left) and 600 \(^\circ \)C (right)

Figure 7 shows simulation results for longer time, when the \(\alpha \) parameter is fitted to the first 24 h of oxidation. In both temperatures (600-650 \(^\circ \)C) the model overestimates the oxide thickness for longer time compared to the experimental data. However, from a life prediction point of view this is acceptable, especially at 600 \(^\circ \)C where the oxide thickness is very thin ,i.e. in the order of 100 nm. At 650 \(^\circ \)C, we were able to use the same \(\alpha \) parameter (black curve) for life prediction of the three alloys containing 20-38 \(\%\)wt Cr after 40 h, showing that the model is generic for different chromia forming alloys. We could not find any long time exposure for FeCrAl alloys with medium/low amount of Cr in the literature at 600-700 \(^\circ \)C to test the predictivity of the model for these alloys. This work can be a motivation for such experiment in the future work.

3.4 Oxidation at 900 \(^\circ \)C and Above

The results in the previous section indicate that Cr oxidation can be a safe method for life prediction of FeCr and FeCrAl alloys at 600-650 \(^\circ \)C. As explained in section 4.1 only pure \(Al_2O_3\) is observed on the surface of FeCrAl alloys at higher temperatures; thus, in our model we used Dictra simulations for metastable solid Al at 900 \(^\circ \)C to reproduce the experimental values for the oxide thickness of these alloys. There is limited experimental work on oxidation of pure Al at 900 \(^\circ \)C in the literature and in those that investigated, higher temperatures (summarized in Ref. 19) the exposed samples actually contained liquid aluminum covered by coarse and dense solid \(\alpha -Al_2O_3\) scale as melting point of pure Al is 660 \(^\circ \)C. Liquid state can effectively change the kinetics of the reaction. This can clearly be seen in Fig. 8 where oxide thickness for “pure Al” from work by Coker[75] is compared with growth of alumina in the FeCrAl alloys.[36,37,63,76]

The size of the grains of \(Al_2O_3\) is measured in Ref. 19 yielding an average of 4 \(\mu \)m at 1150-1550 \(^\circ \)C after short time exposures, e.g. 1–12 h. Based on these results we estimated the grain size at 900 \(^\circ \)C to be 0.5-1 \(\mu \)m after 24 h exposure. Figure 8 shows Dictra simulations results, compared with experimental data for a Fe5Al23Cr alloy[76], Fe5Al21Cr alloy (with trace of reactive elements) from Refs. 36, 37 and a Fe6Al18Cr model alloy[63].

The solid curve shows the results from a simulation for Al-O with \(\alpha \) (activation energy for grain boundary diffusion) equal to 0.524 and grain size of 0.5 \(\mu \)m. Using these values the oxide thickness after first 24 h exposure for the Fe5Al23Cr alloy from[76] was reproduced by treating \(\alpha \) as a fitting parameter. It should be noted that Huntz et al.’s data[76] are for three different Fe5Al23Cr alloys: high purity alloy, Y-doped and C+Y doped alloys. All sets of their data are shown in Fig. 8 while only the oxide thickness for the high purity alloy after 24 h is used for fitting. The fitted value for \(\alpha \) (0.524) was used for longer time exposures shown in Fig. 8. It can be seen that although using this parameter overestimates the oxide thickness for all the alloys after 168 h, the results are reliable from a life prediction point of view.

Fig. 8
figure 8

Dictra simulation results for oxidation of Al at 900 \(^\circ \)C compared with experimental data for different FeCrAl alloys[36,37,63,76], and “pure Al”[75]

4 Discussion: Strategies for High-Temperature Corrosion Simulations of Fe-Based Alloys

Our results show the potential of the equilibrium calculations to predict the phases in high-temperature oxidation and using the Dictra simulations to predict the kinetics of oxidation. Even though the results are for the binary systems, the predictions are relevant for higher-order systems too. In this section we suggest two strategies to improve and generalize modelling of the oxidation in alloys. In the first one, the oxide’s grain boundary composition and kinetics are assumed to depend on the bulk properties, while in the second one the oxide’s grain boundary is suggested to be modelled in an independent way.

4.1 Modelling the Dependence of the Activation Energy of the Oxide’s Grain Boundary Diffusion on the Bulk Composition

For many cases of practical interest the growth rate of oxide scales is dominated by grain boundary diffusion. If so-called type A kinetics[30] is assumed, the contribution from grain boundary diffusion can be taken into account using Eq 7.

According to Klotsman et al.[77], cited by Stechauner and Kozeschnik[78], \(\alpha \) in Eq 6 is typically in the range of \(0.3-0.4\).

In the previous section we showed that by fitting this parameter for a pure element oxide (\(Cr_2O_3\) and \(Al_2O_3\)) good agreement can be obtained with the oxidation kinetics of the alloys. We also showed a dependence of the \(\alpha \) parameter on the bulk composition and by treating \(\alpha \) as a fitting parameter for reproducing the experimental results for the oxide thickness, i.e. mass gain, a generic model can be developed that relates the alloy composition to the oxide growth rate. A possible trend found in the fitted values in relation to the bulk’s chemical composition may be explained by doping the oxide and thereby changing the diffusion properties. However, this needs to be investigated experimentally. The results may be used in any case in the future to formulate a mathematical expression for the possible trend found in the \(\alpha \) values in relation to the bulk composition.

4.2 A Model to Handle the Oxide’s Grain Boundary Diffusion Using an Extra Sublattice Model

Compound Energy Formalism (CEF)[79] is used to model the solution phases in the Calphad databases. In CEF, each phase is modeled by one or more sublattices based on the structure of the atoms. For example in a phase represented by the formula (A,B)\(_k\)(D,E,F)\(_l\), the first parenthesis show mixture of A and B atoms in one sublattice and the second parenthesis show mixture of D, E and F in the other one. The coefficients k and l are the stoichiometric coefficients and one mole of formula units thus contains \(k+l\) moles of atoms. For the interstitial species, defects in the form of vacancies (VA) are added to the sublattice that contains them, e.g. (A,B)\(_k\)(C,VA)\(_l\) where C is an interstitial species.

Diffusion in oxides is intrinsically more complex than diffusion in metallic solid solutions and typically exhibits a very strong dependence on the defect structure. To tackle this, Hallström et al.[21] chose to model diffusion separately on each sublattice such that the total effective bulk mobility is obtained by summing over all sublattices. The general expression for the mobility of element k used by Hallström et al. can be found in Ref. 23 (see also Ref. 26) shown in Equ. 8:

$$\begin{aligned} M_k=\frac{\sum _s a_s\sum _j b_{kj}y_j^s z_{Va}^sM_j^s}{\sum _s a_s\sum _j b_{kj}y_j^s} \end{aligned}$$
(8)

where

  • \(M_j^s\) is the mobility of species j on sublattice s

  • \(a_s\) is the number of sites on sublattice s

  • \(b_{kj}\) is the number of moles of element k per mole formula unit of species j

  • \(y_j^s\) is the site fraction of species j on sublattice s

  • \(z_{Va}^s\) is identical to \(y_{Va}^s\) if vacancy is a defined species on sublattice s, otherwise it is set to one.

This approach has been quite successful to model diffusion in oxides.[21,23,26] A more qualified approach is to make separate assessments of the grain boundary diffusivities rather than rely on the bulk diffusivities and the empirical parameter \(\alpha \), though this is difficult due to the scarcity of data. Stechauner and Kozeschnik[78] performed such assessments for grain boundary and pipe self diffusion in Al, Fe and Ni. Segregation to the grain boundary may make the situation even more complex. A case in point is carbon diffusion through oxide scales, see e.g. Ref. 25, where the carbon does not even dissolve in the lattice.

The discussion above suggests that another way to handle grain boundary diffusion would be to introduce a “grain boundary sublattice”. The sublattice description for a phase with constituents \(A,B,C,\dotsc \) that normally has two sublattices may then be as:

$$\begin{aligned} \left( A,B,\dotsc \right) _a\left( C,\dotsc \right) _b\left( A,B,C,\dotsc \right) _\varepsilon \end{aligned}$$

where a and b are the number of the sites on the first and second sublattice, respectively. The third sublattice would be the grain boundary sublattice where the number of sites \(\varepsilon \propto \delta /d\). Such a treatment could potentially handle both grain boundary segregation and grain boundary diffusion and does not require any change in formalism or computer code. The grain boundary is sometimes pragmatically handled as a separate phase, but from a thermodynamic (and computational) point of view this is less attractive since the prescribed presence of such a “phase” reduces the degrees of freedom of the system.

5 Conclusion

In the present work we suggest new strategies to model high-temperature corrosion of Fe-based alloys using the Calphad approach. It is shown that kinetic-based simulations for oxidation of Al and Cr can be an accurate representative for life prediction of protective layers in these alloys at different temperatures. In addition, some hypotheses are developed for different corrosion mechanisms of FeCr and FeCrAl alloys based on the equilibrium calculations and kinetic simulation results. It is suggested that addition of Al to the bulk affects the activation energy for grain boundary diffusion in the oxide.

We also suggest two new methods to model the composition dependence of the grain boundary diffusion in the oxide which dominates the growth rate of oxide scales in the case of Fe-based alloys.