An Experimental and Numerical Study of the Free Surface in an Uphill Teeming Ingot Casting Process

The generation of nonmetallic inclusions due to mold flux entrapment is detrimental for the final steel quality in the uphill teeming process. Herein, physical and mathematical modeling using a full‐scale mold is conducted to investigate the change of the hump height and the surface height at the free surface during filling. Moreover, the Weber number is also studied to analyze the tendency of mold flux entrapment at the free surface. The results show that the Reynolds stress model can predict the hump height more accurately compared with the k − ε‐based turbulence models. In addition, the result shows that a lower filling rate can result in a lower hump height as well as less fluctuations of the hump and surface heights. Furthermore, the Weber number for the low flow rate case is lower than 12.3, which indicates a small or no risk for mold flux entrapment. The present findings suggest that a reduction in velocity in the mold entrance region are an effective means of reducing the risk of mold flux entrapment.


Introduction
Although the fraction of steel produced in the world by ingot casting was decreased to 3.5% in 2018, it still comprised 59 million metric tons, including about 565 000 metric tons in Sweden. [1] Ingot casting is still frequently used by the manufacturer of lowalloy steel and steel for special applications, which only can be produced by this process. Specific examples are high-carbonchromium-bearing steel, thick plates, seamless tubes, forgings bars, and wire rods. [2] There are two methods for the filling of ingots, namely the top pouring method and the bottom teeming method (which is also called the uphill teeming method). In the top pouring method, molten steel is more exposed to the surrounding air which can result in serious reoxidation problems. During pouring, the stream of molten steel impinges on the surface of the melt. Thus, this can bring mold powder and reoxidation products back into the bulk, which results in macroinclusions. [3] In addition, the metal which splashes on the mold wall will result in defects on the ingots. Overall, this makes the top pouring method nonsuitable for high-quality steels, compared with the uphill teeming method. [4] In the latter case, the liquid steel flows down into a trumpet and passes through a horizontal runner to a vertical runner located at the bottom of the mold. This reduces the molten steel exposure to air, the entrapment of mold powder, and the occurrence of splashing on the mold wall.
Nonmetallic inclusions pose a severe problem because they act as stress raisers in cast steels, which lower mechanical properties. [5] According to Thomas and Zhang, [3,5] the nonmetallic inclusions can be divided into endogenous and exogenous inclusions. Endogenous inclusions are mainly formed due to deoxidation products or precipitates during cooling. Typical deoxidation products are Al 2 O 3 for aluminum-killed steels, SiO 2 for silicon-killed steels, or MnO in steel-containing manganese. [6] Exogenous inclusions are formed due to chemical reactions and mechanical interactions of molten steel with its surroundings. These exogenous inclusions represent a very serious problem compared with endogenous inclusions, as they have larger sizes. This, in turn, causes harmful effects to the steel's mechanical properties.
Eriksson et al. [7] conducted experiments with and without the use of mold powder during the filling of a mold for the production of an aluminum-killed-bearing steel grade. The results showed that a large number of new types of inclusions containing MnO were formed along with more alumina inclusions when no mold powder was used compared with the case with the use of mold powder during the short filling time of the ingot. It can be concluded that reoxidation by air will result in more as well as new types of inclusions. Thus, the addition of mold powder at the surface of molten steel is important, as the mold flux protects molten steel from atmospheric oxidation and absorbs inclusions from molten steel. However, mold flux entrapment is a serious problem when adding mold powder. Zhang et al. [8] investigated the larger inclusions in plain-carbon steel ingots cast when using the uphill teeming method. The researchers found that almost 59% of the ingot inclusions were formed due to reoxidation, whereas the large-size inclusions contained O, Na, Al, Si, K, Mn, S, and Ca. The latter represents 8% of the total inclusions and the composition indicates that they originate from mold flux. Large inclusions have a larger negative effect on the fatigue life than the smaller ones. [2] Furthermore, they also cause inferior DOI: 10.1002/srin.201900609 The generation of nonmetallic inclusions due to mold flux entrapment is detrimental for the final steel quality in the uphill teeming process. Herein, physical and mathematical modeling using a full-scale mold is conducted to investigate the change of the hump height and the surface height at the free surface during filling. Moreover, the Weber number is also studied to analyze the tendency of mold flux entrapment at the free surface. The results show that the Reynolds stress model can predict the hump height more accurately compared with the k À ε-based turbulence models. In addition, the result shows that a lower filling rate can result in a lower hump height as well as less fluctuations of the hump and surface heights. Furthermore, the Weber number for the low flow rate case is lower than 12.3, which indicates a small or no risk for mold flux entrapment. The present findings suggest that a reduction in velocity in the mold entrance region are an effective means of reducing the risk of mold flux entrapment.
surface appearances, poor polishability, and reduced resistance to corrosion. [9] In addition, ingot casting is the last step in the steel-making process. Therefore, these inclusions are hard to remove from liquid steel, as the separation time is limited. Thus, many investigations [10][11][12][13] have been performed to find ways to eliminate or reduce the formation of inclusions that originate from the mold flux.
Eriksson et al. [14] made a numerical simulation to study the flow pattern during filling with different opening angles of the inlet. The results showed that a 25 angle of the inlet nozzle can reduce the risk of slag entrapments, leading to an almost flat rising surface and low surface velocities. Ragnarsson [12] also investigated the opening angle in an ingot casting process where both water model experiments and numerical simulations were performed with a large range of opening angles. Their findings indicated that the inlet angle of 5 used in the industry was a good alternative with respect to the inclusion removal. Also, Tan et al. [13] conducted mathematical simulations for an uphill teeming process during the initial filling stage. The researchers found that the hump height varies a lot during the initial filling stage, due to air entrapments in the trumpet. Recently, Yin et al. [15] studied the influence of three different angled runners on the air entrapment in an ingot casting process. The results indicated that the 30 angled runner results in a much more stable increase in hump height during the initial filling stage.
The use of swirling flows in the filling operation was investigated by some researchers due to the increased cross section and reduced axial velocity, which makes the flow into the mold more clam. [16] Both mathematical and physical modeling with the utilization of swirl blades were used by Hallgren et al. [10] Their study showed that the unevenness of flow velocity was decreased in the vertical runner when swirl was employed. Hallgren et al. [11] also expanded the work to include the mold. The results revealed that a calmer filling was obtained with the use of a swirl blade compared with the case without a swirl blade. Tan et al. [17] also conducted numerical simulations with the implementation of a swirl blade during the filling operation. A more comprehensive model (considering the trumpet, runner, and mold) was simulated in the researchers' work. The results induced that more calm filling with less fluctuations could be achieved in the mold with the implementation of swirl blades.
Most of the research so far focusing on hump height has been done using numerical simulations. Only Hallgren et al. [18] performed physical model experiments for the validation of the predictions. However, the system was simplified as it did not consider the trumpet and horizontal runner. In addition, the mold was simplified with a cylinder of 150 mm diameter and the flow was fixed at a 0.5 m s À1 velocity at the entrance of the vertical runner. This is a very small flow rate and a lower volume size for the mold compared with the real industrial system. Therefore, validations and verifications for a real industrial size of the ingot filling process are highly needed to improve the understanding on how mold flux entrapments occur in the ingot casting process.
In the present work a water model, which is widely implemented in the field of process metallurgy, was built to validate the numerical simulation results, focusing on the hump height. The water model consists of the runner and the mold, which was made according to the former work [13] based on real industrial dimensions. Specifically, velocities in the horizontal runner ranging from 1.05 to 0.5 m s À1 were studied for a runner and mold system. It should be noted that the trumpet was not included in this work. Also the hump height generated at the free surface was recorded using a camera. In addition, the experimental results have been used to validate the present numerical results, focusing on the formation of the hump height.

