Modeling Inclusion Formation during Solidification of Steel: A Review

Abstract

The formation of nonmetallic inclusions in the solidification process can essentially influence the properties of steels. Computational simulation provides an effective and valuable method to study the process due to the difficulty of online investigation. This paper reviews the modeling work of inclusion formation during the solidification of steel. Microsegregation and inclusion formation thermodynamics and kinetics are first introduced, which are the fundamentals to simulate the phenomenon in the solidification process. Next, the thermodynamic and kinetic models coupled with microsegregation dedicated to inclusion formation are briefly described and summarized before the development and future expectations are discussed.

1. Introduction

Nonmetallic inclusions are generally considered to be detrimental to the properties of steels such as ductility, fatigue, strength, and corrosion. Many efforts have been made in the last few decades to achieve a lower amount of nonmetallic inclusions in the steel matrix and to control their size and chemical composition through optimizing steelmaking technologies, such as tundish and protected slags. This evolution led to so-called “clean steel production” [1]. In parallel, new tools such as computational thermodynamics or more sophisticated material analysis methods have become available; consequently, knowledge of the relationships between nonmetallic inclusions and the microstructure and mechanical properties of steels has increased [2]. In the 1980s, Takamura and Mizoguchi [3,4] introduced the concept of “oxides metallurgy” in steels, illustrating that the finely dispersed oxides could act as heterogeneous nuclei for other—less harmful—precipitates and for intragranular (acicular) ferrite, which may contribute to the improved mechanical properties of steel. Considering the aspects of steel cleanness and the utilization of nonmetallic inclusions, the concept of inclusion engineering was further proposed, which is explained in Figure 1. Key objectives, on the one hand, include modifying harmful inclusions into harmless particles and, on the other hand, to produce inclusions with an adjusted composition, structure, size, and number density to optimize the microstructure [5,6,7].

In steelmaking, the first inclusion populations form during deoxidation. The high content of dissolved oxygen is precipitated as oxides by the addition of oxygen affine elements such as aluminum, manganese, or silicon. This process is well understood and the formed inclusions can be partly separated out later into the ladle slag. The control of fluid flow and slag compositions during ladle treatment are important. The residual oxide inclusions and the inclusions generated in the casting process—also being sulfides and nitrides—will mostly remain in the solid steel. These inclusions are usually small, but may also play a significant role in determining the quality of the steel products. It is therefore important to study inclusion formation during the cooling and solidification process. Note that the aforementioned inclusions are termed as endogenous inclusions according to their origin. On the other hand, exogenous inclusions result from reoxidation, the entrainment of oxide bifilms [8,9], and slag and ceramic corrosion. In the most continuous casting steels, the endogenous inclusions are dominant. Hence this work concentrates on the endogenous inclusion formation.

Microsegregation is a fundamental aspect to be considered and inclusion formation is not an exception. Microsegregation results from uneven partition in the solid and liquid steel at the dendritic scale. The further diffusion of solutes influences their distributions in the phases. The enriched concentrations can lead to the growth and transformation of pre-existing inclusions and the nucleation of new inclusions. In addition, this phenomenon results in the formation of defects during the casting process (e.g., hot tearing) and negatively affects product quality (inhomogeneous microstructure) [10,11].

The online control of inclusion formation during steel solidification is still extremely difficult. The increasing development of computer science and computational thermodynamics offers a powerful and valuable tool to simulate inclusion formation and microsegregation. At the beginning of the 1990s, Matsumiya [12] presented an overview of the mathematical analysis of chemical compositional changes of nonmetallic inclusions during the solidification of steels. The commonly applied microsegregation models and the coupled inclusion formation thermodynamic models were reviewed. Based on that work, this paper aims at summarizing the coupled models on inclusion formation during steel solidification including both thermodynamics and kinetics. First, however, the popular microsegregation models and fundamentals on inclusion formation are briefly introduced, and the recent developments and future tasks on the proposed topics are highlighted.

2. Fundamentals

When simulating the formation of inclusions during solidification, fundamental theories and sub-models are necessary. As the fundamental input, the models evaluating segregated concentrations of solutes were selectively introduced. Then, general formation thermodynamics and kinetics of the inclusions were reviewed based on former reports.

2.1. Microsegregation

Due to the importance of microsegregation, this research topic has been widely investigated. Kraft and Chang [13] summarized a variety of modeling work on microsegregation despite ongoing development. The existing modeling approaches consist of numerical and analytical solutions. Analytical models were constructed based on several simplifying assumptions, such as dendrite geometry and diffusion in matrix. Numerical models provide more freedom to physical phenomenon while they require more efforts on both computational mathematics and coding. Compared with the numerical models, the simplifications of the analytical solution also bring the limitations: fixed partition and diffusion coefficient, assumed advance velocity, simple treatment on dendrite coarsening, unified density and no consideration on solid phase transformation and microstructure. However, the coupled model of inclusion formation during solidification aims at purely calculating the concentrations of solutes in the steel matrix. Hence, a relatively simple and easy method of handling microsegregation models is preferable, while more elaborate models (e.g., two-dimensional model [14]) and software products (e.g., DICTRA^® [15] and IDS^® [16]) exist that are dedicated to complex phenomena such as microstructure evolution and phase transformation [17,18,19]. Meanwhile, analytical models were developed to overcome some of the limitations. In this section, the microsegregation models widely coupled to calculate inclusion formation are briefly described.

2.1.1. Lever Rule

The Lever Rule assumes the complete diffusion of solutes in both liquid and solid. At a specific solid fraction ($f_S$), the interfacial concentrations of solute in solid ($C_S^{\ast }$) and liquid ($C_L^{\ast }$) are equal to those in solid ($C_S$) and liquid ($C_L$) far away from the interface. A general mass balance can be given as Equation (1). The redistribution of solutes in the solid and liquid phases is described by the partition coefficient ($k=C_S^{\ast }/C_L^{\ast }$). Consequently, the concentrations in the residual liquid are obtained from Equation (2).

During solidification, the complete diffusion of solutes can hardly be reached, especially in the solid phase. The microsegregation calculated by the Lever Rule is therefore underestimated. For fast diffusion elements in steel such as carbon, the Lever Rule predictions can be close to the real situation. However, this method cannot avoid that the inclusion formations in the solidification process are decreased and postponed, or even missed.

$$C_S{f}_S+C_L{f}_L=C_0$$

(1)

$$C_L=\frac{C_0}{1+k{f}_S-f_S}$$

(2)

2.1.2. Scheil Model

A more practical model on microsegregation was proposed by Scheil [20], which was also derived by Gulliver [21]. In contrast to the Lever Rule, it assumes no diffusion in solid and successful mixing in liquid. With the interfacial equilibrium, the solute enrichments can be calculated with Equation (3), which is the differential form of the Scheil Model:

$$C_L\left(1-k\right)d{f}_S=\left(1-f_S\right)d{C}_L$$

(3)

$$C_L=C_0{\left(1-f_S\right)}^{k-1}$$

(4)

Note that, in this case, the equilibrium partition coefficient changes with the proceeding solidification. In most subsequent applications, the partition coefficient was assumed as constant for simplification and the absence of local values. Furthermore, the integrated form of the Scheil Model was obtained in Equation (4). Besides the concentrations in the residual liquid, the concentration profiles in solid are also available. Due to the lack of diffusion in the solid, the compositions of the formed solid phase remain unchanged.

Compared to the applied conditions of the Lever Rule, the Scheil Model is more appropriate for substitutional solutes with low diffusivity. In contrast, it overestimates the microsegregation for interstitial solutes such as carbon and nitrogen, which diffuse quickly in steel. The interfacial concentrations are infinite when the solid fraction approaches one, which also limits the application of the Scheil Model, because the final concentrations and solidus temperature are important expectations. To overcome the aforementioned limitations, an improved Scheil Model that considered the back diffusion was proposed [22], which is also termed as the partial equilibrium method [23]. In the partial equilibrium method, the perfect diffusion of interstitial solutes and no diffusion of substitutional solutes in solid steel are accounted. This simple, but powerful scheme was believed to be suitable to multicomponent alloys [23].

2.1.3. Brody–Flemings Model and Clyne–Kurz Model

Brody and Flemings [24] proposed a model based on the analysis of the Scheil Model [20,21]. In the model, the finite diffusion in the solid steel is accounted for by introducing back diffusion. When assuming a parabolic thickening of dendrite, Equation (5) is achieved for estimating the solute enrichments in the residual liquid. In Equations (5) and (6), $C_L$ and $C_0$ are the concentrations in the residual liquid and the initial value, respectively; $\mathsf{\alpha }$ is the so called back diffusion coefficient as given in Equation (6); and $D_S$ is the diffusion coefficient in solid. Note that the partition and diffusion coefficients are assumed as constants in the equations:

$$C_L=C_0{\left(1-\left(1-2\alpha k\right)f_S\right)}^{\frac{k-1}{1-2\alpha k}}$$

(5)

with

$$\alpha =\frac{4D_S{t}_f}{{\lambda }^2}$$

(6)

where $\alpha$ is equal to 0 and 0.5 in the Scheil Model and Lever Rule, respectively, which is found from Equation (5). It is easy to understand that there is no diffusion in the solid, as assumed in the Scheil Model. It is not reasonable to achieve the Lever Rule ($\alpha$ = 0.5); however, as back diffusion should be infinite in a well-mixed solid. As a result, this model is confined to the limited solid diffusion solutes.

To overcome the limitation of the Brody–Flemings Model, Clyne and Kurz [25] replaced the back-diffusion coefficient $\alpha$ with $\Omega$, as given in Equation (7). With this mathematical treatment, when $\alpha$ is equal to zero and infinite in the Brody-Flemings Model approaches to the Scheil Model and Lever Rule, respectively. Later the improved model (Clyne–Kurz Model) has been widely applied in microsegregation prediction:

$$\Omega =\alpha \left\{1-exp\left(-\frac{1}{\alpha }\right)\right\}-\frac{1}{2}exp\left(-\frac{1}{2\alpha }\right)$$

(7)

2.1.4. Ohnaka Model

Ohnaka [26] introduced a columnar dendrite diffusion model where one-dimensional diffusion in the triangle area is considered as an approximation for three-dimensional diffusion. The analytical solution in differential form is given in Equation (8). Similarly, by providing the constant partition coefficients and diffusion coefficients, Equation (9) is obtained by integrating Equation (8):

$$\frac{d{C}_L}{C_L}=\frac{\left(1-k\right)d{f}_S}{\{1-\left(1-\frac{\beta k}{1+\beta }\right)f_S\}}$$

(8)

where $\beta$ is equal to $2\alpha$ and $4\alpha$ for plate and columnar dendrite models, respectively; $\alpha$ is the back-diffusion coefficient given by the former models.

$$\frac{C_L}{C_0}={\left(1-\Gamma ·f_S\right)}^{\left(k-1\right)/\Gamma }, with \Gamma =1-\frac{\beta k}{1+\beta }$$

(9)

However, partition coefficients actually depend on the concentrations of various chemical components and temperature rather than being constant. Diffusion coefficients are also strongly influenced by temperature; therefore, local partition coefficients and diffusion coefficients for different compositions and temperatures are desired for predicting microsegregation in various steels. The present authors in [27] modified the differential equation (Equation (9)) to the difference equation (Equation (10)). In this way, the changes of the partition and diffusion coefficients were taken into consideration. Local partition coefficients and diffusion coefficients were calculated at each solidification step, but within the increase of solid fraction by $\mathrm{\Delta }{f}_s$, they were assumed to be constants. In the proposed model, with the help of the thermodynamic library ChemApp [28], the non-equilibrium solidification temperature was also reasonably predicted.

$$C_L^{+}=C_L{\{\frac{1-\Gamma \left(f_S\right)·f_S}{1-\Gamma \left(f_S\right)·(f_S+\mathrm{\Delta }{f}_S)}\}}^{\frac{1-k\left(f_S\right)}{\Gamma \left(f_S\right)}}, with \Gamma \left(f_S\right)=1-\frac{4\alpha \left(f_S\right)k\left(f_S\right)}{1+4\alpha \left(f_S\right)}$$

(10)

where $C_L^{+}$ and $C_L$ are the concentrations of solutes in the residual liquid at solid fractions of $f_s$ and $f_s+\mathrm{\Delta }{f}_s$, respectively; $k\left(f_S\right)$ and $\alpha \left(f_S\right)$ are the local partition coefficient and back diffusion coefficient at the solid fraction $f_S$.

2.1.5. Ueshima Model

Ueshima et al. [29] applied a finite difference method to model the solute distribution in both solid and liquid phases during the solidification of steel and a hexagonal columnar dendrite shape was assumed. The local equilibrium at the transformation interfaces existed and the redistribution of the solutes depended on the partition coefficient. Providing that the solutes only diffuse in one dimension, the model solved the diffusion equation (Equation (11)) in the triangle transvers cross section of the dendrite, which was numerically discretized. Solving the achieved difference equations, the concentrations in the analyzed region were tracked during and after solidification. To calculate the inclusion formation during solidification, the Ueshima Model is especially useful when the precipitation not only in the residual liquid, but also in the solid phase need to be considered. Meanwhile, the influence of peritectic transformation is possibly accounted for:

$$\frac{\partial C}{\partial t}=D·\frac{{\partial }^{2}C}{\partial x^2}$$

(11)

where $C$ is the concentration of solutes; $t$ is time; $D$ is the diffusion coefficient in local phase; and $x$ is the diffusion distance.

The abovementioned microsegregation models offer the predicted solute concentrations and temperatures that are the primary input for simulating inclusion formation during solidification. Hence, in addition to the flexibility to be coupled, reasonable and qualified predictions are also desirable, which is the “driving force” for the continuous improvement on an easily handled model. Considering both aspects and the requirements, models were selectively applied.

2.2. Thermodynamics of Inclusion Formation

In metallurgical processes, thermodynamics are mainly concerned with the state change of a system influenced by energy motion [30]. With the help of energy difference, the possibility and extent of chemical reactions are defined. Concerning an inclusion as a new phase in a steel matrix, its stability can be evaluated using thermodynamics.

The formation reaction of simple stoichiometric inclusion is generally described using Equation (12). Here $\left[P\right]$ and $\left[Q\right]$ are the formed elements of inclusion $P_x{Q}_y$, which are dissolved in liquid steel where $x$ and $y$ are the atom numbers in the molecule. The Gibbs free energy change is the most popular thermodynamic criteria. At a given temperature, the Gibbs energy change ($\mathrm{\Delta }G$) for the reaction is given by Equation (13).

$$x\left[P\right]+y\left[Q\right]=P_x{Q}_y$$

(12)

$$\mathrm{\Delta }G=\mathrm{\Delta }{G}^0+RTln\left(\frac{a_{P_x{Q}_y}}{a_P{}^x{a}_Q{}^y}\right)$$

(13)

where $\mathrm{\Delta }{G}^0$ is the standard Gibbs energy change, which is a function of temperature; $R$ is the gas constant; and $a_i$ is the activity of species $i$. For formation of common inclusions, the empirical expressions of standard Gibbs energy change in liquid iron are available.

When

• $\mathrm{\Delta }G<0$, the reaction can happen in the right direction and the inclusion is stable.
• $\mathrm{\Delta }G>0$, the reaction proceeds towards the left and means that the inclusion $P_x{Q}_y$ will not precipitate.
• $\mathrm{\Delta }G=0$, the reaction reaches the equilibrium state, where Equation (14) is achieved.

  $$\mathrm{\Delta }{G}^0=-RTln\left(\frac{a_{P_x{Q}_y}{}^{eq}}{{\left(a_P{}^{eq}\right)}^x{\left(a_Q{}^{eq}\right)}^y}\right)$$

  (14)

  $$\mathrm{\Delta }G=RTln\left[\frac{{\left(a_P{}^{eq}\right)}^x{\left(a_Q{}^{eq}\right)}^y}{a_P{}^x{a}_Q{}^y}\right]\approx -RTln\left(\frac{K}{K^{eq}}\right)$$

  (15)

where the superscript $eq$ means equilibrium.

When assuming that inclusion $P_x{Q}_y$ is a pure solid phase, its activity is equal to one ($a_{P_x{Q}_y}$ = $a_{P_x{Q}_y}{}^{eq}$ = 1). Consequently, Equation (13) can be written into Equation (15). In the dilute solution, the Gibbs energy change is also estimated by the ratio of concentration product ($K=C_P{}^x{C}_Q{}^y$) and solubility product ($K^{eq}={\left(C_P{}^{eq}\right)}^x{\left(C_Q{}^{eq}\right)}^y$), which is termed supersaturation ($S_{sat}$). This means that when the supersaturation is larger than one, the inclusion will be stable. Simultaneously, thermodynamics decide the chemical driving force for the inclusion formation as displayed in Figure 2. It was also found that supersaturation promoted the proceeding of the reaction until it reached the equilibrium state, while it gradually decreased due to the consumption of solutes. Correspondingly, the absolute value of free energy change approached zero. To some extent, this driving force was the link between thermodynamics and kinetics. The detailed application of the chemical driving force to inclusion nucleation and growth is discussed later.

