Wilhelm-Busch-Str. 4,D-30167 Hanover, GERMANY

baake@ewh.uni-hannover.de, nacke@ewh.uni-hannover.de

UMBRASHKO A., JAKOVICS A.

Laboratory for Mathematical Modeling of Environmental and Technological Processes, University of Latvia, Zellu str. 8, LV-1002 Riga, LATVIA

andrey@modlab.lv, ajakov@latnet.lv

Abstract: Experimental results show that heat and mass transfer processes in the turbulent melt flow of induction furnaces are significantly influenced by low-frequency large scale oscillations of the main flow eddies. Large Eddy Simulation (LES) of the turbulent melt flow in induction crucible furnace carries out with good conformity the transient three-dimensional oscillations of the dominating toroidal flow eddies. This 3D transient model offers new possibilities for simulation of heat and mass transfer processes in induction furnaces.

1. Introduction

Fluid flow in industrial metallurgical installations has become a subject of numerical modeling many years ago. The variety of geometries, boundary conditions, fluid forcing factors and a whole range of flow Reynolds numbers provide a good challenge for numerical investigations. The absence of an universal and always reliable modeling approach together with a wide choice of available non-universal turbulence schemes turns it into a non-trivial problem, at least while up to now direct numerical simulation of the turbulent melt flow remains practically inapplicable. As a wide spread example can be mentioned the melting of alloys in induction furnaces. The flow pattern in these installations is formed by the influence of electromagnetic forces and usually takes form of several dominating vortices. Flow patterns obtained with two-dimensional solvers based on Reynolds Averaged Navier-Stokes (RANS) equations usually are in good agreement with estimated and measured time-averaged values [1-3]. Resulting spatial distribution of the temperature and alloys compound concentration depends strongly on the heat and mass exchange between these vortices. Numerical investigations show that two-dimensional turbulence models fail to describe correctly the heat and mass transfer processes between the main vortices. There was developed an engineering approach for this problem described in [1], but for more generic and therefore more flexible solution it is necessary to investigate advanced simulation methods. Today it is possible to run un-stationary three-dimensional numerical calculations of fluid dynamic problems using advanced turbulent models with higher grid resolution requirements and get reliable results in reasonable time. Concluding all these preconditions the calculations presented in this paper were based on Large Eddy Simulation (LES) method, which can be described as a compromise between the solving of RANS equations and Direct Numerical Simulation (DNS). The main idea is that main flow structure is resolved directly, while only small eddies, which size is comparable with grid size, are modeled additionally. Finer meshing is required and, consequently, more computational resources than for two parameter turbulence models, e.g. k-ε turbulence model, but less than it is necessary for the DNS.

2. Low-frequency oscillations

Most of the experimental and numerical investigations were carried out with a model induction crucible furnace (ICF), which has a radius of 158 mm and a height of 756 mm, where the inductor height is 570 mm (Fig.1). Wood's metal, which has a melting point of 72°C, and a dynamic viscosity of 4.2∙10^-3 kg/m∙s, a density of 9,700 kg/m³ and a conductivity of 1∙10⁶S/m was used as a model melt. Experimental velocity measurements carried out using a potential probe show presence of low-frequency flow oscillations in this model furnace (Fig.2). Most intensive of them have a characteristic period about 8-12 seconds depending on inductor current I. These oscillations can be described very simple as periodical up and down movements of the toroidal eddies of the mean melt flow, where the axis symmetry of the melt flow does not longer exist. The main oscillation frequency f increases in dependence on the time-averaged velocity: f ~ v_ch ~ I_ind

Figure 2. Measured axial velocity oscillations at half-height of the inductor near the crucible wall (left) and at half radius (right)

These flow instabilities have a large characteristic time and space scale and therefore they can not be considered as turbulence. So it is doubtfully, that stationary calculations with help of any existing turbulence modeling approach could correctly simulate scalars transport processes in such flow formation. This statement was verified with several widespread turbulence models. All of them do not predict the high intensity of heat exchange in the discussed region. As expected, unstationary 2D calculations using k-ε turbulence model also didn't lead to acceptable results. The absence of azimuthal velocity in 2D analysis limits the dynamics of the large eddies in incompressible fluid. Also 3D calculations using k-ε turbulence model lead to not acceptable results because the low-frequency oscillations are dumped by the high turbulent viscosity. Transient 3D calculations appeared to be the only possible way to resolve the given problem.