Experimental Section
To investigate the hump height and surface height inside the mold, a water model was built. The dimension of the mold is shown in Figure 1, where the dimension is taken from Tan. [13] However, the bottom of the mold is modified to be a flat surface due to restrictions during manufacturing. A flange is fixed at the bottom of the mold, which is used to tightly connect the runner. In addition, the mold is made by plexiglass and is transparent, to visualize the hump height formation inside the mold. Three holes are placed at the side wall of the mold, used as outlets. The function of the holes is to help the system reach a semisteady state, which is significantly useful for visualizing the hump height formation.
The whole water model experiment sketch is shown in Figure 2. The water from the reservoir was pumped into the horizontal runner. The flow meter was connected between the horizontal runner and the pump, which was used to control the volume flow rate of water to obtain the desired inlet velocity. The vertical runner and horizontal runner were fixed to the mold, where the connecting junction was made by 3D printing. Thereafter, the water was transported through the runner to the mold and out of the mold through open holes at the sidewall. Thus, a semisteady state can be reached using this system, when the difference of the volume flow rates between the inlet and the outlet approaches zero. In addition, some holes were plugged to investigate the hump height at three different surface heights. Moreover, a camera was placed at the side of the mold to record the hump height and surface height at a semisteady state.

Computational Domain
The computational domains for the present simulation are shown in Figure 3, i-iii. To save simulation time, three different computational domains were built based on the height of the holes. The computational domain was separated into two parts: the inner part, which is above the vertical runner and as is shown as a small cylinder in Figure 3, and the outer part, which represents the rest of the computational domain. The inner part is the location for the hump height formation, where a finer mesh was generated. The outer part has a more coarse mesh compared with the inner part. This can help decrease the total element size of the mesh, contributing to a reduction in simulation time. Three turbulence models (standard k À ε, realizable k À ε, and Reynolds stress turbulence models) were implemented in the present simulation.