2.3. Kinetics of Inclusion Formation

On the basis of thermodynamics, kinetics defines the rate of the chemical reaction. Specific to inclusion formation, the evolution of size and number density are described using kinetics. In this way, the size distribution of inclusions can be studied and controlled. Furthermore, inclusion composition and amount are simultaneously achievable.

2.3.1. Nucleation

Inclusion can homogeneously nucleate in the melt or heterogeneously on the existing matrix, which are accordingly termed as homogeneous nucleation and heterogeneous nucleation. Classical nucleation theory [31,32,33,34,35] is widely used and illustrated as valid to investigate precipitation related topics. The development and detailed description of classical nucleation can be referred to in [36].

Homogeneous Nucleation

When assuming a spherical nucleus with a radius of $r$ generating, the free energy change of the system is given in Equation (16). In Equation (16), the first term describes the Gibbs energy change caused by the chemical reaction of nucleus formation. $\mathrm{\Delta }{G}_V$ is the volume energy change of inclusion formation, which is calculated by the ratio of molar Gibbs energy change and the molar volume of inclusion. The second term is the energy obstacle resulting from the new interface formation. Since $\mathrm{\Delta }{G}_V$ and ${\sigma }_{inL}$ (interfacial energy of inclusion and liquid steel) are constants under the current condition, the critical radius for possibly stable inclusion nuclei ($r^{\ast }$) corresponding to free energy change ($\mathrm{\Delta }{G}_{hom}^{\ast }$) are obtained through differentiating as given in Equations (17) and (18). Furthermore, it was found that when $r<r^{\ast }$, the nucleus dissolved into liquid to minimize the system free energy; when $r>r^{\ast }$, the nucleus tended to grow up and become stable. For specific inclusions with a certain radius, the driving force of nucleation was dependent on the formation Gibbs energy change ($\mathrm{\Delta }G$), or supersaturation of comprised elements ($\frac{K}{K^{eq}}$) as given in Equation (15). In this manner, the thermodynamics and kinetics of inclusion formation are connected.

$$\mathrm{\Delta }{G}_{hom}=\frac{4\pi r^3}{3}\mathrm{\Delta }{G}_V+4\pi r^2{\sigma }_{inL}$$

(16)

$$r^{\ast }=-\frac{2{\sigma }_{inL}}{\mathrm{\Delta }{G}_V}$$

(17)

$$\mathrm{\Delta }{G}_{hom}^{\ast }=\frac{16\pi {\sigma }_{inL}{}^3}{3\mathrm{\Delta }{G}_V^2}$$

(18)

As for the rate of nucleation, Volmer and Weber [31] first proposed an expression, and Becker and Döring [32] further improved it, which has formed the basis for almost all subsequent treatments, as described in Equation (19) [36]:

$$I=I_{A}exp[-\frac{\mathrm{\Delta }{G}_{hom}^{\ast }}{k_{b}T}]$$

(19)

In Equation (19), $I_A$ is a frequency factor that is the product of the number of nucleation sites, the atom or molecule diffusion frequency across to the liquid and inclusion embryo interface, and the probability of the particle successfully adsorbing on the embryo. $\mathrm{\Delta }{G}_{hom}^{\ast }$ is the maximum Gibbs energy change for the homogeneous nucleation; $T$ is temperature; and $k_b$ is the Boltzmann constant. For estimating the frequency factor, Turnbull and Fisher [33] proposed an expression, as given in Equation (20):

$$I_A=\frac{N_A{k}_{b}T}{h}exp[-\frac{Q_D}{RT}],$$

(20)

where $N_A$ is the Avogadro constant; $h$ is the Planck constant; and $Q_D$ is the activation energy for diffusion.

Turpin and Elliott [37] applied the above method and estimated the frequency factor (s⁻¹∙m⁻³) with the pertinent data [38] for several oxides in an iron melt: Al₂O₃, 10³²; FeO∙Al₂O₃, 10³¹; SiO₂, 10³⁴; FeO, 10³⁶. Rocabois et al. [39] suggested that the factor ranged from 10³⁵ to 10⁴⁵. Turkdogan [40] and Babu et al. [41] took a value of 10³³ for oxides in their calculations. It was believed that this frequency factor could be considered as a constant during the calculations [37,39,40,41].

On the other hand, based on Equations (18) and (19), nucleation rate is strongly influenced by the critical Gibbs energy change. Meanwhile, interfacial energy also plays an important role (${\sigma }_{inL}$), and is calculated by Equation (21) [42]. The interfacial energy of inclusion and pure liquid as well as contact angle can both be measured by the sessile drop method [43] and calculated by mathematical models together with phase diagrams [44]. The two methods for achieving interfacial energies are normally in binary or ternary systems. The multi-components and multi-phases in liquid steel influence the values. The referred value of interfacial energies between common inclusions and steel are summarized in

Table 1. The interfacial energies between inclusions and liquid Fe-based melt.

Inclusion Types	Interfacial Energies (j/m2)
Al2O3	1.5 [37]; 1.8 [45]; 2.0 [44]; 2.27 [46]; 1.32–0.777 ln(1 + 40$C_O$) 1 [47]
Ti3O5	1.0 [48]; 1.32–0.777 ln(1 + 40$C_O$) [47]
SiO2	1.4 [44]; 1.47 ± 0.23 [49]; 1.7 [50]; 0.9 [37,51]
MnO	0.6 [44]; 1.45 [52]; 1.45 ± 0.23 [49]; 1.2 [51]
CaO	1.5 [44]; 1.7 [50]
MgO	1.2–1.8 [51]
FeO	0.18 [51]; 0.3 [37]; 1.32–0.777 ln(1 + 40$C_O$) 1 [47]
MnS	0.7 [44]; 0.2–1.0 [45]
TiN	0.3 [44]
AlN	1.0 [44]

¹ $C_O$ is the oxygen concentration.

.

$${\sigma }_{inL}={\sigma }_{in}-{\sigma }_{L}cos\phi ,$$

(21)

where ${\sigma }_{in}$, ${\sigma }_L$ are interfacial energies of inclusions, liquid with vapor, respectively, and $\phi$ is the contact angle between inclusions and liquid.

Heterogeneous Nucleation

For simplification and using less uncertain parameters, the aforementioned homogeneous nucleation theory was applied in most simulations. In practice, heterogeneous nucleation is the dominant nucleation format due to the existence of impurity particles and boundaries. Compared with homogeneous nucleation, the smaller energy obstacle of heterogeneous nucleation earns its popularity. Assuming that a sphere inclusion nucleates on a flat surface with contact angle $\theta$, the system Gibbs energy change of this heterogeneous nucleation ($\mathrm{\Delta }{G}_{het}$) is derived in Equations (22) and (23):

$$\mathrm{\Delta }{G}_{het}=\left(\frac{4\pi r^3}{3}\mathrm{\Delta }{G}_V+4\pi r^2{\sigma }_{inL}\right)·f\left(\theta \right)$$

(22)

$$f\left(\theta \right)=\frac{\left(2+cos\theta \right){\left(1-cos\theta \right)}^2}{4}$$

(23)

By differentiating Equation (22), the critical free energy change $\mathrm{\Delta }{G}_{het}^{\ast }$ can be obtained as given in Equation (24), which indicates that heterogeneous nucleation is much easier than homogeneous nucleation. The corresponding critical radius was the same with homogeneous nucleation (Equation (17)).

$$\mathrm{\Delta }{G}_{het}^{\ast }=\frac{16\pi {\sigma }_{inL}}{3\mathrm{\Delta }{G}_V^2}·f\left(\theta \right)=\mathrm{\Delta }{G}_{hom}^{\ast }·f\left(\theta \right)$$

(24)

As for the heterogeneous nucleation rate, it had a similar form as that of homogeneous nucleation as given in Equation (25) [53]. The contact angle $\theta$ needs to be defined when using heterogeneous nucleation, which can be quite challenging to decide since it varies for different cases.

$$I=I_{B}f{\left(\theta \right)}^{1/6}exp[-\frac{\mathrm{\Delta }{G}_{het}^{\ast }}{k_{b}T}]$$

(25)

where $I_B$ is the frequency factor and similar with $I_A$.

2.3.2. Growth

