Numerical simulation of DC casting of large-size rare earth magnesium alloy ingot under low-frequency electromagnetic field

For studying the changes of macro-physical field in the casting process of large-scale rare earth magnesium alloy, through the numerical simulation method, a two-dimensional axisymmetric multi-physical field coupling model was established by using the multi-physical simulation software COMSOL Multiphysics. The changes of temperature field, flow field, Lorentz force, and liquid fraction of large-size rare earth magnesium alloy with diameter of 750 mm under different electromagnetic parameters (magnetic field frequency and current intensity) in steady state of direct-chill (DC) casting were studied. The results reveal that using a magnetic field can reduce the temperature gradient and greatly accelerate the melt flow, the depth of the sump is reduced by about 50 mm. As the current intensity rises, the flow rate in the melt becomes accelerated, the sump depth becomes shallower, while the melt area with a liquid fraction of 0.5 to 0.63 increases. The Lorentz force rises as the magnetic field frequency increases, but the skin depth of the magnetic field decreases from 64.9 to 36.4 mm.


Introduction
Recently, magnesium alloys have been used in automotive lightweight materials, aerospace materials, automotive wheels, and other applications because of their low density, high specific strength, and specific stiffness [1,2]. However, traditional magnesium alloys still have some significant disadvantages, such as low absolute strength, poor mechanical properties at high temperature, low plasticity at room temperature, poor processing deformation ability, and corrosion resistance. Rare earth magnesium alloys have higher strength, better corrosion resistance, and heat resistance than traditional magnesium alloys, which makes alloys with high rare-earth content become a research hotspot in recent years [3]. Nevertheless, the shrinkage coefficient of high rare earth alloy is large, and it is very easy to form shrinkage holes and microcracks in the solidification process, so it is easy to have cracks and other problems in the preparation of high rare earth alloy. Therefore, the procedure for preparing rare earth magnesium alloy must be chosen carefully. The traditional casting methods include semi-continuous casting and die casting, among which the semi-continuous casting method is mostly used in the manufacture of magnesium alloy billets. It has high production efficiency and is suitable for batch production, a high degree of mechanization and low energy consumption, and other advantages [4,5]. The traditional DC casting often presents coarse dendrite, columnar structure, and component segregation, especially the large-size billet, which leads to poor deformation and mechanical properties [6,7]. Consequently, it is particularly important to use certain means to refine the grain and reduce the defects in the casting process. Many studies have shown that applying external fields such as magnetic fields can effectively reduce ingot defects and improve ingot quality.
Zhang et al. [8] established a model describing the interaction of multiple physical fields in the process of conventional DC casting and low-frequency electromagnetic casting (LFEC) with commercial software ANSYS and FLUENT. It was found that the low-frequency electromagnetic field had a significant impact on the velocity distribution, temperature distribution, and sump depth in the process of DC casting. Ma et al. [9] found that axial magnetic force leads to melt convection, radial magnetic force leads to melt vibration, and Joule heat generated by pulsed magnetic field is concentrated near the melt surface during the pulse action period. Chen and Shen [10] carried out numerical simulation research on solidification characteristics under pulse magnetic field, harmonic magnetic field, and absence of magnetic field. Hatic et al. [11] used a meshless numerical model of direct chill casting under low-frequency electromagnetic field to study the effects of low-frequency electromagnetic force on temperature, liquid fraction, and fluid flow under different current intensities and frequencies. Jia et al. [12] established a transient two-dimensional axisymmetric mathematical model of the coupling of pulsed electromagnetic field with fluid flow and solidification by using COMSOL Multiphysics software, and they simulated and discussed influence of magnetic field on the DC casting process of AZ80 magnesium alloy. The results show that the Lorentz force and melt convection are significantly enhanced with the increase of current intensity, and the increase of frequency has the opposite effect. Duan et al. [13] created a two-dimensional axisymmetric model by using the finite element method and discussed the influence of coil connection mode on the DC casting process of magnesium alloy. It indicates that the inverse pulse magnetic field can be used to regulate the range and intensity of melt convection by modifying the coil.
There are many researches on low-frequency electromagnetic semi-continuous casting process, but there are few researches on DC casting process of large-size rare earth magnesium alloy, and DC casting is carried out at high temperature, magnetic field, temperature field, and velocity field are difficult to measure, and the cost of experiment is high. With the advancement of computing technology, numerical simulation methods can now accurately simulate the DC casting process in low-frequency magnetic fields; hence so, the DC casting process of large-size rare earth magnesium alloy under low-frequency electromagnetic field was studied by using a two-dimensional axisymmetric mathematical model, the effects of magnetic field, frequency and current intensity on melt flow, Lorentz force, temperature, and solidification properties were discussed. The appropriate process parameters can be obtained by changing relevant parameters, which can also play a good guiding role in the process formulation of actual production. Figure 1 shows a schematic diagram of electromagnetic semi-continuous casting. Melt is injected into the shunt plate from the inlet and solidified after passing through the primary cooling zone of the crystallizer, air cooling zone, and secondary cooling zone. The two-dimensional axisymmetric geometric model is used in this research to simplify the calculation, and the infinite remote in the magnetic field is replaced by the infinite field, as shown in Fig. 2. Table 1 shows the geometric model dimensions.