Numerical Assumptions
The following assumptions were made in the present simulations: 1) The fluid behaves as an incompressible Newtonian fluid.
2) Heat transfer is not considered. 3) Chemical reactions are neglected. 4) Material properties are assumed to be constant. 5) No initial entrapment of gas occurs at the trumpet and runner system.

Governing Equations
The governing equation can be expressed in the following generic form according to our assumptions [19] ∂ðρΦÞ where ρ is density, Φ is the conserved property, t is time, u is the mean velocity vector, Γ is the diffusion coefficient, and S Φ is the source term. The various transport equations are summarized in Table 1.

Standard k À ε Model
The standard k À ε model, which was proposed by Launder and Spalding, [20] is a two-transport equation model based on the turbulence kinetic energy (k) and its dissipation rate ε.
The assumption for the standard k À ε model is that the molecular viscosity is not considered and that the flow is assumed to be fully turbulent. [21] Based on the numerical assumption earlier, the transport equation may be written as follows:  where ρ is the density, μ is the dynamic viscosity, and μ i (i ¼ x, y, z) is the mean velocity in the direction of x, y, and z. Furthermore, G k represents the generation of turbulence kinetic energy caused by mean velocity gradients. The turbulence viscosity μ k is calculated using the following equation

Realizable k À ε Model
The realizable k À ε model was proposed by Shih et al. [22] It satisfies some mathematical constraints with respect to Reynolds stress, which are more consistent with the physics of turbulent flows compared with the standard k À ε model. [21] The turbulent kinetic energy (k) transport equation is the same as the standard k À ε model. The transport equations for the dissipation rate, ε, are shown below. where where the model constants are the following: C 1ε ¼ 1.44,

Reynolds Stress Model
For the Reynolds stress model (RSM), the hypothesis of an isotropic eddy viscosity is abandoned. In addition, the Reynolds are utilized in the transport equation to close the Reynolds-averaged Navier-Stokes equation, together with an equation calculating the dissipation rate. It means that seven equations are used for a 3D simulation. The transport equations for the transport of Reynolds stress can be written as follows [21] ∂ ∂t ρu 0 The terms describing the local time derivative and convection are presented on the left-hand side, whereas the terms representing turbulence diffusion, molecular diffusion, stress production, pressure strain, dissipation, and production by system rotation are given on the right-hand side. Two new terms (pressure strain and rotation terms) are inserted compared with the transport equation for turbulent kinetic energy.

Volume of Fluid Model
Both the gas and liquid phases are considered in the present work. Thus, the volume of fluid (VOF) method, which was developed by Hirt and Nichols, [23] was implemented to track the gas/liquid interface. The properties of air and water are shown in Table 2. The time-dependent volume fraction of fluid, F, is governed by the following equation Air phase only 1 Liquid phase only 0 < F < 1 At the interphase (13) Table 2. Physical properties of air and water. [15] Density