In addition to number density, inclusion content also depends on the particle growth rate. Three mechanisms—diffusion controlled growth, collisions, and coarsening—contribute to the growing up of inclusions [54,55,56]. After an inclusion is thermodynamically stable and the supersaturation satisfies the condition of nucleation, the nucleus begins to grow. The growth is initially promoted by constituents diffusing towards the particle and chemical reaction. In liquid steel, collisions of individual particles lead to further size enlargement with a reducing number density. Coarsening, referred to as Ostwald ripening [57], is caused by larger inclusions growing via the consumption of smaller particles.

Diffusion Controlled Growth

To describe the diffusion controlled growth and dissolution, the diffusion equation of Fick’s second law should be solved under the boundary conditions. It is also needed to satisfy the mass balance equation. Aaron et al. [58] derived the exact solution to calculate the growth and dissolution rate of spherical and planar particles. In addition, three approximated solutions on precipitate size change were compared with the exact solution. The details could be in the original work [56]. Also there might be other solutions and empirical equations. The present work concerned on the applied solutions on inclusion growing in the solidifying steel.

One of the most frequently used expressions to evaluate the diffusion controlled growth rate of a spherical particle was derived by Zener, as given in Equation (26) [59]. The detailed derivation of this equation can be found in the original publication [59].

$$\frac{dr}{dt}=\frac{D_L}{r}\frac{C_L-C_{inL}}{C_{in}-C_{inL}}$$

(26)

where $dr/dt$ is the growth rate of the particle with a radius of $r$; $D_L$ is the solute diffusion coefficient in the liquid steel; and $C_L$, $C_{in}$, and $C_{inL}$ are the solute concentrations in liquid steel, inclusion and at the inclusion-liquid steel interface, respectively.

From Equation (26), it was found that the driving force was mainly dependent on the solute concentration difference in liquid and at the inclusion/liquid interface. For calculating the interfacial concentrations, it was assumed that a thermodynamic equilibrium exists at the interface as expressed in Equation (27). In this equation, the superscripts $P$ and $Q$ represent the formed elements of the inclusions; and $K^{eq}$ are the solubility products of the inclusion under current conditions.

$$C_{inL}^P·C_{inL}^Q=K^{eq}$$

(27)

One further assumption to solve interfacial concentrations has the following possibilities: mass balance according to the stoichiometric formula as given in Equation (28) [60,61]; considering the diffusion of the formed elements as defined in Equation (29) [62]; assuming the ratio of the diffusion fluxes for the formed elements through the inclusion/liquid interface are equal to the stoichiometric ratio as derived in Equation (30) [63,64]. From Equations (27) and (28), it was found that the interfacial concentrations were possibly equal to the equilibrium concentrations. When using Equations (28) or (29), the selection of the controlled element is involved, for instance, oxygen is commonly considered as the controlled element for oxide growth [60,61,62]. For Equation (30), the growth of inclusions is controlled by the diffusions of both solutes. For inclusions with more than two elements, similar equations can be constructed from any two of the elements.

$$\frac{C_L^P-C_{inL}^P}{C_L^Q-C_{inL}^Q}=\frac{x{M}^P}{y{M}^Q}$$

(28)

$$\frac{C_L^P-C_{inL}^P}{C_L^Q-C_{inL}^Q}=\frac{x{M}^P}{y{M}^Q}·\sqrt{\frac{D_L^Q}{D_L^P}}$$

(29)

$$\frac{C_L^P-C_{inL}^P}{C_L^Q-C_{inL}^Q}=\frac{x{M}^P}{y{M}^Q}·\frac{D_L^Q}{D_L^P}$$

(30)

where $M^P$ and $M^Q$ are molar weights of elements $P$ and $Q$; $D_L^P$ and $D_L^Q$ are the liquid diffusion coefficients for the elements $P$ and $Q$.

In contrast, another expression to estimate the growth rate of a spherical particle does not consider the interfacial phenomenon given in Equation (31). The derivation of this mechanism can be found in [65,66,67]. It was found that the driving force for growth is the difference between liquid concentration ($C_L^P$) and equilibrium value ($C_{eq}^P$) instead of interfacial concentration. This simplification has gained in popularity; however, growth is only controlled by the element $P$ in this situation.

$$r\frac{dr}{dt}=\frac{M_{in}}{100·M^P}\frac{{\rho }_{Fe}}{{\rho }_{in}}{D}_P\left(C_L^P-C_{eq}^P\right)$$

(31)

where $M$ is the molar weigh; $\rho$ is the density; $P$ stands for the controlled solutes; and $in$ and $Fe$ mean inclusion and liquid steel, respectively.

Turkdogan [40] proposed an inclusion growth model based on the work of Ham [68], as given in Equation (32). It was assumed that the number of growing inclusions was fixed and each one had its own sphere diffusion zone with radius ${\overline{r}}_0$. This model was derived by solving Fick’s diffusion law under the assumption of a pseudo-steady state, and a detailed description on the formulating process can be found in the original publication [40].

$$\overline{r}={\overline{r}}_0{\left(\frac{C_0-C_L}{C_{in}}\right)}^{\frac{1}{3}}$$

(32)

where $\overline{r}$ is the oxide radius in after growing; ${\overline{r}}_0$ is the radius of reactant diffusion zone, which is defined by the number density of growing inclusions; $C_0$ is the initial concentration of the solute; $C_L$ and $C_{in}$ are concentrations in bulk melt and inclusions, respectively. In this case, oxygen was selected as the diffusion-controlled solute.

Collisions

The collision growth of inclusions in liquid steel even during the solidification process should be taken into account. The radius of particles generated by collisions is usually calculated using the unchangeable total volume and the decreasing number density. According to the theory of collisions, the collision frequency ($N_{ij}$, m⁻³∙s⁻¹) can be calculated using Equation (33) [54]:

$$N_{ij}=\beta \left(r_i, r_j\right)·n_i·n_j$$

(33)

where $\beta \left(r_i, r_j\right)$ (m³∙s⁻¹) is a function of collision frequency of particles with radius of $r_i$; and $r_j$, $n_i$ and $n_j$ are the corresponding number densities of these two group particles.

Normally there are three types of collisions contributing to the growth of inclusions in liquid steel and their collision frequency functions are expressed as Equations (34)–(36) [69,70,71]:

Brownian motion:

$${\beta }_B\left(r_i, r_j\right)=\frac{2k_{b}T}{3\mu }·\left(\frac{1}{r_i}+\frac{1}{r_j}\right)·\left(r_i+r_j\right)$$

(34)

Stokes collision:

$${\beta }_S\left(r_i,r_j\right)=\frac{2\pi g\left({\rho }_{Fe}-{\rho }_{in}\right)}{9\mu }·\left|r_i-r_j\right|·{\left(r_i+r_j\right)}^3$$

(35)

Turbulent collision:

$${\beta }_T\left(r_i,r_j\right)=1.3{\alpha }_T{\pi }^{1/2}·{\left(\epsilon /v_k\right)}^{1/2}·{\left(r_i+r_j\right)}^3$$

(36)

where ${\beta }_B\left(r_i, r_j\right)$, ${\beta }_S\left(r_i, r_j\right)$, and ${\beta }_T\left(r_i, r_j\right)$ are Brownian motion, Stokes, and turbulent collision frequency functions, respectively, for the particles with radius of $r_i$ and $r_j$; $k_b$ is the Boltzmann constant; $T$ is temperature; $\mu$ is the dynamic viscosity of liquid steel; $\pi$ is circumference ratio; $g$ is the gravitational acceleration; ${\rho }_{Fe}$ and ${\rho }_{in}$ are the densities of liquid steel and inclusion; ${\alpha }_T$ is the turbulent coagulation coefficient; $\epsilon$ is the turbulent dissipation rate; $v_k$ is the kinematic viscosity of the melt.

Then, the total collision frequency function can be obtained:

$$\beta \left(r_i, r_j\right)={\beta }_B\left(r_i, r_j\right)+{\beta }_S\left(r_i,r_j\right)+{\beta }_T\left(r_i,r_j\right)$$

(37)

Note that among the three formats of collision, Brownian motion and Stokes collisions are fundamental steps that can happen without liquid flow, and turbulent collisions decide the intensity of collisions in most cases with liquid flow. Based on Equations (33)–(37), the growth of particles resulting from collisions was considered.

Coarsening

Coarsening is derived based on the reduction of interfacial energy. This process is realized through the shrinkage of smaller particles and growth of larger ones. Greenwood theoretically analyzed this process and a change rate of particle size was formulated [72]. Based on the theory of Greenwood, Lifshitz, and Slysov proposed the equation (shown as Equation (38)) to estimate the mean radius change [56,73]. Coarsening is particularly important when the formation of inclusions reaches equilibrium.