Flow field and temperature field
In this study, a single region volume-average model is adopted, which applies the governing equation to all regions in the solidification process. All models are integrated, and there is no obvious distinction between liquid region, mushy region, and solid region. Instead, according to the energy distribution obtained by solving the model equation, each region is implicitly defined [14], and the equation is as follows: Mass conservation equation: Momentum conservation equation: (1) t + ∇ ⋅ ( U) = 0 and T are the laminar viscosity and turbulent viscosity of the liquid respectively, U is melt flow velocity, S m is the momentum source term, including thermal buoyancy ( S i ) [15], Darcy source term ( S ), and outfield source term ( F m ), the formula is as follows: where is the coefficient of thermal expansion, T 0 is the reference temperature, S and F m are given in Eq. (17) and Eq. (26), respectively. Energy conservation equation: where C eff is the equivalent specific heat, including the specific heat and the release of latent heat, Q is the Joule heat.
In the process of DC casting, the turbulence model is often used to deal with the natural convection in the melt and the forced convection of the melt under the action of electromagnetic stirring. The development from turbulence to laminar  flow is a gradual transition process. The movement of solidification interface of ingot is treated by wall function method. In this study, the fluid flow is described using the standard k-ε turbulence model, and the viscous flow near the wall is analyzed using the wall function approach. The turbulence model equation is as follows: Turbulent kinetic energy k: Turbulent flow energy dissipation rate : Source term: The turbulent viscosity coefficient is T [16]: where S k and S are modified source terms of turbulent kinetic energy equation and dissipation rate equation, respectively.
where A is the parameter of paste zone, its value is 10 5 kg∕(m 3 ⋅ s) ; is a little amount set to keep the denominator from going to zero. Table 2 contains the constants used in the k-ε model.
The solidification of alloy is realized in a certain temperature range, from liquid zone to mushy zone to solid zone. The single region volume model is used in this paper, and the whole region is regarded as a governing equation. The equivalent specific heat method is used to solve the energy equation in the solidification process in this simulation. Assuming the latent heat is completely released in the mushy region: where L is the latent heat of crystallization. In this paper, Gaussian curve is used to describe the release of crystallization latent heat in the mushy region [17]: where T m is the melting point and ΔT is half the transition temperature of the mushy area, then the equivalent specific heat is as follows: where C p is the specific heat. f s and f l represent the solid fraction and liquid fraction, respectively. They are a function of temperature and are calculated using the lever law, which has the following formula: where T l and T s are the liquidus temperature and solidus temperature, respectively. The melt flow in the mushy region is often treated by adding a source term to the momentum equation, the change term of the melt is regarded as the momentum source in the equation of porous medium, and the change term of the melt is considered according to the Darcy's theorem [18]: where K is the permeability of porous media, U s is the casting speed. The calculation stage of the model can be separated into four parts, with the permeability of each section provided in Table 3 [8].
In the table, f * s is the solid rate of dendrite overlapping. In this paper, permeability can be regarded as a function of solid rate. K 0 is the initial permeability, which is linked to the internal structural parameters of the ingot like the secondary dendrite arm spacing (DAS) [19]. Calculation formula:

Electromagnetic field
In the process of electromagnetic casting, the electromagnetic field equation is based on Maxwell's equation and Ohm's law, the differential equation is as follows [20]: Ampere's law: Faraday's law: Gauss law of electric field: Gauss law of magnetic field: where B is the magnetic induction intensity, H is the magnetic field strength, B = 0 H , 0 is relative permeability; D is the electric induction intensity, E is the electric field strength, D = 0 E ; 0 is the dielectric constant, J is the current density, 0 is the resistivity, and t is the time.
To solve the above equation, the constitutive equation of J is given by Ohm's law: where is the conductivity. By deriving the above equation, the formula of electromagnetic force is obtained: Joule heat during electromagnetic casting: Since low-frequency electromagnetic casting produces low-frequency time-varying electromagnetic field, and the period of electromagnetic field is much shorter than that of melt momentum response, the electromagnetic force and electromagnetic heat in the process of low-frequency electromagnetic casting are [21]:

Thermophysical properties and boundary conditions
In this simulation, the temperature fields and flow fields were simulated using the commercial multi-physical simulation software COMSOL Multiphysics. Thermal conductivity, specific heat, and other important parameters in the simulation process were calculated by Jmatpro. The changes in thermal conductivity and specific heat with temperature are shown in Fig. 3a, b, respectively. Table 4 illustrates other physical parameters and the initial conditions of this simulation. In addition, the control variable method is adopted in this paper. When a certain parameter is changed in the simulation process, other parameters remain unchanged. Figure 4 depicts the boundary conditions of the flow field and thermal field used in the simulation procedure, the top surface and mold surface are regarded as static walls, for heat transfer conditions, the Cauchy boundary conditions are used, as follows: The heat transfer coefficient of the air cooling area is set to be constant: 15 W/(m K). The secondary cooling area and air cooling area are regarded as moving walls, for the secondary cooling zone, use the following formula to calculate the heat transfer coefficient [22]: where Q c is the flow rate of cooling water, T water is the temperature for cooling water.
The central part is the axis of symmetry, and the adiabatic boundary condition is adopted.
Boundary conditions at the inlet: the pouring temperature is 953 K, velocity of inlet melt: where R s is the ingot radius, and R z is the inlet radius. Turbulent kinetic energy at inlet: Energy dissipation rate of turbulent flow at inlet: where I n is turbulence intensity: where R e is Reynolds number: The boundary condition of the magnetic field is relatively simple. Current is applied to multiturn coil, and magnetic insulation is set at the axis of symmetry and the periphery of the infinite field.  Specific heat, C p Figure 3b Liquidus temperature,

Simulation procedures and model assumptions
The temperature field and flow field of semi-continuous casting are simulated in this paper using the commercial finite element software COMSOL Multiphysics. Through the physical field interface of fluid heat transfer, input the casting parameters and set the boundary conditions, the coupling interface adopts non-isothermal flow, the coupling of temperature field and flow field can be realized, and the temperature field and flow field distribution of melt can be calculated. When adding the magnetic field, the coupling interface adopts Lorentz force coupling, by changing the magnetic field conditions, the analysis of melt temperature, flow field, and solidification process under different magnetic fields can be achieved. Since this study involves the coupling between multiple physical fields, some reasonable assumptions are used to simplify the problem: 1. The melt is treated as an incompressible fluid and the Boussinesq model is used to calculate the thermal buoyancy in the melt [23]. 2. This model does not calculate the meniscus and the solute field. 3. Only steady-state simulation is considered in this simulation procedure. 4. There is no consideration for the formation of a solid shell. 5. Melt flow has little effect on magnetic field, because the magnetic Reynolds number R m < < 1 ( R m = 1 U 0 L 0 = 0.028 , where 1 and σ are permeability and electric conductivity, respectively, U 0 and L 0 are characteristic velocity and characteristic length, respec-tively). As a result, in the equation J = (E + U × B) , U × B is ignored [17].