Boundary Conditions
The inlet is located at the horizontal runner and a constant velocity was applied. A parameter study was conducted to study the influence of the velocity of the inlet on the filling performance. Specifically, the following velocities were tested, 0.5, 0.6, 0.7, 0.8, 0.9, and 1.05 m s, À1 for both simulations and experiments. The pressure outlet condition with a zero-gauge pressure was adopted at the top and the side of the mold, which is shown in Figure 2. In addition, a no-slip condition was applied at the walls. Finally, a standard wall function was applied for the boundary condition at the walls.

Solution Method
All simulations in the present work were done using the commercial computational fluid dynamics software ANSYS Fluent 19.3. The number of cells used in the simulations was determined by conducting a mesh sensitivity study, as described later. Moreover, the PRESTO! discretization method was applied to discretize pressure. The coupled scheme was used to solve the pressure-velocity coupling. The second-order upwind scheme was adopted to solve the momentum equations. Also, the second-order upwind schemes were used to calculate the turbulent kinetic energy and turbulent dissipation rate. In addition, the Courant number of VOF was fixed at 0.25. Finally, the transient calculations were done up to several seconds, until the mass flux of the liquid phase between the inlet and the outlet reached a semisteady state condition.

Mesh Sensitivity Study
Because the discretization error is influenced by the calculation grid size, [24,25] mesh sensitivity was studied with different grid sizes in the present work. Also, the grid convergence (GCI) method [26] was adopted to calculate numerical uncertainty. Mesh sensitivity evaluations were done in the present model with different grid sizes. Typically, the values of velocities are reported as a solution ϕ k for different grids sizes, when the simulations reach a steady-state solution. It should be pointed out that the present simulations do not reach a steady-state condition. Instead, it is assumed that a semi-steady state solution can represent a steady-state solution. Figure 4 shows the normalized mass flux of the liquid phase, which was calculated from Equation (14), for the case with a velocity of 1.05 m s À1 at the inlet.
where Q inlet is the mass flux of the liquid phase at the inlet and Q outlet is the mass flux of the liquid phase at the outlet. It is apparent that the time around 40 s can be stated as the semisteady state, as the normalized mass flux fluctuates around the value of zero. This, in turn, corresponds to the values of AE0.2. Table 3 shows the velocities in the fixed point (0, 0.022, 0) at different grid sizes of the mesh, where the velocities are extracted during several times over a 35 s period. Thereafter, an averaged value is determined. By substituting the velocity ϕ into Equation (15)- (20), the numerical uncertainty value GCI for a fine mesh was calculated as 1.39% in the current simulation. It shows a low numerical uncertainty value for the GCI evaluation method. Therefore, the fine mesh is adopted in the present simulation. It should be noted that only the GCI value for the fine mesh can be calculated using the above equations. It should be mentioned that the above calculation was done for the computational domain i as shown in Figure 3. In the present work, three computational domains were used in the simulations. Thus, the same calculations were also done for the computational domains of ii and iii. Consequently, the grid sizes of 398 343 and 428 106 were determined for the computational domains of ii and iii, respectively.    (20) where ΔV i is the volume, ε 32 ¼ ϕ 3 À ϕ 2 , ε 21 ¼ ϕ 2 À ϕ 1 and r ¼ h coarse h fine .