$${\overline{r}}^3={\overline{r}}_0^3+\frac{4}{9}\frac{2{\sigma }_{inL}{V}_{in}{C}_{0}D}{R\left(C_{in}^i-C_0\right)}·t$$

(38)

where ${\overline{r}}_0$ and $\overline{r}$ are the mean radius before coarsening and at time $t$, respectively; ${\sigma }_{inL}$ is the interfacial energy between inclusion and liquid steel; $V_{in}$ is the molar volume of inclusion; $C_0$ and $C_{in}^i$ are the concentrations of the controlled solute $i$ at initial state and in inclusion, respectively; $D$ is the diffusion coefficient of solute $i$ in the matrix; and $R$ is the gas constant.

2.3.3. Dissolution

When the formed inclusion is thermodynamically unstable, it starts to dissolve. Considering dissolution as a diffusion-controlled process, Whelan [74] derived the following expression to calculate the dissolution rate as expressed in Equations (39) and (40):

$$\frac{dr}{dt}=-\frac{{\alpha }_{d}D}{2r}-\frac{{\alpha }_d}{2}\sqrt{\frac{D}{\pi t}}$$

(39)

with

$${\alpha }_d=2\frac{C_{inL}-C_L}{C_{in}-C_{inL}}$$

(40)

$$\frac{dr}{dt}=-\frac{{\alpha }_{d}D}{2r}$$

(41)

where $D$ is the diffusion coefficient of the solute in the matrix; and $t$ is the time for dissolution.

During the dissolving of particles, the elements diffuse from the inclusion/liquid interface towards a liquid. If the transit item in Equation (39) is neglected, Equation (41) is derived [75]. Compared with diffusion controlled growth, dissolution is believed to be an inverse process. Note that when putting Equation (40) into Equation (41), it becomes the reverse process of the growth suggested by Zener (Equation (26)).

2.3.4. Behavior of Inclusions at the Solidification Interface

The behavior of inclusions at the solidification interface influences their final compositions and size distribution, and particles can be pushed in the residual liquid or engulfed by the solid phase. The pushed particles are able to transform or grow due to the enriched solutes. In contrast, the engulfed inclusions change little in the solid phase. To investigate this topic, a number of models [76,77,78,79,80] have been developed based on the force balance on the inclusion at the advancing liquid/solid interface. Meanwhile, Confocal Scanning Laser Microscopy (CSLM) has also been applied for in situ observation of the behavior of particles [47,81,82,83,84]. Most models define a critical solidification velocity above which inclusions are engulfed. Wang et al. [84] reviewed the representative critical velocities modeled by different authors and their validity was compared with CSLM experimental results. They found that the models from Stefanescu et al. [85] and Pötschke and Rogge [86] well predicted the pushing and engulfment of the regular liquid inclusions while the critical velocity of irregular Al₂O₃ was underestimated. As a widely applied model [47,87], the critical velocity for particle engulfment ($V_{cr}$) proposed by Stefanescu and Catalina [85] is given in Equation (42):

$$V_{cr}={\left(\frac{\mathrm{\Delta }{\gamma }_0{a}_0^2}{3\eta kR}\right)}^{1/2}$$

(42)

with

$$\mathrm{\Delta }{\gamma }_0={\gamma }_{PS}-{\gamma }_{PL}$$

(43)

where ${\gamma }_{PS}$ and ${\gamma }_{PL}$ are the interfacial energies of particle/solid steel and particle/liquid steel, respectively; $a_0$ is the atomic distance; $\eta$ is the viscosity of liquid steel; $k$ is the ratio of thermal conductivity of particle to that of liquid steel; and $R$ is the particle radius.

In addition to critical velocity, Wu and Nakae [78] derived criteria for particle pushing and engulfment by considering only the interfacial energy balance as presented in Equation (44). With this criteria, the particles are pushed when ${\theta }_{PLS}$ is larger than 90° [87].

$$\cos {\theta }_{PLS}=\frac{{\gamma }_{PL}-{\gamma }_{PS}}{{\gamma }_{SL}}$$

(44)

where ${\gamma }_{SL}$ is the interfacial energy of liquid/solid steel.

On simulating inclusion formation during solidification, Yamada and Matusmiya [88] accounted that particles at the solidification front were trapped by the solid without pushing. The engulfed inclusion content was estimated by the change of solid fractions when assuming that the particles were distributed homogeneously. For more dedicated work, the application of the aforementioned critical velocity models is quite promising.

3. Models on Inclusion Formation

On the basis of the analysis of microsegregation and inclusion formation, considerable efforts have been made to develop a coupled model on the changes of inclusions during the solidification process. With a thermodynamic model, the stability, compositions, constituents, and types of inclusions can be possibly achieved; changes of inclusions during the cooling and solidification process are simulated; and the influences of solute concentrations and enrichment on the formation of inclusions are predicted. Based on formation thermodynamics, kinetic models are able to evaluate the evolution rate of inclusions, and the size distribution and number density are achievable. The effects of cooling conditions and the concentrations of formed elements on inclusion size and amount can be investigated and controlled. In this section, different coupled thermodynamic and kinetic models are briefly reviewed.

3.1. Thermodynamic Models

Many simulations have been performed to predict the precipitations in residual liquid steel based on the calculated solute enrichments [12]. With the development of alloy steels, it is quite desirable to analyze the formation of various inclusions simultaneously. In the 1980s, the first thermodynamic model to simulate the compositional changes of inclusions was reported by Yamada and Matsumiya [88], which coupled the SOLAGSMIX [89] and Clyne–Kurz Model [25]. SOLAGSMIX is a Gibbs energy minimization program that can calculate thermodynamic equilibrium for multicomponent systems. At that time, this program was still in the infancy of ChemSage [90], thus extra data of standard formation free energies for nonmetallic inclusions, the activities coefficients of species in molten steel, compositions of liquid oxides and liquid temperature had to be introduced into SOLAGSMIX. The basic assumptions of this coupled model were as follows: (1) solute enrichments in the residual liquid steel during solidification were estimated by the Clyne–Kurz Model; (2) the existence of an equilibrium between the segregated solutes and inclusion phases in the residual liquid steel at each solidification step; (3) the formed inclusions were distributed homogeneously in the residual liquid steel; (4) the inclusions were trapped by the solidification interface without pushing out, and the inclusions in the solid were inert in future solidification steps. In addition, the values of partition coefficients and diffusion coefficients were needed for the microsegregation calculation. Using the proposed model, the formation process of calcium oxides and sulfide during the solidification of hydrogen-induced-crack resistant steel were analyzed. In this case, Ca was added to control the sulfide shape in the steel. One of the calculated results is shown in Figure 3, which shows the compositional evolution of all inclusion types and the stabilities of the possible complex oxides. Before achieving the solid fraction of 0.5, the amount of CaS increased gradually by consuming CaO. The mass fraction of the various inclusions changed little when the solid fraction ranged from 0.5 to 0.9. At end of the solidification process, CaO became unstable and transformed to CaS due to the strong segregation of S and the corresponding decreased temperature. The released oxygen reacted with Al to form Al₂O₃. Based on the calculation, the formation of CaS can suppress the formation of MnS and the shape can be controlled. In a later study, the authors of Reference [91] proposed a more generalized model through coupling the Clyne-Kurz Model with ThermoCalc [92] where the basic assumptions were similar with the former model. The thermodynamic equilibrium such as phase stability and liquid temperature were calculated using different databases, while the reasonability of microsegregation prediction was enhanced using the local partition coefficient, which was achieved from the equilibrium calculation at each solidification step. The formation of inclusions in stainless steel was calculated as a case study.

Based on the microsegregation model [29], Ueshima et al. [93,94] simulated the behavior of MnS formation during the solidification of resulfurized free-cutting steel. In the calculation, it was assumed that MnS started to crystallize in liquid or precipitate in solid when the corresponding concentration products of Mn and S exceeded the equilibrium solubility. To consider the distribution of Mn and S, extra fine nodes in the MnS precipitating area were divided. Part of the predicted results is shown in Figure 4. From Figure 4, it was found that most MnS crystallized in the interdendritic region and a small amount of MnS precipitate in the dendrite. Comparing the calculated results with the unidirectional solidification tests, the distribution of MnS and solutes Mn and S were well predicted, which suggested that the rate of both crystallization and precipitation were controlled by the diffusion of Mn (the dotted line in Figure 4).