3. Numerical modeling

LES and two-equation RANS approach are using different turbulence modeling techniques, therefore, the calculated subgrid turbulent viscosity distribution has principally different character. The highest values of the LES subgrid viscosity are in the zone of eddies interaction in the near-wall region. Another unsimilarity is that turbulent viscosity in case of LES is one order of magnitude less than it is predicted by k-ε model. This difference increases if we improve the spatial resolution of the numerical grid. The smaller become cells, the less energy contain eddies, which are modeled with subgrid turbulent viscosity (Fig. 3). Therefore, further mesh refining will lead to the situation when eddies of all considerable scales are resolved directly. In LES this tendency is obvious when examine the expression for subgrid turbulent viscosity m_sub in Smagorinsky-Lilly scheme used in FLUENT [4-6]:

, ,

where d is the distance from the closest wall, k = 0.42, C_s = 0.1 is Smagorinsky constant and V_c is volume of the computational cell. The averaged radius or cell size directly influences this expression. Kinetic energy of turbulence predicted with k-ε model doesn't depend on spatial resolution, and averaged velocity field doesn't change significantly with mesh refinement. Hence, velocity gradients, which are responsible for modeling of turbulent parameters, also remain the same. Certainly, if we use smaller cell sizes, then additional flow details may appear and influence the resulting turbulence field. But this dependence is not as obvious as in LES. For example, increasing number of mesh elements from ¼ of million to 3 ½ millions the subgrid viscosity calculated with LES becomes one order of magnitude less, but turbulent viscosity of k-ε model remains almost the same (Fig.3).

In k-ε simulation, due to high predicted turbulent viscosity, the viscous forces are comparable with inertial ones and, therefore, all possible large-scale pulsations are damped. Calculating with LES the subgrid viscosity is at least ten times smaller even on relatively coarse mesh. As result, oscillations, which arise from any disturbing factors, develop in time and reach significant amplitude.

Transient calculations started with uninitialized flow field v(r)≡0. After several seconds flow pattern comes to the symmetric state, which resembles results of 2D stationary modeling (Fig.4B), and then begins to oscillate (Fig.4C). There were chosen three control points for velocity along the crucible radius at half-height of the inductor, where the time-averaged velocity is approximately zero. It is possible to analyze the oscillations and compare the time-depended behavior of the velocity components with the experimental results. Velocity oscillations obtained with transient numerical simulation using LES turbulence model are shown on the Fig.5. The Smagorinsky-Lilly subgrid-scale model on tetrahedral mesh with 3.5∙10⁶ elements and with time step 10 ms is taken as a base variant. Other cases differ in only one parameter. In case of 0.4∙10⁶ elements mesh were calculated 130 seconds, in other cases only about 60 seconds. Simulated oscillations amplitudes are in good agreement with those from measurements and don’t depend on number of elements. In good agreement with the experiment the maximum of the oscillation intensity is near the wall and decreases in the center of the furnace (Fig.2,3). Table 1 lists the parameters of different numerical investigations and resulting characteristic values of low-frequency oscillations. The oscillations intensity in the given point is estimated by the time average of axial velocity square.

Figure 4. Flow pattern after 2, 6 and 10 seconds of calculations

Number of grid elements, x10⁶	Subgrid turbulence scheme	Total calculated flow time, s; [time step], ms	Oscillations intensity <V_z²>, cm²/s²	Low-frequency oscillations period, seconds
0.4	S-L	130 [10]	88	14
3.5	S-L	60 [10]	110	12
3.5	S-L	60 [5]	95	8.5
3.5	RNG	60 [10]	60	10.5
Experimental data [7]		56 (measured)	122	9

Table 1. Comparison between calculated and measured low-frequency oscillations intensity and characteristic period between the main flow vortices