A Comparison of Turbulence Models for Numerical Modeling
The air and the liquid steel interfaces experienced a very complicated pattern during the 3D filling process. With a continuous supply of liquid steel from the bottom of the mold, a hump will be formed during the filling process. A higher height of the hump can result in a higher mold of flux entrapments. [11] The higher hump formation is due to strong teeming velocity at the bottom of the mold. Consequently, a strong turbulence can be generated at the interface, contributing to a large contact area between the liquid steel and the mold flux. This may cause mold flux entrapments. To model this process, one of the important issues is the prediction of the interface between the gas and the liquid phase. This can be reflected from the interface of water and air. To ensure that the results are accurate and trustworthy, the standard k À ε, realizable k À ε, and Reynolds stress turbulence models were validated using the physical modeling results. The interface between the air and liquid phase is shown in Figure 5, which represents a cross-section height plot for the free surface. Figure 5a shows the gas/liquid interface behavior for the standard k À ε model, ii for the realizable k À ε model, and iii for the RSM. Figure 5d shows the interface result from the cold model experiment. It should be mentioned that both the simulation and experimental results were performed using a 1.05 m s À1 inlet velocity, which was equal to a mass flow of 23 kg s À1 for a real industrial process. [13] It can be noticed that the gas and liquid interface outside the hump is almost flat for the predictions using the standard k À ε and realizable k À ε turbulence models. However, the predictions using the RSM turbulence model show a much more wavy liquid/gas interface. The predictions using the RSM turbulence model show a similar phenomenon that was found in the water model experiments. The hump height was also recorded for both simulations and water model experiments, the definition of the hump height is shown in Figure 6. The averaged hump height and its standard deviation are shown in Table 4. It should be noted that the water experiments were conducted for several minutes to reach a semisteady state. The similarity between the water model experiment and simulation is defined as follows.
where R experiment is the result for the water model experiment, R simulation is the result for the simulation. In the present work, R experiment and R simulation refer to the hump height and the  This indicates that there is significant hump height fluctuations, even when it has reached the semisteady state. However, the predictions using the standard k À ε and realizable k À ε model show a very low deviation of the hump height, which are 0 and 1.31 mm, respectively. This means that hump height is much more stable according to the predictions using the standard k À ε and the realizable k À ε turbulence models. The stable hump height can result in a lower turbulence generation at the gas/liquid interface, which can be an explanation for the absence of surface waves in the simulations when using the standard k À ε and realizable k À ε turbulence models. The deviation of the hump height for the case with RSM is larger, which contributes to higher fluctuations in hump height. This result agrees better with the results from the water model experiments. In addition, the averaged hump height is also different between the three turbulence models. The largest difference of the hump height is found between the realizable k À ε turbulence model and the experiments. However, the RSM and the standard k À ε models show an 88% agreement compared with the water model experiments according to Equation (21). Considering the deviation of hump height, the realizable k À ε and standard k À ε turbulence model predictions have a large difference compared with the experimental results. Specifically, from the experimental results, they have a 16% and 0% agreement with the experimental results for the standard k À ε turbulence model and the realizable k À ε turbulence model, respectively.
It should be noted that the above simulation results for the three turbulence models were simulated using an inlet velocity of 1.05 m s À1 , which was calculated from the industrial filling rate. [13] To determine the best turbulence model in the present work, the range of the inlet velocities should be enlarged. However, the present simulations are very time-consuming. Specifically, almost 20 days were needed for each parameter study, when the simulation reached a semisteady state. Therefore, only the realizable k À ε and the RSM models were tested in the present work. The main reason for choosing the realizable k À ε model is that its turbulence model was widely used for numerical simulations during the filling of the ingot casting process by many researchers. [10,11,13,18] Thus, it is of great interest to include the realizable k À ε turbulence model in the comparison. Figure 7 shows the hump height variation with different inlet velocities for the RSM and realizable k À ε models as well as for the experimental results.
It is clearly shown that the RSM model correlates best to the experimentally determined hump height. Furthermore, the realizable k À ε model has the largest difference (only about 1% agreement with the experimental results) when looking at the dynamics of the hump. It is shown that the k À ε models predict small and none fluctuations. This is in disagreement with the experimental results. Also, the RSM model gives a better prediction, showing the dynamic surface. Given the above results, the RSM model is recommended to be used for numerical treatment of the current system. Note that while the qualitative behavior is correct, there is no perfect quantitative match, especially, for higher flow rates.