The research group at Institut de Recherche de la Sidérurgie (IRSID) used the same technique as Yamada [88,91] and suggested a model based on Chemical Equilibrium Calculation for the Steel Industry (CEQCSI) [95], an in-house developed software, and the Clyne–Kurz Model [96,97]. In the model, both the stoichiometric and complex solution inclusions could be considered based on the thermodynamic equilibrium calculation with CEQCSI. In one application, the compositions of oxides in semi-killed steel from different industrial processes were well predicted when compared with the experimental results [96]. It was suggested that the contents of alloy elements Ca, Al, and Mg should be well controlled to avoid the formation of harmful inclusions such as alumina and spinel. In another case, the precipitation of the (Mn, Fe, Cr)S solution phase during the solidification of high carbon steels were calculated in [97], and both compositions and amounts showed good agreement with the experiments.

Choudhary and Ghosh [98] described a methodology to predict the formation of inclusions during cooling and solidification. In the cooling process, inclusion changes were calculated using the Equilib module of FactSage [99]. A sequential calculation was performed by coupling the Clyne–Kurz Model and FactSage. In this manner, the segregated solute concentrations (estimated by the microsegregation model) were input into FactSage to predict the inclusion evolutions during solidification. Note that the consumption of the formed inclusion formation was accounted for when inputting the segregated concentrations into FactSage. The inclusion changes of a low carbon Si–Mn killed steel were calculated. The types of inclusion at the end of solidification were the complex liquid inclusion (MnO–SiO₂–Al₂O₃), SiO₂ and MnS. Figure 5 displayed the mass fraction and composition changes of the Liquid inclusion. It was found that the liquid inclusion continuously precipitated during the solidification process (left axis). Its compositions varied with the solid fraction (right axis). The content of Al₂O₃ in the complex inclusion decreased with steel solidification, which was attributed to the consumption of pure alumina formation and the subsequent less segregation. The predicted inclusions match the measured ones from the industrial samples well.

To simulate inclusion behavior during casting and solidification, researchers at the Helsinki University of Technology (now Aalto University) combined InterDendritic Solidification (IDS) [16] software with the thermodynamic library ChemApp [28] (ICA [16,100]). IDS is a more elaborate model for solidification and phase transformation when compared with simple mathematical models. This program was constructed based on a thermodynamic substitutional solution model, a magnetic ordering model, and Fick’s diffusion law, and has a similar diffusion geometry as the Ueshima Model [29]. IDS contains its own database so it can provide solidification-related thermophysical properties such as enthalpy and specific heat. ChemApp could bridge the self-programmed model and databases in the FactSage databank [99]. In these cases, ChemApp calculates the thermodynamic equilibrium for inclusion formation. Figure 6 shows an example of the calculation of inclusion changes during casting and solidification in high carbon steel. In this example, the Ca treatment was expected to modify the hard alumina to soft calcium aluminates, which showed that the formation of various inclusions (including sulfides and complex oxides) could be predicted. The components of liquid slag oxide were mainly CaO and Al₂O₃. The liquid slag phase and CaS were stable at the beginning of solidification. With decreasing temperature, the slag phase transformed into various calcium aluminates before finally transforming to corundum (Al₂O₃). Accompanying these changes was a gradual increase in the amount of CaS. At the end of solidification, the residual liquid contained a high sulfide content due to sulfur enrichment. Based on the well corresponding industrial experience, it was suggested that the modeling work can offer indicative calculations on inclusion formation during solidification. Holappa et al. applied the model to calculate inclusion changes in a Ca treatment Al-killed steel and the predictions were in line with the experimental results [101,102].

Along the line of the former coupling models, the present authors [103] proposed a thermodynamic model coupling microsegregation and inclusion formation using one ChemSage [90] datafile. The thermodynamic equilibrium was calculated using ChemApp [26] to determine the liquid temperature, solute partition coefficients at the solidification interface, and inclusion formation in the residual liquid. The solute enrichment was predicted by the step-wised Ohnaka Model [27], as described in Section 1. At each solidification step, the composition was transferred to ChemApp, which was used to determine the necessary thermodynamic information such as partition coefficients and phase stability. The formed particles at the solidifying front were trapped without pushing as suggested by Yamada and Matsumiya [88]. The logarithm of the coupling microsegregation and inclusion formation was tested through an overall mass balance. The inclusions changed in three Al-Ti alloyed steels were calculated and compared with laboratory experimental results. The inclusion types and compositions were well predicted. The simulated inclusion formation process could indirectly explain the formation of heterogeneous inclusion types.

On the other hand, because of its strong databases, commercial software [92,99,104,105] with cooling modules has become a popular thermodynamic tool for simulating inclusion formation during solidification. FactSage [99] accounts for both the Lever-Rule and Scheil Model [20,21], which are valid for specific solidification processes. MTDATA [104] also includes the Scheil cooling process., and Thermo-Calc [92] and Matcalc [105] coupled the improved Scheil Model, which considered the back diffusion of interstitial solutes such as C and N. Compared with other software, however, Matcalc concentrates more on precipitation kinetics and microstructure evolution after solidification. Using these models, almost all thermodynamic information of inclusion formation was achievable. Note that predictions from the software were highly dependent on the database. In addition to thermodynamics, more reasonable and flexible models that consider kinetics such as the diffusion geometries in the solidification process and nucleation and growth of inclusions are desirable, which is why several simple in-house models were developed.

3.2. Kinetic Models

While the present work concentrated on solidification, investigations into the inclusion formation kinetics in the liquid process such as deoxidation and welding have also contributed many efforts to the development of the modeling work. Hence, the kinetic models on inclusion formation in both the liquid and solidification processes were reviewed in this section.

3.2.1. The Liquid Process

In addition to chemistry, the number density and size of inclusions are also important aspects. In 1966, Turpin and Elliot [37] first applied classical nucleation theory to investigate oxide nucleation in ternary steel melts. In the application, only homogeneous nucleation was considered to avoid complications caused by the introduction of heterogeneous substrates, and the nucleation frequency factor (pre-exponent) for a variety of oxides was estimated. On the basis of the critical nucleation Gibbs energy change, the critical nucleation concentration of oxygen for supersaturation (with the equilibrium concentrations of alloys) in different ternary systems were calculated. Using this method, the effects of super-cooling and interfacial tension on the nucleation of possible nuclei were investigated. At the same time, corresponding cooling and solidification experiments were designed and the experimental results supported the nucleation theory. It also indicated that the interfacial tensions between the oxides and liquid iron were the main limitation on applying the theory to experiments.

In the same year as Turpin and Elliot [37], Turkdogan [40] analyzed the kinetics of nucleation, growth and the flotation of oxide inclusions in liquid melt. It was assumed that the nuclei resulted from the homogeneous nucleation of deoxidation products. The growth of inclusions was controlled by the solute diffusion and the growth rate was derived by Equation (32). Since the existing equilibrium at the inclusion/steel interface and the flux of reactants are equal, only one reactant (oxygen) needs to be considered in the calculation. Furthermore, oxide flotation was accounted for using the Stokes law. With practical consideration, this approach was applied to study deoxidation efficiency and the removal of inclusions by assuming different number densities of the growing inclusions. The calculated results showed that a critical number density could be achieved to reach the highest deoxidation efficiency.

Later, Mathew et al. [106] proposed conceptions that considered not only nucleation and growth, but also Stokes and gradient collisions. Their analysis concluded that the collisions of the inclusions were the main reason for growth in deoxidation products. The importance of the interfacial energy between inclusions and steel were addressed.

Based on the aforementioned concepts, Babu et al. [41] extended the application of classical nucleation theory (Equation (19)) and diffusion-controlled growth (Equation (26)) to weld metal deoxidation. In their model, they further applied overall kinetics to describe transformation extent (ζ) as given in Equation (45). Using this model, the Time–Temperature–Transformation (TTT) curves for various oxides were calculated as an example (Figure 7), and showed that the reaction kinetics of Al₂O₃ were faster than that of SiO₂ when approaching the extent of 0.9, which was attributed to the higher nucleation rate. Additional calculations and analyses were concerned with the influence of oxygen content, deoxidizing element concentrations, and temperature on the inclusion characters. In a subsequent publication [107], the researchers coupled thermodynamics and kinetics as well as weld cooling curves to simulate the inclusion formation. The calculated results from the proposed model (including the composition, size, number density, and oxidation sequence) were verified with the experimental results. The preliminary work on the coupled heat transfer, fluid flow, and inclusion model was also discussed. As a continuous work presented in [60,62,108], the research group published a more elaborated model accounting for growth, dissolution, collision, and coarsening of inclusions in the weld pool where Al₂O₃ was selected as an example. The calculated size distributions agreed with the experimental results, which indicates that this kind of fundamental model could be used to simulate inclusion formation.