Grid and model validation
In the finite element simulation, the mesh quality affects the convergence of the model. If the mesh is too large, the convergence is poor, the mesh is too small, and the calculation time is long; therefore, the quality of the grid is an important factor affecting the success of simulation. Because the simulated transition region (the region where the phase transition occurs) requires a high degree of discretization, thus, this paper uses the built-in adaptive meshing technology to solve the model. Figure 5 shows the grid used in this simulation. To test the independence of the grid, three kinds of grids with the number of 15,891 (mesh 1), 23,781 (mesh 2), and 47,154 (mesh 3) are used to calculate the temperature distribution without magnetic field. And draw the temperature change along path 3. As can be observed in Fig. 6b, changing the grid has little impact on the simulation results. The simulation is tested by the casting experiment of Bao et al. [24] to ensure its dependability, and the distribution of sump in DC casting process and LFEC process was compared; the measured values in the figure are the values reported in the literature, and the casting parameters and ingot specifications used in the simulation are consistent with the experiment. Figure 7 is a comparison diagram of sump distribution between simulation result and experimental result, as can be observed that the distribution of sump in the experiment and simulation is very similar. However, due to the ideal two-dimensional axisymmetric model used in the simulation, the sump is completely symmetrical, which is somewhat different from the actual situation, it shows that this simulation is still reasonable to a certain extent.