Hump Height and Surface Height
Ragnarsson and Sichen [27] studied the generation of inclusions by mold powder in an ingot casting process. The results showed that the sizes of inclusions, which were generated due to mold flux entrapments, were in the range of 10-30 μm and they were solid at typical casting temperatures. These solid inclusions along   with their composition indicated that the mold powder was entrained into the liquid metal before the liquid slag layer was formed. In addition, Hallgren et al. [11] found that a hump was indeed formed in the plant trials. The hump showed an intensive interaction with the mold powder as well as more spreading of the mold powder during the initial filling stage compared with the later filling stage. The interaction will lead to reoxidation by strong oxidizers such as aluminum in liquid steel, which reacts with FeO inside the mold powder. Thus, it is of great interest to optimize the hump height during the initial filling stage, to minimize the number of the nonmetallic inclusions formed due to reactions with the casting powder. Figure 6 shows the definition of the hump height H, surface height L, and the height difference ΔH.
A parameter study of the filling conditions of an ingot casting process has been conducted. The effect of the inlet velocities on the filling was studied, to find the relation between the hump height and the inlet velocities. In this work, a parameter study of the velocity was made, where the velocities at the inlet in each case was fixed to one of the following values: 0.5, 0.6, 0.7, 0.8, 0.9, and 1.05 m s À1 . A 1.05 m s À1 inlet velocity corresponds to the 23 kg s À1 teeming speed used in industrial uphill teeming processes. The influence of inlet velocity on hump height is evaluated in the present work. It should be noticed that the hump height formation created when using the six studied inlet velocities was measured at the first hole, which is outlet 1 shown in Figure 2. The recorded averaged hump heights from the parameter study are shown in Figure 8. It is clearly shown that the hump height from the simulation has a 94% agreement compared with the experimental results. Typically, it reveals that the averaged hump height sharply decreases as the inlet velocity decreases. This is because the lower inlet velocity results in a lower momentum of the liquid phase at the entrance of the mold. Therefore, a lower hump height is formed at the gas/liquid interface. In addition, the deviation shows a decreasing trend as the inlet velocity is decreased. This indicates that higher velocities can result in higher fluctuations of the hump height, which contribute to a strong interaction with the mold powder. From the results in a previous study, [15] it can be explained that the fluctuations of the hump height are believed to be caused by the unevenness of the velocity in the vertical runner. Specifically, a higher velocity was formed at the left-hand side and lower velocity was formed at the right-hand side of the vertical runner.
To determine the interface of the gas/liquid outside the hump, the surface height was evaluated. The method is presented in the Experimental Section. Figure 9 shows the surface height variation with the inlet velocity for both simulations and experimental results. It can be noticed that the surface height decreases as the inlet velocity decreases. This is due to the smaller mass flow rate with filling time at the inlet. Moreover, it is also shown that the fluctuations of the surface height increase as the inlet velocity increases. It reveals that the surface height fluctuates more intensively at higher inlet speeds. It means that a higher fluctuated gas/liquid interface is formed outside the hump. For higher inlet velocities, larger fluctuations of the hump and surface heights are formed. This may cause mold flux entrapments. Overall, it can be concluded that a lower velocity can result in a lower hump height as well as lower fluctuations of the hump height and surface height. This can intensively decrease the mold flux entrapments. As a consequence, a lower number of nonmetallic inclusions generate during the filling stage of the ingot casting process.

Height Difference
It is also of great interest to investigate the hump height at the later stages of the filling process, to know the mold flux entrapments. To have a good comparison of the hump, the height difference (ΔH) was recorded instead of hump height. The surface height is different during filling, which results in a large difference in the hump height. In this work, the height difference was measured at a semiequilibrium state using three The parameter study of inlet velocities was carried out using the following values: 0.5, 0.8, and 1.05 m s À1 . Figure 10, 11, and 12 show the height differences at three opening holes with different inlet velocities. Considering the height difference, it can be seen that the simulation results have an almost 93% agreement with the experimental results. It appears that the height difference decreases with an increased teeming time. This is due to more water being present in the mold, which can act as flow resistance for the incoming water. Also, due to the fact that kinetic energy is transferred to a gravitational potential energy with the increased surface height. In addition, the case with an inlet velocity of 1.05 m s À1 shows a larger standard deviation compared with the other cases. The reason is that the higher velocity results in higher fluctuations of the hump height. Also, it can be noticed that the standard deviation of the height difference decreases with an increased teeming time. It can be concluded that the free surface becomes more stable as the filling proceeds.