$$\zeta =1-exp\left\{-I\left(\frac{8\pi }{15}\right){\left(2D_L^O\frac{C_L-C_{inL}}{C_{in}-C_{inL}}\right)}^{\frac{1}{2}}{t}^{\frac{5}{2}}\right\}$$

(45)

where $I$ is the homogeneous nucleation number density; $D_L^O$ is the diffusion coefficient of oxygen in liquid melt; $t$ is the time for inclusion formation; $C_L$, $C_{in}$ and $C_{inL}$ are the solute concentrations in liquid steel, inclusion and at inclusion/liquid steel interface, respectively; and the solute is oxygen in this case.

Zhang et al. [109,110] established a model considering Ostwald ripening and collision growth instead of diffusion-controlled growth in the deoxidizing process. In the model, the pseudomolecule of the inclusions was assumed as the basic unit and clusters of the pseudomolecule existed before nucleation. The size distribution and evolution at different formation stages were predicted. Later, Zhang and Lee [70] improved the aforementioned model by considering more details on Ostwald ripening and various collisions. At the same time, a numerical method was introduced to reduce the load of the enormous computation. In the following work, a similar mathematical model was proposed by Lei et al. [71]. In their model, the deoxidation products were divided into embryos and inclusion particles. The two parts (with corresponding equations) were separately solved to speed up the calculation, and their predictions on inclusion size distribution were consistent with the experimental results. In addition, the influence of diffusion coefficients and turbulent energy dissipation rate were also evaluated.

3.2.2. During Solidification

Apart from the formations in the melt, inclusions precipitating during solidification have received a large amount of attention, especially after introducing the concept of oxides metallurgy [3,4]. Goto et al. [111,112,113] described a coupled model of oxide growth and microsegregation for studying the precipitations during solidification as described in Figure 8. In the model, Equation (31) was applied to estimate the growth of oxides and the Ohnaka Model [26] was used to predict solute enrichment in the residual liquid. Figure 8 shows that oxides were assumed to form in the interdendritic liquid and grow with the driving force of segregated and equilibrium concentration difference. The consumption of reactants was calculated with a local mass balance. Using the presented model, the effects of the cooling rate on the mean size evolution and supersaturation for oxides were investigated. Note that in the calculation, the number density of oxides was set as a constant based on the experimental results. It was suggested that a higher cooling rate enhanced the supersaturation and frequency of oxide formation, but reduced their size.

In the following study, Ma and Janke [114] predicted inclusion growth through mass balance and also calculated microsegregation with the Ohnaka Model [26]. Based on the former models, Liu et al. [115] applied the Ueshima Model [29] to predict the solute concentration changes in both solid and liquid steel. They calculated inclusion growth by mass balance, while the oxygen content in the solid was accounted, which had been omitted by Ma and Janke [114]. In the models proposed by both Ma et al. [114] and Liu et al. [115], the solute concentrations were assumed to reach equilibrium after the formation of the inclusions. Providing the constant number densities, the influences of cooling rate on the growth of inclusion with various initial radii were studied and it was found that the size of the secondary oxides was greatly affected by cooling rate. Using a similar method, Yang et al. [116] studied TiN growth on the pre-existing MgAl₂O₄ oxide in the solidification process. The initial size of the oxide was assumed, and it was found that the larger size of the oxide limited the growth extent of periphery TiN, which further reduced the proportion of the complex inclusion.

Suzuki et al. [56] proposed a similar model with Goto et al. [111] on inclusion growth in stainless steel, while the Ueshima Model was used for microsegregation prediction. In the same work, the solidification temperature range was divided into 10 regions and Ostwald ripening (Equation (46)) was applied to calculate the inclusion growth. In the calculation, the number of particles was assumed as proportional to the liquid volume and the nucleation rate was set as constant. The initial particle radius was assumed as 1.3 μm. After comparison with the experimental results, it was suggested that the growth of inclusions formed during solidification was controlled by diffusion coalescence.

Osio et al. [45] proposed a model to investigate the effects of solidification on inclusion formation and growth in low carbon steel welds where it was assumed that diffusion controlling the growth of oxides in the deoxidizing process was the primary mechanism for inclusion growth during solidification. In the model, the growth model from Turkdogan [40] (Equation (32)) was simplified and applied. At each solidification step, the nucleation of inclusions was calculated using homogeneous nucleation theory (Equation (19)). Solute enrichments were evaluated by the Scheil Model. The size of particles and their corresponding number densities were tracked to determine the size distribution. The particles at the solidification front were assumed to be rejected into the residual liquid. The influence of local solidification time and solute content on Al₂O₃ formation was studied using the model, and Figure 9 displays the predicted size distribution of Al₂O₃ under different O and Al contents. It was found that both number density and size as well as size range increased with a higher Al content; the increasing oxygen content also resulted in a larger size and number density, which was in agreement with the experimental results [41,117]. A longer local solidification time promoted the growth of particles and reduced the number density. Note that the oxides formed before solidification were assumed to be removed by the weld pool in the modeling process, while the particles with an initial size and number density could be accounted.

Rocabois et al. [39,118] combined classical nucleation theory with microsegregation to describe the formation process of titanium nitride (TiN). In the model, the homogeneous nucleation theory was applied; thermodynamic equilibrium was calculated using CEQCSI [95]; microsegregation was calculated with the Lever Rule; and a mixed control of diffusion and interfacial reaction for inclusion growth was assumed as given in Equations (46) and (47). With Equation (46), the interfacial concentrations and flux could be solved and the particles at the solidification front could be treated as total rejection or engulfment. Using the presented model, the size distribution of TiN was obtained, and the predicted amount of evolution fit well with the experimental results. Next, the model was extended to one complex solution for oxides by Lehmann et al. [119], which enabled the calculation of composition changes and size evolution. The oxide formation in an Al–Ti alloyed low-carbon steel was calculated. Figure 10 shows the evolution of the oxide size distribution. The main components of this complex oxide are Ti₂O₃, Al₂O₃, SiO₂, and MnO and shows that at 1492 °C (which is the initial stage of the inclusion formation), the most numerous inclusions were always the smallest. With growth, particle size with a peak number of densities obviously increased. When the temperature decreased from 1491 °C to 1484 °C, the size distribution shape remained due to decreases of the supersaturation and nucleation rate, and the inclusions can continue to enlarge.

$$\begin{array}{ll}J& =\frac{D_L^{Ti}}{r}\frac{{\rho }^{Fe}}{100M^{Ti}}\left({\left[\%Ti\right]}_L-{\left[\%Ti\right]}_{inL}\right)\\ & =\frac{D_L^N}{r}\frac{{\rho }^{Fe}}{100M^N}\left({\left[\%N\right]}_L-{\left[\%N\right]}_{inL}\right)\\ & =k_c\left(a_{inL}^{Ti}·a_{inL}^N-K_{eq}^{TiN}\right)\end{array}$$

(46)

$$4\pi r^2{M}^{TiN}Jdt=\frac{4}{3}\pi {\rho }^{TiN}d\left(r^3\right)$$

(47)

where $J$ is the molar flux; $D_L^{Ti}$ and $D_L^N$ are the diffusion coefficients in liquid steel of $Ti$ and $N$, respectively; ${\rho }^{Fe}$ and ${\rho }^{TiN}$ are the densities of liquid steel and $TiN$; $M^{Ti}$, $M^N$, and $M^{TiN}$ are the molar weights of $Ti$, $N$, and $TiN$, respectively; $\left[\%Ti\right]$ and $\left[\%N\right]$ are the concentrations of $Ti$ and $N$ referred to a 1% dilute solution; $a_{inL}^{Ti}$ and $a_{inL}^N$ are the activities of $Ti$ and $N$ at the interphase of inclusion and liquid steel respectively; $k_c$ is the kinetic constant; $K_{eq}^{TiN}$ is the solubility product of $TiN$ in liquid steel; $t$ is time; and $r$ is the radius of the particle.

