Multiphasic Study of Fluid-Dynamics and the Thermal Behavior of a Steel Ladle during Bottom Gas Injection Using the Eulerian Model

: In this work, the ﬂuid dynamic and thermal behavior of steel was analyzed during argon gas stirring in a 140-t reﬁning ladle. The Eulerian multiphase mathematical model was used in conjunction with the discrete ordinates (DO) thermal radiation model in a steel-slag-argon system. The model was validated by particle image velocimetry (PIV) and the analysis of the opening of the oil layer in a physical scale model. The effect of Al 2 O 3 and Mg-C as a refractory in the walls was studied, and the Ranz-Marshall and Tomiyama models were compared to determine the heat exchange coefﬁcient. The results indicated that there were no signiﬁcant differences between these heat exchange models; likewise, the radiation heat transfer model adequately simulated the thermal behavior according to plant measurements, ﬁnding a thermal homogenization time of the steel of 2.5 min for a gas ﬂow of 0.45 Nm 3 · min − 1 . Finally, both types of refractory kept the temperature of the steel within the ranges recommended in the plant; however, the use of Al 2 O 3 had better heat retention, which would favor reﬁning operations.


Introduction
Secondary refining is the process of making steel where it is chemically adjusted to obtain the desired quality. Synthetic slag is commonly used for the removal or reduction of unwanted elements and compounds. The enhancement of chemical reactions between slag and steel compounds occurs using the injection of argon gas [1,2]. Various studies have shown that the fluid dynamics structure [3][4][5][6][7], the mass transfer [8,9], and the chemical reaction rates [10] are affected by the flow of injected gas. In addition, gas stirring helps to eliminate temperature gradients in steel [11,12] generated mainly by exothermic reactions, the heat losses through the ladle walls, electrode reheating, and the ladle residence time [13]. Computational fluid dynamics (CFD) is a technique that helps to analyze the effect of variables such as temperature in the metallurgical process and is currently one of the most used techniques [12,[14][15][16][17]. The thermal control of secondary refining is an important process because it affects the deoxidation and desulfurization processes, as well as the rate of chemical reactions [18]. The use of numerical and experimental modeling in secondary refining has been used to study the heat loss of the steel through the walls of the ladle [19], during the holding period in the furnace [20][21][22][23], the casting to the tundish (teeming) [24,25], as well as the thermal mixing of steel [26]. Castillejos et al. [27] studied the process of thermal homogenization of steel using gas agitation and concluded that after 3 min, the temperature is completely homogeneous.
It has also been shown that the slag layer, in addition to modifying the fluid dynamics structure [9], works as an insulating layer that prevents heat losses in the surface mainly due to the radiation mechanism [14]. In the present work, the CFD technique was used for the complete study of the fluid dynamics and thermal behavior of steel in a metallurgical ladle during argon gas stirring. The Euler-Euler approach was used with drag and non-drag forces using the Eulerian multiphase model, and the Ranz-Marshall [28] and Tomiyama [29] models were compared to calculate the heat exchange coefficient. The fluid dynamics results were validated by means of a physical scale model using the particle image velocimetry (PIV) technique and analysis of the opening of the oil layer. For the calculation of heat losses by radiation, the discrete ordinates (DO) radiation model was used, and the results were compared with those obtained with different methodologies. Finally, the effects of two refractories commonly used in the industry on the variation of the global temperature of steel were compared.

Physical Modeling
The experimentation was carried out on a 1:7-scale acrylic model of the ladle furnace using water, synthetic oil Mobil DTE 24 ® (ExxonMobil, Mexico City, Mexico) and air to emulate steel, slag, and argon, respectively. To carry out the stirring, air was injected by means of a nozzle located at 2/3 of the radius bottom, supplied from a Goni ® compressor (Industrial de Herramientas S.A. de C.V., Mexico City, Mexico) with a capacity of 2 Hp. The thickness of the oil layer was considered 3% of the water height to preserve the slag/steel ratio of the plant system, while the air flow was controlled by a flow meter Cole Parmer ® (Cole-Parmer Instrument Company, Vernon Hills, IL, USA) with a working range of 0-5 Nl·min −1 . To calculate the scaled gas flow, the modified Froude number [30] was employed, described by Equation (1): where, Q is the gas flow, Nm 3 ·s −1 ; ρ is the density, kg·m −3 ; d represents the gas inlet diameter, m; and H is the height of the liquid, m. The air flow used was 2.576 × 10 −5 Nm 3 ·s −1 and corresponds to an industrial argon flow of 0.5 Nm 3 ·min −1 . To obtain the fluid dynamics structure, the particle image velocimetry technique was used through a Dantec Dynamics ® (Dantec Dynamics A/S, Skovlunde, Denmark) time-resolved camera for image correlation and flow field determination. Image capture was performed once fluid dynamics quasi-stability had been reached. For the illumination of the 50 µm polyamide particles, a laser with a time and pulse energy of 100 ns and 10 mJ, respectively, was used. Image capture was carried out with a 3260 fps high-speed camera with a 1 Mp resolution, located at 90 • from the illuminated laser plane. The capture rate was 1000 frames during a 2 s period. A schematic of the experimental setup and dimensions is shown in Figure 1.

Boundary Conditions and Material Properties
For the simulation of the industrial process, a 140-t ladle furnace stirred by an argon gas bottom for a time of 15 min was considered. Gas flow rates of 0.3, 0.45, and 0.5 Nm 3 ·min −1 were injected. For the numerical simulation, the following considerations were taken into account: the calculations were developed in a transitory state, the injected gas flow was constant, the walls had a non-slip condition, the gravitational acceleration was exerted only in the negative direction of the y-axis, the bath temperature was homogeneous at the start of the simulation, and the thickness of the slag layer was uniform at the beginning of the calculations. Table 1 shows the parameters and properties of the materials considered in the experimental and numerical techniques.  The viscosity of the slag is a function of temperature and composition, and, in this work, a slag composed of 55% CaO, 30% Al 2 O 3 , 7.5% MgO, 7.5% SiO 2 was used. The calculation of the viscosity, µ slag , was performed using Equation (2) obtained from data reported in the literature [31]: The dimensions of the ladle are shown in Figure 2a, as well as the boundary conditions and the distribution of the phases considered. For the analysis of the temperature into the steel during the simulation, 15 monitoring points were established, homogeneously distributed throughout the volume of the steel, and thus guaranteeing the global evolution of the temperature, as shown in Figure 2b. On the other hand, Table 2 shows the position of each monitoring sensor in cylindrical coordinates. Figure 2c shows the employed mesh, and the sensitivity technique was chosen to optimize its refinement based on the axial velocity of the steel at the middle ladle height, as shown in Figure 3. According to the results, no differences between the mesh of 1,247,337 elements and that of 1,500,000 elements were found; the first being selected to perform the calculations and thus to reduce the computational cost. Two different heat fluxes were employed, which were associated with the two refractory types most used in the plant for metallurgical ladles, the Al 2 O 3 and Mg-C heat fluxes of −6400 and −13,750 W·m −2 , respectively [14].   The solutions of the continuity, momentum, and energy equations were performed using the Eulerian multiphase model [32], included in the Ansys Fluent 15.0 ® software (Ansys Inc., Canonsburg, PA, USA). This model solved the equations and performed the mass and energy balances using a pressure-based solver, considering a convergence criteria 1 × 10 −4 for all variables. A Dell Precision T1700 ® Intel Xeon E3-1241 v3-3.5 GHz, 16 GB RAM (Intel Corporation, Mountain View, CA, USA) was employed.

Momentum Equation
The force balance in system was solved by Equation (4) [32]: Here, p is static pressure, Pa; ρ → g is force due to gravitational acceleration, m·s −2 ; → F li f t,q and → F td,q are lift force and turbulent dispersion force for phase q, respectively (in N); and → R pq is the interaction force between phases. The strain-stress tensor, τ q , is defined by Equation (5): where, µ q and λ q are the shear and bulk viscosities for phase q, respectively [32].

k-ε Realizable Turbulence Model
Turbulence was modeled using the k-ε realizable model [33], and it is described by Equations (6) and (7): where, k and ε are kinetic turbulent energy, J·kg −1 , and the turbulent dissipation rate, m 2 ·s −3 , respectively. µ t is the turbulent viscosity, m 2 ·s −1 . σ k and σ ε are the Prandtl number for k and ε whose values are 1.0 and 1.2, respectively [33]. G k and G b represent the generation of turbulence kinetic energy due to the mean velocity gradients and buoyancy, respectively. C 1ε and C 2 are constants, and their corresponding values are 1.44 and 1.9 [33]. C 1 was calculated according to Equations (8)- (10): Turbulent viscosity was determined by Equation (11): C µ was calculated using Equations (12)- (19): Here, u represents the mean velocity component. S is a scalar measure of the deformation tensor. Ω ij is the mean rate-of-rotation tensor along the moving reference frame, and ω k represents angular velocity.

Troshko-Hassan Turbulence Interaction
The influence of the disperse phase (argon and slag) on the turbulence equations was modeled using the Troshko-Hassan model [34]. For continuous phase q, the source terms were calculated by Equations (20) and (21): where, → U is the phase-weighted velocity of the phase; C ke and C td are constants whose values are 0.75 and 0.45, respectively; and the term K pq is the interphase momentum exchange coefficient, and it is related with the momentum equation (Equation (4)) according to Equation (22): → v p and → v q are the velocities of each phase, and τ p is the characteristic time of induced turbulence defined by Equation (23): Here, C V M is the virtual mass coefficient, and C D represents the drag coefficient. To calculate the kinematic turbulent viscosity (ν) of the dispersed phase, Equation (24) is used:

Symmetric Drag Model
To determine the interphase momentum exchange coefficient, K pq , the symmetric model was employed [32]. This model uses Equations (25) and (26) to calculate properties such as density and viscosity: Particle diameter was determinate, using the diameters of phase q and p, by Equation (27): K pq was calculated by Equations (28)- (31): Here, Re expresses the dimensionless Reynolds number. A i is the interfacial area, m 2 ; f is the drag function, and it is calculated according to the criteria described in Equation (32) [33]:

Lift Model
The lift effect was included using the Tomiyama lift model [29]. According to the model, the lift force can be determined by Equation (33): where, C l is the lift coefficient, which is calculated according to Equations (34)- (38): Eo represents the dimensionless Eötvös number; this number relates buoyancy and surface tension forces. Eo is a modified Eötvös number, and it is based on the longer axis of the bubble, d h . σ is surface tension, N·m −1 ; g is gravity acceleration, m·s −2 ; and d b is bubble diameter, m.

Virtual Mass
The virtual mass force exerted on a rising bubble is calculated by Equation (39) [35]: where, the phase material time derivative is formulated as Equation (40):

Turbulent Dispersion Force
The force due to the turbulent dispersion effect was modeled using Burns et al.'s model [36]. Here, the dispersion scalar, D tq , was calculated with Equation (41): where, D q and D p are the scalar dispersion of the continuous and the disperse phase, respectively. Thus, the turbulent dispersion force was determined by Equation (42): where, σ pq is the dispersion Prandtl number, and its value is 0.9; the C TD constant value was set to 1 [36].

Energy Conservation
The energy equation from the Eulerian multiphase model is described by Equation (43) [33]: where, h q is the specific enthalpy for phase q; → q q is the heat flux, W·m −2 ; and Q pq represents the intensity of the heat exchange between phases.

Heat Exchange Coefficient
The volumetric energy transfer rate between phases, Q pq , in W is defined by Equation (44) [33]: where, h pq is the volumetric heat transfer coefficient between phases p and q; A i represents the interfacial area, m 2 ; and T p and T q are phases p and q in K, respectively. The volumetric heat transfer coefficient was calculated from the Nusselt number, Nu, using Equation (45) [33]: where,κ q is the thermal conductivity of phase q, W·m −1 ·K −1 , and d p is the bubble diameter, m. There are several models to determine h pq ; however, two models were considered in this study. The Ranz-Marshall model [28] proposes Equation (46)  Here, Re and Pr are the Reynolds and Prandtl dimensionless numbers, respectively. The Prandtl number was computed according to Equation (47): where, C p q is the heat capacity of phase q, J·kg −1 ·K −1 .

Discrete Ordinates Radiation Model
The discrete ordinates radiation model [37][38][39] was employed to solve the radiation equation. This model uses a finite number of discrete solid angles associated to the vector direction, and it is described by Equation (49): It is applicable for an emitting, absorbing, and scattering medium. → r , → s , and → s are the position, direction, and scattering direction vector, respectively. s, n, and a represent the path length, the refractive index, and the absorption coefficient (0.9 and 0.83 for steel and slag, respectively). σ s is the scattering coefficient, and σ is the Stefan-Boltzmann constant, 5.669 × 10 −8 W·m −2 ·K −4 . The radiation intensity and the local temperature are defined by I and T, respectively. Φ is a phase function, and Ω is the solid angle. The emissivity was set as 0.3 [40].

Fluid Dynamics Validation
For the validation, the results of velocity magnitude and slag layer opening obtained by numerical simulation were compared with those obtained by the PIV technique, as shown in Figure 4a,b, respectively. The effects of all drag and non-drag forces were included in the simulation according to the conclusions established in previous studies [6,41]. The highest velocity was reached in the plume region because of the rising air bubbles. The transfer of the momentum between the bubbles and the water produced a movement of the liquid from the bottom towards the free surface, producing a recirculation pattern as observed in points 1 and . The location of the recirculation pattern and the magnitudes of the water velocities predicted adequately fit those obtained experimentally. Figure 4c,d compares the openings of the oil layer at the same experimental and numerical calculation time,~15 s. The opening was primarily located in the region of the broken free surface by the plume; however, the convective movement and momentum exchange at the interface prolonged the opening beyond the spout. It was observed that the shape of the oil aperture was very similar in both cases, and an oil aperture difference of~3.7% was determined using the ImageJ ® analyzer software (National Institutes of Health, Bethesda, MD, USA).

Behavior of the Gas Flow Rate and the Slag Eye Aperture
The effect of the injected argon flow, 0.3 and 0.5 Nm 3 ·min −1 , on the fluid dynamics structure of the steel is shown in Figure 5 for the industrial ladle. The velocity of the plume was slow using the gas flow of 0.3 Nm 3 ·min −1 because of the ferrostatic pressure exerted by the steel that delayed the ascending gas plume (Figure 5a). In both cases, there was the formation of a recirculation pattern, but it was located in different positions. With low injection, the momentum was less, and the drag effect was increased, causing a recirculation in the middle of the ladle at point 1 . With a flow of 0.5 Nm 3 ·min −1 (Figure 5b), the recirculation position was located close to the free surface caused by a higher momentum exchange between argon and steel, point . The injected gas flow influenced the shape and size of the opening slag eye. A flow of 0.3 Nm 3 ·min −1 (Figure 5a) had a slightly smaller opening than when using a high gas flow (Figure 5b). This difference was attributed to a lower impact force at the steel-slag interface, which led to a reduction in slag displacement. The quantitative difference of the aperture areas was~1.1% and was similar to that reported in the literature [14].

Thermal Behavior of Steel in the Ladle
The Ranz-Marshall model [28] has been used to calculate the thermal behavior of steel in a satisfactory way [11]; however, in this work, the models of Ranz-Marshall and Tomiyama [29] were compared. The latter was used because it was used in the solution for the lift forces in the fluid dynamics calculation. Figure 6 shows the steel temperature predicted by the Ranz-Marshall and Tomiyama models in a 15 min process period. The injected gas flow was 0.45 Nm 3 ·min −1 , and a boundary condition of heat loss through the walls of −6400 W·m 2 was used, omitting the radiation mechanism for this case. During the first 5 min, there was no homogeneity in the temperature of the steel because of the transfer of heat in the plume and the walls. However, after 5 min, the temperature was the same in all the sensors, and they dropped to a temperature of approximately 1868 K. This final temperature did not show significant differences between the models since both reached a value of 1868 K. In the industrial process, the usage of argon gas generates the opening of the slag layer, causing the steel to be exposed to the atmosphere, favoring the transfer of heat through the radiation mechanism according to Farrera et al. [14]. This is shown in Figure 7, where the effect of the incorporation of the radiation model is compared. The numerical calculation using the DO model yielded a difference of 9 K in the final temperature of the steel when it was not used. The initial drop in temperature at 30 s was more pronounced when the DO model was used because of the removal of heat by the free surface, and, subsequently, it descended at a constant rate. The foregoing demonstrates the importance of incorporating a radiation model for the analysis of the thermal behavior of steel in a secondary refining ladle.   Figure 8 also shows the comparison of the plant data [14] with the results obtained in the mathematical model under the same conditions with the Mg-C refractory. In this, it can be seen that the temperature difference at the end of the calculation time was only 1 K and that both slopes are similar, allowing the model to be used with confidence to predict the thermal evolution of the steel. In the industry, a final temperature range of the steel is recommended (delimited by the dotted lines), and both numerical models using different refractory materials are kept within the range.
The temperature contours of the thermal evolution of the steel using an injection flux of 0.45 Nm 3 ·min −1 are shown in Figure 9. Initially, the temperature of the steel was at 1883 K according to the initial condition (Figure 9a). At 30 s (Figure 9b), the effect of radiation caused a decrease in temperature near the free surface of the steel. At 1 min (Figure 9c), the steel remained hotter in the lower zone of the ladle, showing a drop in temperature as it approached the free surface. The influence of the plume was also noticeable when modifying the contour in the upward zone, as indicated in point 1 . At 1.5 and 2 min-represented by Figure 9d,e, respectively-temperature gradients were significantly reduced, reaching the total thermal homogenization of the steel at 2.5 min, as shown in Figure 9f. This thermal homogenization time was similar to that reported by Castillejos et al. [27] where they established that convective mechanisms are dominant in the homogenization of temperature. However, they did not include the effect of free surface radiation, and the correspondence in this homogenization time was due to the fact that the authors [27] employed a greater through-wall cooling rate of −12,500 W·m −2 .

Conclusions
In this work, the fluid dynamics and thermal evolution of a multiphase system were studied using the Eulerian model with drag and non-drag forces and the Ranz-Marshall, Tomiyama, and DO models during the gas stirring process in the ladle furnace. According to the results, the following conclusions were drawn:

•
The selected multiphase models correctly predicted the fluid dynamics of a threephase system (water-oil-air). The velocity of the water and the shape and area of the opening of the oil layer predicted mathematically had good agreement with the experimental data, through the PIV technique. • It was observed that there were no significant differences between the results calculated for the coefficients of heat exchange between the Ranz-Mạrshall and Tomiyama models. However, the Tomiyama model was chosen because of its affinity with the lift model for fluid dynamics simulation. • Despite both refractories retaining enough heat to keep the steel temperature within industrial limits, the Al 2 O 3 one had better heat retention, which would favor refining operations. However, for the desulfurization process, the use of basic refractories such as Mg-C is recommended to favor the reactions involved, and the operational cost is lower.

•
The inclusion of the radiation mechanism contributes to the precision of the simulation of the thermal evolution of the steel into the ladle furnace. The DO radiation model correctly predicted heat loss through the free surface, reaching steel thermal homogenization at 2.5 min. Funding: This research received no external funding.

Data Availability Statement:
This study not report any data.