Apparently we deal with a nonlinear system of Navier-Stokes equations, which has a solution tending to the limiting cycle. This means, that it is possible to estimate, with a good probability, the flow direction in a given point after relatively long time period. The numerical investigations show, that characteristic period of this main cycle doesn't change significantly when we improve the spatial or temporal resolution. Its length is determined by electromagnetic forces intensity and domain geometry; therefore, the initial conditions have no strong influence.

The discrete particle tracing approach has been carried out to investigate convective scalar transport mechanism in considered flow. In the very beginning virtual particles are placed on the top of computational domain. These particles are assumed to have the same density as fluid and this leads to the expectation that their path will coincide with the streamlines of the flow. When the flow in closed domain without inlets and outlets is stationary, the streamlines also are closed and particle trajectories should be looped. Then it is improbable that particle will penetrate into the neighboring flow region if the turbulent transfer is neglected. Therefore, transport processes between the main flow eddies generally would have diffusive character in steady-state flow. In this case scalar exchange intensity will strongly depend on the semi-empiric turbulent parameters like turbulent viscosity and turbulent Prandtl (δ_t=c_pμ_t/λ_t) or Schmidt (Sc_t=μ_t/ρD_t) number. Latter parameters magnitude often depends both on the type of fluid and on the type of the flow and has to be determined experimentally. The flexibility of using such approach for different industrial installations is rather low.

The subgrid viscosity model seems to have more universal character. Transient simulations, with a small time step and appropriate meshing, allow resolving the wide range of flow formations involved in scalar transport. Particle trajectories, traced in such unstationary simulations, show, that estimated convective mass exchange between the main flow eddies is quite intensive (Fig.6A,B).

Four particles were launched simultaneously at z = z_max, r = r_max/2, δφ =π/2 and traced 20 seconds (about 5-6 eddy turnover times). The choice of starting position was caused by conformity with industrial alloying process, when additional components are added on the melt surface. The typical tracing result states, that only one particle didn't change the host eddy, but others didn't stay in one eddy longer than two or three turnover times. Therefore we can make the conclusion, that convective transport mechanism plays significant role in the heat and mass exchange between the main flow eddies. The same tracing procedure was used with the averaged velocity field from transient LES calculations. As expected, all particles rotated in the initial eddy with relatively small azimuthal drift (Fig.6C). Probably, if we would take longer time for averaging, the trajectories would tend to those in steady-state flow pattern.

4. Conclusions

Large Eddy Simulation turbulence model was applied for three-dimensional transient calculations and it proved to be a very promising tool for numerical simulation of complex turbulent flows. Produced results for velocity field and flow oscillations are in good agreement with experimental data. An extensive numerical analysis was performed to study how spatial and temporal resolution together with subgrid turbulence model choice influence solution stability and reliability. Particle tracing in transient flow showed the high efficiency of convective transfer due to simulated low-frequency oscillations. These results state, that, due to constant increasing of computational power, Large Eddy Simulation could become an universal tool for practical engineering applications.

Further numerical research in this area may include calculations with deformed surface of the melt in order to study its influence on the flow regime. Also the numerical simulation and optimization of the melting processes in industrial installations could be a logical continuation of the presented work.

5. References

1. Nacke B., Baake E., Jakovics A., Umbrashko A.: Heat and mass transfer in turbulent flows with several recirculated flow eddies. Fourth International PAMIR Conference, Presquile de Giens, France, 2000, p. 71-76

2. Jones. R.M., Harvey A.D., Acharya S.: Two-equation turbulence modeling for impeller stirred tanks. Transactions of the ASME, Vol. 123, Sep. 2001, p. 640-648

3. Bojarevics V., Pericleous K., Cross M.: Modeling turbulent flow dynamics in AC magnetic field. The 3^rd International Symposium on Electromagnetic Processing of Materials, Nagoya, Japan, 2000, p. 85-90

4. Fluent 5 User’s guide. Fluent Inc. 1999

5. Smagorinsky J.S.: General circulation model of the atmosphere. Monthly Weather Review. p.91-99, 164, 1963.

6. Lilly D.K.: On the application of the eddy viscosity concept in the inertial subrange of turbulence. NCAR Manuscript, 123, 1966.

7. Baake, E.: Grenzleistungs- und Aufkohlungsverhalten von Induktions-Tiegelöfen. VDI-Verlag Düsseldorf 1994.