Based on previous work, You et al. [64] proposed a comprehensive model on the formation of MnS during the solidification of steel. The model coupled the formation kinetics of MnS with the step-wise Ohnaka Model, which was linked to a thermodynamic database [27]. Homogeneous nucleation (Equation (27)) and diffusion controlled growth (Equations (34) and (38)) were applied to calculate the formation of MnS. Particle size distribution (PSD) [120] and Particle size grouping (PSG) [121] methods were used to record the size evolution. The collisions of particles in the residual liquid steel were accounted for by inducing a collision factor that considered the normal mechanisms of Brownian motion, Stokes collisions, and turbulent collisions. The collision factor was later calibrated by the experimental results. The particles were assumed to be trapped by the solid phase and the trapped amount was proportional to the step value of the solid fraction [88]. The submerged split chill tensile (SSCT) experiment was used to simulate the solidification process and MnS formation [122,123,124]. The inclusions in the samples were measured using automated scanning electron microscopy/energy dispersive X-ray spectroscopy (SEM/EDS) analysis. With the calibrated model, the evolution of MnS size distribution was predicted as shown in Figure 11 [64] where it was found that the entire distribution shifted to a larger size direction and became flatter as solidification proceeded. This can be attributed to the growth and collision reducing the particle number density. The effects of cooling rate and solute contents on the size distribution of MnS were studied using the model. The good agreement of the predictions with the experimental results indicated the validity of the present model. Figure 12 displays the influence of the cooling rate. Both calculated and experimental results showed that the particle size increased with the decreasing cooling rate. In addition, the total number increased as cooling strengthened. It was suggested that finer particles with a higher number density were achievable by faster cooling, which is beneficial to microstructure optimization.

Descotes et al. [125,126] presented a modeling study on TiN generation and growth during the solidification of a maraging steel (Figure 13). In the model, the heterogeneous and athermal nucleation [127] of TiN on the formed oxide was assumed. The critical supersaturation for TiN nucleation (Equation (48)) was derived based on classical nucleation theory (Equations (25) and (26)). The growth of TiN on pre-existing oxides was calculated using the method suggested by Rocabois et al. [39,118] (Equation (46)), which considered both interfacial reaction and reactant transportation. It was considered that TiN always nucleated on the oxide. A log-normal size distribution of the sphere oxide was assumed and generated by a mathematical method. Solute enrichment was estimated using the Lever Rule (Equation (2)). The particles at the solidification front were engulfed as in the assumption by Yamada and Matsumiya [88]. Once engulfed by the solid, the particles were inert. Figure 14 shows the size distributions of the initial oxide and TiN in the solid and liquid phases at the final stage of solidification. It suggests that the nucleation and growth of TiN happened intensively at the late stages of solidification according to the considerable difference of the distributions in the liquid and solid, which was attributed to the segregated concentrations of Ti and N, and the resultant high supersaturation. It appeared that all the oxide particles were used as nucleation sites for TiN. The predicted maximum particle sizes were found qualitatively in accordance with the industrial observations. Using the model, the effects of the initial oxide number density, N content, and total solidification time were studied. The results offered several unattended trends and understanding of the TiN formation. It was concluded that the final inclusion size increased with the initial N content while this effect was reduced by the prolonged local solidification time.

$$S_c=exp(\frac{2{\gamma }_{TiN/steel}\frac{M_{TiN}}{{\rho }_{TiN}}\sin \theta }{RT\varphi })$$

(48)

$$\sin \theta \approx \frac{\varphi }{r_{nucl}^{TiN}}$$

(49)

where $S_c$ is the critical supersaturation of nucleation; ${\gamma }_{TiN/steel}$ is the interfacial energy between TiN and the liquid or solid steel; $M_{TiN}$ and ${\rho }_{TiN}$ are the molar mass and density of TiN, respectively; $\varphi$ is the radius of pre-existing oxide; $R$ is the gas constant; $T$ is the local temperature; and $\theta$ is the contact angle between TiN and the pre-existing oxide, as estimated by Equation (49). In Equation (49), $r_{nucl}^{TiN}$ is the radius of the nuclei, which can be calculated using Equation (25).

4. Summary and Outlook

Based on the fundamental principles of inclusion formation, a variety of models were developed and applied. This work paid special attention to inclusion formation during solidification.

Table 2. List of thermodynamic models on inclusion formation during solidification.

Process	Author	Year	Reference	Inclusion Stability	Microsegregation		Annotations
					Model	Temperature
Solidification	Yamada	1990	[88]	SOLGASMIX	Clyne-Kurz	Based on Fe-C phase diagram	Additional thermodynamic data
	Yamada	1991	[91]	Thermocalc	Clyne-Kurz	Multi-components	Unified thermodynamics
	Ueshima	1990	[93]	Empirical	Ueshima	Based on Fe-C phase diagram	MnS form in liquid and solid
	Wintz	1995	[96,97]	CEQCSI	Clyne-Kurz	Multi-components	Solution phase is possible
	Choudhary	2009	[98]	FactSage	Clyne-Kurz	Based on Fe-C phase diagram	Methodology
	Nurmi	2010	[100]	ChemApp	IDS	Multi-components	Solution phase is possible
	You	2016	[103]	ChemApp	Ohnaka	Multi-components	Unified thermodynamics

summarizes the coupled thermodynamic model. The widely applied microsegregation models were combined with inclusion formation thermodynamics. Using thermodynamic databases to evaluate inclusion stability is preferable due to its outstanding advantage, while simple empirical equilibrium equation offered an alternative and easy handling method. Note that various commercial software with strong thermodynamic databases is not listed in Table 2, despite being powerful tools to simulate inclusion formation thermodynamics. On inclusion formation kinetics, simulations in the liquid process such as deoxidation and welding greatly promote the development of the modeling work, and were also reviewed. Nucleation theory, growth, collision and coarsening were applied to describe the behavior of the inclusions. In the meantime, process characteristics such as heat input in welding and fluid flow in deoxidation could further elaborate the models. Based on the work, inclusion formation during solidification was simulated and the related models listed in

Table 3. List of kinetic models on inclusion formation during solidification.

Process	Author	Year	Reference	Inclusion Stability	Number	Size		Microsegregation Model	Annotations
						Growth	Collision
Solidification	Goto	1994	[111,112,113]	Empirical 1	Constant	Diffusion	-	Ohnaka	Mean size
	Osio	1996	[45]	Empirical	CN 2	Diffusion	-	Scheil	Size distribution
	Ma	1998	[114]	Empirical	Constant	Mass balance	-	Ohnaka	Mean size
	Rocabois	1999	[39,118]	CEQCSI	CN	Diffusion and reaction	-	Lever Rule	Size distribution
	Lehmann	2001	[119]	CEQCSI	CN	Diffusion and reaction	-	Lever Rule	Size distribution (Solution phase)
	Suzuki	2001	[56]	Empirical	Constant	Diffusion		Ueshima	Mean size
	Liu	2002	[115]	Empirical	Constant	Mass balance	-	Ueshima	Mean size
	Descotes	2013	[126]	Empirical	CN	Diffusion and reaction	-	Lever Rule	Heterogeneous nucleation
	You	2017	[64]	Empirical or ChemApp	CN	Diffusion	Yes 3	Ohnaka	Size distribution

¹ Empirical = empirical free energy equation; ² CN = Classical Nucleation; ³ Yes indicates the item was considered.

. Among the models, particle size evolution could be described using classical nucleation and growth theory, or in a simple way using mean size by assuming the constant number density. Similar to thermodynamic models, different microsegregation models were applied to predict solute concentrations. The kinetic models were more comprehensive by considering both thermodynamics and kinetics given that most of them focus on single inclusion formation. Though these models have made tremendous contributions to controlling and understanding inclusion formation, further developments are still necessary and expected. In the future, work on the following aspects is suggested:

• For both microsegregation and inclusion formation simulations, links to thermodynamic databases offered a new development space. Meanwhile the unified thermodynamic parameters were achievable.
• In addition to the nucleation and growth of a single phase, modeling work on the competitive formation of various inclusions was appreciated to the multi-alloy steels. Another challenging aspect is the heterogeneous nucleation on existing inclusions. Most oxides are generated before solidification and their compositions and size distributions are prerequisite. The subsequently formed inclusions could heterogeneously nucleate on the oxides or other surfaces.
• During the solidification process, the behavior of particles at the solidifying front is necessary for a dedicated inclusion model. The collision of particles is one challenge due to the complex fluid field.
• Aside from the inclusions formed in the liquid, the precipitations in the solid phase also play an important role in the microstructure and properties of steel. In particular, carbides, sulfides, and nitrides, whose precipitation is strongly influenced by microsegregation, are expected to be considered. Furthermore, coarsening also influences the size distribution.
• The melting experiments and inclusion measurements were primary on improving and supporting the calculations.