Tendency of Mold Powder Entrapments
The tendency for mold flux entrapment can be estimated by using the Weber number. The Weber number is reported as follows [4] We ¼ u 2 steel · ρ steel ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi γgðρ steel À ρ slag Þ q (22) where u steel is the tangential steel velocity in the free surface, ρ steel is the density of steel, γ is the interfacial tension, and g is gravity.
The Weber number has been investigated by Xiao et al. [28] who developed both a physical model and a numerical model to determine the Weber number when oil drops were entrapped into a water bath. They found that oil drops entrained into water when the Weber number exceeded a value of 12.3. They also reported that slag will be entrained into liquid steel when the Weber number is higher than 12.3. According to these results, there is no risk for mold flux entrapments during ingot filling when the Weber number is lower than 12.3. It should be noted that the interfacial tension varies between molten steel and slag during ingot filling. Specifically, the interfacial tension decreases sharply in the beginning, due to intensive mass transfer. Thereafter, the value increases when mass transfer slows down and approaches a stable value when equilibrium is reached. [29,30] In the present calculation, the interfacial tension values were fixed as 0.5 N m À1 , [31] 1 N m À1 , [32] 1.5 N m À1 , [33] and 2 N m -1 . [4] In addition, the tangential velocity has a range at the free surface, where the maximum value is reported in Table 5. It should be noted that the simulations were implemented using liquid steel instead of water, to predict the tendency of the mold flux entrapments. The density for liquid steel is set as 6900 kg m À3 and the density for mold flux is defined as 2500 kg m À3 . [4] first hole second hole third hole   Figure 13 shows the Weber number as a function of velocity for different interfacial tension values between slag and molten steel. It can be noticed that a higher value of interfacial tension (2 N m À1 ) can result in a lower increase in the Weber number. In other words, if the velocity at the free surface is fixed, a higher interfacial tension can result in a lower tendency for mold flux entrapment. Considering the range of the tangential velocities in the different cases, the Weber number for the case with an inlet velocity of 0.5 m s À1 is lower than 12.3 at all tested values of interfacial tension. This reveals that there is no risk for mold flux entrapment for the case with an inlet velocity of 0.5 m s À1 . When the inlet velocity increases, the tendency of the mold flux entrapment increases dramatically. For an inlet velocity of 1.05 m s À1 , the Weber number exceeds the critical value of 12.3 for most of the velocity range. It means that the case with a 1.05 m s À1 velocity has a higher risk for mold flux entrapment compared with the other cases. From the present findings, it recommended that the industry should lower the filling velocity as much as possible. Specifically, a filling velocity of 0.5 m s À1 is a critical value to limit the mold flux entrapments, based on the results from the present setup. To maintain the same casting mass flow rate while decreasing the inlet velocity, it is necessary to modify the process geometry. For example, the angled runner [15] can be used to reduce the initial filling speeds at the mold entrance.

Conclusions
An experimental and numerical study of the variation of hump height during the filling stage of the uphill teeming process was performed. In addition, the Weber number was also used to estimate the tendency for mold flux entrapment. Both the physical and mathematical model results confirm that lower filling speeds can result in a lower risk of mold flux entrapment. Overall, the following main specific conclusions were obtained: 1) For the casting filling speed of 1.05 m s À1 (it correlates to a 23 kg s À1 mass flow of molten steel in the industrial process), the k À ε based turbulence model predictions show a poor agreement (16% and 0% for the realizable and standard k À ε mode, respectively) to water model results, considering the dynamics of the hump. However, the Reynolds stress turbulent model was found to be more accurate (88%) for predicting hump height and free surface dynamics during the filling process.
2) The hump height and surface height were both determined from the physical model and simulation results. The numerical simulations showed a 93% agreement with the experimental results, considering the height difference. Also, it was found that a lower inlet velocity can result in a lower hump height and lower fluctuations of the hump height. 3) In the current setup, the Weber number was found to increase with an increased filling speed, when the interfacial energy is fixed. The use of a 0.5 m s À1 inlet velocity resulted in no risks for mold flux entrapment for all tested interfacial tension values.