Lorentz force
In engineering, the depth at which the current density is reduced to about 37% of the surface is called skin depth [25]. Figures 8 and 9 are the Lorentz force distribution diagram and skin depth variation diagram, respectively. Figure 10 displays the variation of Lorentz force along path 3; as depicted in the figure, the direction of Lorentz force points to the lower-left corner. It illustrates that Lorentz force has axial component and radial component; the axial component promotes the stirring of the melt and the radial component has constraints on the melt. The Lorentz force at the edge of the ingot rises as the current increases, and the skin depth remains unchanged. With the increase of the magnetic field frequency, the Lorentz force at the edge of the ingot also increases, and the Lorentz force has a smaller range of action. The larger the frequency is, the smaller the Lorentz force is near 50 mm from the edge of the ingot; this is caused by the skin effect of electromagnetic field. With the increase of frequency, the included angle of Lorentz force decreases, the stirring ability increases, and the constraint decreases. Figure 11 demonstrates the flow field changes during DC casting and LFEC casting. Figure 12 shows the velocity distribution of melt after flowing out of the shunt plate. After injecting the high-temperature melt onto the shunt plate, a counterclockwise vortex arises inside the shunt plate due to the obstruction of the wall of the shunt plate, and the liquid at the edge of the vortex flows out of the shunt plate under the action of inertia, forms the solidification front after cooling in the mold, and forms a clockwise vortex under the joint action of the solidification front and thermal buoyancy; the remaining melt flows downward under inertia. The overall speed of melt improves after the low-frequency magnetic field is applied, and the skin effect of the magnetic field makes the vortex outside the shunt plate close to the mold. Under the action of the electromagnetic field, the melt has forced convection and a smaller vortex in the opposite direction is generated below the large vortex, which is brought on by the Lorentz force on the melt stirring.  Figure 13 displays the distribution of flow field in melt under different magnetic field frequencies and different current intensities. It can be observed that as the magnetic field frequency rises, the flow field in the melt moves to the edge of the ingot as a whole. As the current intensity rises, it can be seen from the arrow on the surface that a small vortex in the opposite direction is gradually formed under the vortex close to the mold; this is because as the increase of current, the increase of Lorentz force, and the axial and radial components of Lorentz force act on the melt in two directions. Figure 14 shows the influence of magnetic field frequency and current intensity on the velocity of melt flowing from the outlet of shunt plate to the edge of billet; it can be seen that at 120 mm away from the outlet of the shunt plate, with increasing frequency, the melt velocity increases at first, and then decreases. At the edge of the billet, it grows as the magnetic field intensity grows stronger; the skin effect of electromagnetic fields is mostly responsible for this. When the frequency of electromagnetic field increases, the skin depth of electromagnetic field decreases, and the area affected by Lorentz force decreases, resulting in uneven melt flow. The greater the current, the greater the magnetic field strength, and the greater the velocity in the melt flow process. Figure 15 depicts the cloud diagram of temperature distribution and liquid fraction distribution during the DC casting and LFEC casting process. It can be shown that when a magnetic field is applied, the sump depth becomes shallower, the mushy zone becomes wider, the temperature distribution in the melt moves upward as a whole, and an upward concave shape appears on the contour line with the temperature of 876.5 K, combined with the above analysis of the flow field; this is due to the formation of thermal convection between the melts due to the stirring of Lorentz force. Figure 16 demonstrates the temperature change along path 1; when a magnetic field is present, it can be seen that the temperature distribution near the mold is lower than that without magnetic field, which shows that electromagnetic field can effectively improve the DC casting process and accelerate the solidification of melt. Lorentz force has a positive effect on DC casting process. Figure 17 displays the cloud diagram of temperature and liquid fraction distribution under four electromagnetic frequencies; the change of frequency has little influence on the temperature distribution as shown in the graph. With the increase of frequency, the overall temperature moves down slightly, but the isoline with the temperature of 867 K moves down obviously, and the melt area with the liquid fraction of 0.5 to 0.63 becomes smaller with the increase of frequency; this is because the strong electromagnetic force at higher frequency mainly acts on the edge, and the effect on the internal high-temperature liquid region is weakened due to the decrease of skin effect. Figure 18 shows the melt temperature distribution and liquid fraction distribution under four current intensities; as the current intensity rises, it can be seen that the concave of isoline 867 K becomes more serious, while the isoline with temperature of 810 K moves up slightly. It shows that the magnetic field stirs the melt at the edge of the mold, and the greater the current, the stronger the stirring ability. The farther away from the mold, the smaller the influence of the magnetic field on the melt. When the temperature is 582 K, the influence of the current intensity on the temperature distribution can be ignored; this is because the melt solid fraction is large in the later stages of solidification, the fluidity is poor, and the influence of current stirring is small. As the current intensity rises, the melt area with liquid fraction of 0.5 to 0.63 increases, and the melt area with liquid fraction greater than 0.63 decreases; this also confirms that the melt flow near the mold is accelerated, the convective heat transfer is increased, the cooling effect is enhanced, and the melt solidifies faster. Figure 19 is the temperature change along path 3 under the action of different magnetic field frequencies and different With the increase in frequency, it can be shown that the central temperature of the ingot increases; this is because the frequency increases, the action range of magnetic field decreases, the melt solidification speed slows down, and the temperature increases. However, with an increase in current, the melt central temperature drops, because the Lorentz force increases as the current increases, while the skin depth of magnetic field remains unchanged, resulting in the accelerated solidification of the central melt and the decrease of temperature. Figure 20 shows the change of convective heat flux at the center of the outlet of the shunt plate and the edge of the ingot. With the increase of current, the intensity of magnetic field increases, Lorentz force strengthens the flow of melt, accelerates the convective heat transfer of melt, and increases the convective heat flux. Under a certain current, due to skin effect, when the magnetic field frequency is about 40 Hz, the convective heat flux is the largest.

Conclusion
In this research, COMSOL Multiphysics is used to create a two-dimensional axisymmetric model. The current intensity and magnetic field frequency on the temperature field, flow field, Lorentz force, and liquid fraction of rare earth magnesium alloy DC casting in the steady-state stage are analyzed. The following conclusions are drawn: 1. The usage of a magnetic field can accelerate the melt flow, thereby shorten the solidification time, significantly move the melt isotherm upward, reduce the temperature gradient, the sump depth is decreased by approximately 50 mm, and the vortex in the melt moves to the mold. 2. The increase of current intensity can increase the Lorentz force; when the current is increased from 80 to 250 A, the Lorentz force at the edge of the ingot is expanded by about 10 times. It also accelerates melt flow, thereby increasing convective heat transfer, and then reduces temperature gradient and the sump depth. 3. The increase of magnetic field frequency can increase the Lorentz force, but due to the existence of the skin effect, the range of magnetic field is reduced and the melt flow distribution is uneven.

Availability of data and material
The datasets used or analyzed during the current study are available from the corresponding author on reasonable request.
Code availability Not applicable.

Declarations
Ethics approval Not applicable.

Consent to participate Not applicable.
Consent for publication Not applicable.

Conflict of interest
The authors declare no competing interests.