Numerical modeling of casting process

Călin Vamoş
(original), metals and alloys casting, Numerical Modeling, paper

Abstract

Authors

C. Vamoș
-Tiberiu Popoviciu Institute of Numerical Analysis, Romanian Academy

V. Soporan

C. Pavai

Keywords

Cite this paper as:

C. Vamoş, C. Pavai, V. Soporan, Numerical modeling of casting process, Rev. Anal. Numér. Théor. Approx., 25 (1996) no. 1-2, 247-255.

References

see the expansion block below.

PDF

https://ictp.acad.ro/jnaat/journal/article/view/1996-vol25-nos1-2-art25/Vamos-Pavai-Soporan

About this paper

Journal

Rev. Anal. Numer. Theor. Approx.

Publisher Name
Romanian Academy
DOI

https://ictp.acad.ro/jnaat/journal/article/view/1996-vol25-nos1-2-art25

Print ISSN
Online ISSN

MR

?

ZBL

?

[1] R.A. Brown, S.H. Davies (editors), Free Boundaries in Viscous Flows, (1991), Springer Verlag, New York.

[2] D.L. Dwoyer, M.Y. Hussaini, R.G. Voigt (editors), Theoretical Approaches to Turbulence, (1985), Springer Verlag, New York.

[3] C.W. Hirt, B.D. Nichols, Volume of fluid (VOF) method for the dynamics of free boundaries, J. Comp. Phys., 39 (1981), pp. 201-225.

[4] K.H. Hoffmann, J. Sprekels (editors), Free Boundary Problems: Theory and Applications, (1990), Longman, New York.

[5] H.J. Lin, W.S. Hwang, Combined fluid flow and heat transfer analysis for the filling of castings AFS Transactions, 96 (1988), pp. 447-458.

[6] G.I. Marciuk, Method of Numerical analysis, (1983), Ed. Academiei, Bucureşti (in Romanian).

[7] Y. Ohtsuka, T. Ono, K. Mizuno, N. Matsubara, Computer simulation system of molten metal flow in diecasting, 56th World Foundry Congress, Osaka, 23-28 September 1990.

[8] T. Piwonka, V. Voller, L. Katgerman (editors), Proceedings of the Sixth International Conference on Modeling of Casting and Welding Processes, held in Palm Coast, Florida, March 21-26, 1993, A Publication of The Minerals, Metals and Materials Society, Printed in USA.

[9] M. Rappaz, M. Özgü, K. Mahin (editors), Proceedings of the Fifth International Conference on Modeling of Casting and Welding Processes, held in Davos, Switzerland, September 16-21, 1990, A Publication of The Minerals, Metals and Materials Society, Printed in USA.

[10] V. Soporan, V. Constantinescu, M. Crişan, Solidification of alloys-theoretical preliminaries (1995) Editura Transilvania, Cluj-Napoca (in Romanian).

[11] V. Soporan, V. Constantinescu, Modeling on a macrostructural level of the alloy solidification process, (1995), Dacia, Cluj-Napoca (in Romanian).

[12] R. S. Stoehr, C. Wang, Coupled heat transfer and fluid flow in the filling of castings, AFS Transactions, 96 (1988), pp. 733-741.

soon

Paper (preprint) in HTML form

jnaat,+Journal+manager,+1996-Vamos-Pavai-Soporan

NUMERICAL MODELING OF CASTING PROCESS

C. VAMOŞ, C. PAVAI, V. SOPORAN

(Cluj-Napoca)

1. INTRODUCTION

The numerical modeling of casting must take into account the physical phenomena occurring from the beginning of mold filling to the final casting. The mold filling with melt is modeled as a viscous flow with free boundary. Then the cooling phase follows, over which a distinct period is the solidification accompanied by the latent heat release. The shrinkage occurs during the cooling, leading to different defects of the casting.
The most difficult problem is the existence of two moving boundaries. The first one is the free boundary of the liquid metal during the mold filling. The second moving boundary is the solidification front. So far no general numerical methods to solve satisfactorily the free boundaries evolution have been elaborated [1], [4]. For casting modeling, specialized numerical models exist [8], [9].
In this paper we present a three-dimensional numerical model based on the SOLA-VOF technique [3]. We have brought some improvements to the evolution of the melt free boundary and to the shrinkage modeling. We have applied this numerical model to the steel ingot casting and the results obtained allow a new interpretation of the casting technology. In our approach the superficial phenomena are neglected.

2. FORMULATION OF THE MATHEMATICAL PROBLEM

The description of the physical phenomena implied into the casting process can be found in [10], [11]. In this section we present the mathematical problem which must be solved numerically.
There are three space domains in which different physical phenomena occur. The bounded domain of the mold D m ⊂ R 3 D m ⊂ R 3 D_(m)subR^(3)\mathscr{D}_{m} \subset \mathbf{R}^{3}Dm⊂R3 is constant in time. The bounded domain filled with liquid metal is denoted by D l ( t ) ⊂ R 3 D l ( t ) ⊂ R 3 D_(l)(t)subR^(3)\mathscr{D}_{l}(t) \subset \mathrm{R}^{3}Dl(t)⊂R3 and that filled with
solidified metal by D s ( t ) ⊂ R 3 D s ( t ) ⊂ R 3 D_(s)(t)subR^(3)\mathscr{D}_{s}(t) \subset \mathbf{R}^{3}Ds(t)⊂R3. These two domains depend on time. At the initial time there is no melt in the mold cavity, so D i ( 0 ) = ∅ D i ( 0 ) = ∅ D_(i)(0)=O/\mathscr{D}_{i}(0)=\varnothingDi(0)=∅ and implicitly D s ( 0 ) = ∅ D s ( 0 ) = ∅ D_(s)(0)=O/\mathscr{D}_{s}(0)=\varnothingDs(0)=∅. The rest of the space contains air under normal conditions of pressure and temperature. Since the density and specific heat of the air is much smaller than the density and specific heat of the metal and the mold material, we neglect the mechanical and thermal phenomena occurring in the mold environment. Thus the modeling is confined only to the three bounded domains defined above.
The phenomena occurring at the separating surfaces between these domains have a particular importance. If we denote by ∂ D m , ∂ D l ∂ D m , ∂ D l delD_(m),delD_(l)\partial \mathscr{D}_{m}, \partial \mathscr{D}_{l}∂Dm,∂Dl and ∂ D s ∂ D s delD_(s)\partial \mathscr{D}_{s}∂Ds the domains boundaries, then the separating surfaces are given by
(2.1) B l s ( t ) = ∂ D l ∩ ∂ D s ; B m ( t ) = ∂ Q ∩ ∂ Q m ; B m s ( t ) = ∂ R m ∩ ∂ R s (2.1) B l s ( t ) = ∂ D l ∩ ∂ D s ; B m ( t ) = ∂ Q ∩ ∂ Q m ; B m s ( t ) = ∂ R m ∩ ∂ R s {:(2.1)B_(ls)(t)=delD_(l)nn delD_(s);quadB_(m)(t)=delQnn delQ_(m);quadB_(ms)(t)=delR_(m)nn delR_(s):}\begin{equation*} \mathscr{B}_{l s}(t)=\partial \mathscr{D}_{l} \cap \partial \mathscr{D}_{s} ; \quad \mathscr{B}_{m}(t)=\partial \mathscr{Q} \cap \partial \mathscr{Q}_{m} ; \quad \mathscr{B}_{m s}(t)=\partial \mathscr{R}_{m} \cap \partial \mathscr{R}_{s} \tag{2.1} \end{equation*}(2.1)Bls(t)=∂Dl∩∂Ds;Bm(t)=∂Q∩∂Qm;Bms(t)=∂Rm∩∂Rs
Obviously these surfaces depend on time. The boundaries of the domains with the environment are
B 1 ( t ) = ∂ Q ∖ ( ∂ D n ∪ ∂ D s ) ; B s ( t ) = ∂ D s ∖ ( ∂ Q ∪ ∂ Q m ) ; (2.2) R m ( t ) = ∂ P m ∖ ( ∂ Q ∪ ∂ D s ) B 1 ( t ) = ∂ Q ∖ ∂ D n ∪ ∂ D s ; B s ( t ) = ∂ D s ∖ ∂ Q ∪ ∂ Q m ; (2.2) R m ( t ) = ∂ P m ∖ ∂ Q ∪ ∂ D s {:[B_(1)(t)=delQ\\(delD_(n)uu delD_(s));B_(s)(t)=delD_(s)\\(delQuu delQ_(m));],[(2.2)R_(m)(t)=delP_(m)\\(delQuu delD_(s))]:}\begin{gather*} \mathscr{B}_{1}(t)=\partial \mathscr{Q} \backslash\left(\partial \mathscr{D}_{n} \cup \partial \mathscr{D}_{s}\right) ; \mathscr{B}_{s}(t)=\partial \mathscr{D}_{s} \backslash\left(\partial \mathscr{Q} \cup \partial \mathscr{Q}_{m}\right) ; \\ \mathscr{R}_{m}(t)=\partial \mathscr{P}_{m} \backslash\left(\partial \mathscr{Q} \cup \partial \mathscr{D}_{s}\right) \tag{2.2} \end{gather*}B1(t)=∂Q∖(∂Dn∪∂Ds);Bs(t)=∂Ds∖(∂Q∪∂Qm);(2.2)Rm(t)=∂Pm∖(∂Q∪∂Ds)
The mechanical state at a point of radius vector x = ( x 1 , x 2 , x 3 ) x = x 1 , x 2 , x 3 x=(x_(1),x_(2),x_(3))\mathbf{x}=\left(x_{1}, x_{2}, x_{3}\right)x=(x1,x2,x3) at the time t t ttt is given by the velocity field v ( x , t ) = ( v 1 , v 2 , v 3 ) v ( x , t ) = v 1 , v 2 , v 3 v(x,t)=(v_(1),v_(2),v_(3))\mathbf{v}(\mathbf{x}, t)=\left(v_{1}, v_{2}, v_{3}\right)v(x,t)=(v1,v2,v3) and the pressure field p ( x , t ) p ( x , t ) p(x,t)p(\mathbf{x}, t)p(x,t). The velocity does not vanish and the pressure is different from the atmospheric pressure p 0 p 0 p_(0)p_{0}p0 only into the liquid metal. Then for x ∈ D l x ∈ D l xinD_(l)\mathrm{x} \in \mathscr{D}_{l}x∈Dl, we have the continuity equation
(2.3) div v = 0 (2.3) div v = 0 {:(2.3)divv=0:}\begin{equation*} \operatorname{div} \mathbf{v}=0 \tag{2.3} \end{equation*}(2.3)divv=0
and the Navier-Stokes equation
(2.4) ∂ v ∂ t + v ⋅ grad v = g − 1 ρ l grad p + v Δ v (2.4) ∂ v ∂ t + v ⋅ grad v = g − 1 ρ l grad p + v Δ v {:(2.4)(delv)/(del t)+v*gradv=g-(1)/(rho_(l))grad p+v Deltav:}\begin{equation*} \frac{\partial \mathbf{v}}{\partial t}+\mathbf{v} \cdot \operatorname{grad} \mathbf{v}=\mathbf{g}-\frac{1}{\rho_{l}} \operatorname{grad} p+v \Delta \mathbf{v} \tag{2.4} \end{equation*}(2.4)∂v∂t+v⋅gradv=g−1ρlgradp+vΔv
where g g ggg is the gravitational acceleration, ρ l ρ l rho_(l)\rho_{l}ρl is the melt density and v v vvv is the melt dynamical viscosity. The thermal state of the three domains is described by the temperature field T ( x , t ) T ( x , t ) T(x,t)T(\mathbf{x}, t)T(x,t) satisfying the Fourier equation
(2.5) ρ c [ ∂ T ∂ t + ( v ⋅ grad T ) ] = k Δ T + q (2.5) ρ c ∂ T ∂ t + ( v ⋅ grad T ) = k Δ T + q {:(2.5)rho c[(del T)/(del t)+(v*grad T)]=k Delta T+q:}\begin{equation*} \rho c\left[\frac{\partial T}{\partial t}+(\mathbf{v} \cdot \operatorname{grad} T)\right]=k \Delta T+q \tag{2.5} \end{equation*}(2.5)ρc[∂T∂t+(v⋅gradT)]=kΔT+q
where ρ ρ rho\rhoρ is the mass density, c c ccc is the specific heat and k k kkk the thermal conductivity of the material in the domain where (2.5) holds. The term q q qqq represents the heat sources per unite of mass and time such as the heat released by chemical reactions.
Equations (2.3) - (2.5) have to be completed by the constitutive relations. We use the simplest assumption that the material constants do not depend on temperature. Their values for each domain are denoted by the indexes: m m mmm for mold, l l lll
for the liquid metal, and s s sss for solidified metal. So we neglect both the convection due to the non-homogeneous heating of the melt and the change of the metal chemical composition on a side and the other of the solidification surface B l s B l s B_(ls)\mathscr{B}_{l s}Bls.
The initial and boundary conditions should be added to the system of equations (2.3)-(2.5). We have initial conditions only for the temperature field
(2.6) T ( x , 0 ) = T m for x ∈ D m (2.6) T ( x , 0 ) = T m  for  x ∈ D m {:(2.6)T(x","0)=T_(m)quad" for "quadxinD_(m):}\begin{equation*} T(\mathbf{x}, 0)=T_{m} \quad \text { for } \quad \mathbf{x} \in \mathscr{D}_{m} \tag{2.6} \end{equation*}(2.6)T(x,0)=Tm for x∈Dm
where T m T m T_(m)T_{m}Tm is the initial temperature of the mold. The condition at the boundary between the liquid metal and the mold wall is the usual one
(2.7) v ( x , t ) = 0 for x ∈ A n . (2.7) v ( x , t ) = 0  for  x ∈ A n .  {:(2.7)v(x","t)=0" for "xinA_(n)". ":}\begin{equation*} \mathbf{v}(\mathbf{x}, t)=0 \text { for } \mathbf{x} \in \mathscr{A}_{n} \text {. } \tag{2.7} \end{equation*}(2.7)v(x,t)=0 for x∈An. 
The solidification of the metal on the surface B l s B l s B_(ls)\mathscr{B}_{l s}Bls is accompanied by the increase of the mass density. So, there is a mass flux from the liquid state to the solid one through the surface B l s B l s B_(ls)\mathscr{B}_{l s}Bls implying the condition
(2.8) v ( x , t ) = V ( x , t ) for x ∈ B l s . (2.8) v ( x , t ) = V ( x , t )  for  x ∈ B l s .  {:(2.8)v(x","t)=V(x","t)" for "xinB_(ls)". ":}\begin{equation*} \mathrm{v}(\mathbf{x}, t)=\mathrm{V}(\mathbf{x}, t) \text { for } \mathbf{x} \in \mathscr{B}_{l s} \text {. } \tag{2.8} \end{equation*}(2.8)v(x,t)=V(x,t) for x∈Bls. 
The velocity V V VVV can be determined in terms of the displacement of B l s B l s B_(ls)\mathscr{B}_{l s}Bls but we shall not use this expression in the following. Condition (2.8) expresses the shrinkage process. The free boundary B l B l B_(l)\mathscr{B}_{l}Bl has two parts. One is the ingate of the melt S ⊂ B l S ⊂ B l SsubB_(l)\mathscr{S} \subset \mathscr{B}_{l}S⊂Bl, for which we have
(2,9) V ( x , t ) = { U ( t ) if t ≤ t ∗ 0 if t > t ∗ (2,9) V ( x , t ) = U ( t )  if  t ≤ t ∗ 0  if  t > t ∗ {:(2,9)V(x","t)={[U(t)," if "t <= t^(**)],[0," if "t > t^(**)]:}:}\mathbf{V}(x, t)= \begin{cases}\mathbf{U}(t) & \text { if } t \leq t^{*} \tag{2,9}\\ 0 & \text { if } t>t^{*}\end{cases}(2,9)V(x,t)={U(t) if t≤t∗0 if t>t∗
where U ( t ) U ( t ) U(t)\mathbf{U}(t)U(t) is a given function of time and represents the pouring velocity and t ∗ t ∗ t^(**)t^{*}t∗ is the filling time of the mold. For x ∈ B , ∖ S x ∈ B , ∖ S xinB,\\S\mathbf{x} \in \mathscr{B}, \backslash \mathscr{S}x∈B,∖S, the condition represents the contimuity of the interaction force and has the form
(2.10) p n i − ρ v ∑ k = 1 3 ( ∂ v i ∂ x k + ∂ v k ∂ x i ) n k = p 0 n i . (2.10) p n i − ρ v ∑ k = 1 3   ∂ v i ∂ x k + ∂ v k ∂ x i n k = p 0 n i . {:(2.10)pn_(i)-rho vsum_(k=1)^(3)((delv_(i))/(delx_(k))+(delv_(k))/(delx_(i)))n_(k)=p_(0)n_(i).:}\begin{equation*} p n_{i}-\rho v \sum_{k=1}^{3}\left(\frac{\partial v_{i}}{\partial x_{k}}+\frac{\partial v_{k}}{\partial x_{i}}\right) n_{k}=p_{0} n_{i} . \tag{2.10} \end{equation*}(2.10)pni−ρv∑k=13(∂vi∂xk+∂vk∂xi)nk=p0ni.
where n = ( n 1 , n 2 , n 3 ) n = n 1 , n 2 , n 3 n=(n_(1),n_(2),n_(3))\mathbf{n}=\left(n_{1}, n_{2}, n_{3}\right)n=(n1,n2,n3) is the exterior normal.
For the temperature field, the conditions on B s s , B m B s s , B m B_(ss),B_(m)\mathscr{B}_{s s}, \mathscr{B}_{m}Bss,Bm and B s m B s m B_(sm)\mathscr{B}_{s m}Bsm represent the equality of the heat fluxes
(2.11) k 1 ∂ T 1 ∂ n − k 2 ∂ T 2 ∂ n = q s (2.11) k 1 ∂ T 1 ∂ n − k 2 ∂ T 2 ∂ n = q s {:(2.11)k_(1)(delT_(1))/(del n)-k_(2)(delT_(2))/(del n)=q_(s):}\begin{equation*} k_{1} \frac{\partial T_{1}}{\partial n}-k_{2} \frac{\partial T_{2}}{\partial n}=q_{s} \tag{2.11} \end{equation*}(2.11)k1∂T1∂n−k2∂T2∂n=qs
where indexes refer to the media separated by the surface with the normal vector n n n\mathbf{n}n. The term q s q s q_(s)q_{s}qs is nonvanishing only for x ∈ B l s x ∈ B l s xinB_(ls)\mathbf{x} \in \mathscr{B}_{l s}x∈Bls where it represents the heat released per time unit and per surface unit by the metal solidification. It can be expressed in
terms of the displacement of B l s B l s B_(ls)\mathscr{B}_{l s}Bls but we shall not use this approach. In addition to (2.11), the continuity of temperature must also be imposed. For x ∈ B m ∪ B s ∪ ( B ∖ S ) x ∈ B m ∪ B s ∪ ( B ∖ S ) xinB_(m)uuB_(s)uu(B\\S)\mathbf{x} \in \mathscr{B}_{m} \cup \mathscr{B}_{s} \cup(\mathscr{B} \backslash \mathscr{S})x∈Bm∪Bs∪(B∖S), the boundary condition expresses the heat flux towards the environment
(2.12) k ∂ T ∂ n = − α ( T − T 0 ) (2.12) k ∂ T ∂ n = − α T − T 0 {:(2.12)k(del T)/(del n)=-alpha(T-T_(0)):}\begin{equation*} k \frac{\partial T}{\partial n}=-\alpha\left(T-T_{0}\right) \tag{2.12} \end{equation*}(2.12)k∂T∂n=−α(T−T0)
where α α alpha\alphaα is the global coefficient of the heat transfer and T 0 T 0 T_(0)T_{0}T0 is the air temperature. Condition (2.12) also holds for S S S\mathscr{S}S if t > t ∗ t > t ∗ t > t^(**)t>t^{*}t>t∗. During the casting filling
(2.13) T ( x , t ) = T p for x ∈ S and t ≤ t ∗ (2.13) T ( x , t ) = T p  for  x ∈ S  and  t ≤ t ∗ {:(2.13)T(x","t)=T_(p)quad" for "quadxinSquad" and "quad t <= t^(**):}\begin{equation*} T(\mathbf{x}, t)=T_{p} \quad \text { for } \quad \mathbf{x} \in \mathscr{S} \quad \text { and } \quad t \leq t^{*} \tag{2.13} \end{equation*}(2.13)T(x,t)=Tp for x∈S and t≤t∗
where T p T p T_(p)T_{p}Tp is the pouring temperature.
The separating surface B l s B l s B_(ls)\mathscr{B}_{l s}Bls between the liquid metal and the solidified one has a very intricate form. The phases are interpenetrated by lengthened formations with dimensions extended over some magnitude orders. The calculation of the shape of this surface in macroscopic volumes, like those of castings, goes beyond the capacity of the existing computers. In order to eliminate the microscopic details of B l s B l s B_(ls)\mathscr{B}_{l s}Bls we use a spatial averaging analogous to that used for turbulent flows [2]. If the microscopic structure of B l s B l s B_(ls)\mathscr{B}_{l s}Bls is homogeneous in the averaging volume, then a new variable representing the volume fraction occupied by the solidified metal is defined as
(2.14) s ¯ ( x , t ) = { 0 if T > T l T 1 − T ( x , t ) T l − T s if T s ≤ T ≤ T l 1 if T < T l (2.14) s ¯ ( x , t ) = 0  if  T > T l T 1 − T ( x , t ) T l − T s  if  T s ≤ T ≤ T l 1  if  T < T l {:(2.14) bar(s)(x","t)={[0," if ",T > T_(l)],[(T_(1)-T(x,t))/(T_(l)-T_(s))," if ",T_(s) <= T <= T_(l)],[1," if ",T < T_(l)]:}:}\bar{s}(\mathbf{x}, t)=\left\{\begin{array}{ccc} 0 & \text { if } & T>T_{l} \tag{2.14}\\ \frac{T_{1}-T(\mathbf{x}, t)}{T_{l}-T_{s}} & \text { if } & T_{s} \leq T \leq T_{l} \\ 1 & \text { if } & T<T_{l} \end{array}\right.(2.14)s¯(x,t)={0 if T>TlT1−T(x,t)Tl−Ts if Ts≤T≤Tl1 if T<Tl
where T l T l T_(l)T_{l}Tl is the liquidus temperature and T s T s T_(s)T_{s}Ts the solidus temperature. By averaging, equations (2.3) - (2.5) keep their form but they refer to averaged fields denoted by v ― , p ¯ v ¯ , p ¯ bar(v), bar(p)\overline{\mathrm{v}}, \bar{p}v―,p¯ and T ¯ T ¯ bar(T)\bar{T}T¯.
The averaging, also induces the change of the constitutive relations. For s ¯ ( x , t ) ∈ ( 0 , 1 ) s ¯ ( x , t ) ∈ ( 0 , 1 ) bar(s)(x,t)in(0,1)\bar{s}(\mathbf{x}, t) \in(0,1)s¯(x,t)∈(0,1), about the point x x x\mathbf{x}x there is a mixture oi liquid and solidified metal called mushy region. The mechanical properties of this mixture are complex but we have adopted the usual simple solution. If s ¯ ≥ 0.3 s ¯ ≥ 0.3 bar(s) >= 0.3\bar{s} \geq 0.3s¯≥0.3, the mixture is considered a rigid solid. For s ¯ < 0.3 s ¯ < 0.3 bar(s) < 0.3\bar{s}<0.3s¯<0.3 the melt may flow through the microscopic structures of the solidified metal and an additional drag is included in an averaged viscosity [7]
v ¯ = v ( 1 − s ¯ ) n v ¯ = v ( 1 − s ¯ ) n bar(v)=(v)/((1-( bar(s)))^(n))\bar{v}=\frac{v}{(1-\bar{s})^{n}}v¯=v(1−s¯)n
where n ∈ [ 5 , 30 ] n ∈ [ 5 , 30 ] n in[5,30]n \in[5,30]n∈[5,30]. It is assumed that the solidified metal does not participate to the flow and the density in the averaged equation (2.4) is equal to ρ l ρ l rho_(l)\rho_{l}ρl not to the mixture density. This interpretation also implies that the mean velocity v ― v ¯ bar(v)\overline{\mathbf{v}}v― is defined by
averaging only the melt velocities, without taking into account the solidified metal.
In the averaged equation (2.5), the material constants have averaged values. For example for mass density we have
ρ ¯ = ( 1 − s ¯ ) ρ l + s ¯ ρ s . ρ ¯ = ( 1 − s ¯ ) ρ l + s ¯ ρ s . bar(rho)=(1- bar(s))rho_(l)+ bar(s)rho_(s).\bar{\rho}=(1-\bar{s}) \rho_{l}+\bar{s} \rho_{s} .ρ¯=(1−s¯)ρl+s¯ρs.
The same expression is used for c ¯ c ¯ bar(c)\bar{c}c¯ and k ¯ k ¯ bar(k)\bar{k}k¯. By averaging the heat source q s q s q_(s)q_{s}qs in (2.11), it is changed into an internal heat source
(2.15) q s = ρ s L ∂ s ¯ ∂ t (2.15) q s = ρ s L ∂ s ¯ ∂ t {:(2.15)q_(s)=rho_(s)L(del( bar(s)))/(del t):}\begin{equation*} q_{s}=\rho_{s} L \frac{\partial \bar{s}}{\partial t} \tag{2.15} \end{equation*}(2.15)qs=ρsL∂s¯∂t
where L L LLL is the latent heat of fusion.
The initial condition (2.6) is not modified by averaging, but the boundary conditions are, because now the averaged domains are defined as
(2.16) D ― l ( t ) = { x ∈ D ∣ ⋃ D s ∣ s ( x , t ) < 0.3 } (2.16) D ¯ l ( t ) = x ∈ D ∣ ⋃ D s ∣ s ( x , t ) < 0.3 {:(2.16) bar(D)_(l)(t)={xinD_(∣)uuuD_(s)∣s(x,t) < 0.3}:}\begin{equation*} \overline{\mathscr{D}}_{l}(t)=\left\{\mathbf{x} \in \mathscr{D}_{\mid} \bigcup \mathscr{D}_{s} \mid s(\mathbf{x}, t)<0.3\right\} \tag{2.16} \end{equation*}(2.16)D―l(t)={x∈D∣⋃Ds∣s(x,t)<0.3}
D ― s ( t ) = { x ∈ Q ∪ D s ∣ s ( x , t ) ≥ 0.3 } . D ¯ s ( t ) = x ∈ Q ∪ D s ∣ s ( x , t ) ≥ 0.3 . bar(D)_(s)(t)={xinQuuD_(s)∣s(x,t) >= 0.3}.\overline{\mathscr{D}}_{s}(t)=\left\{\mathbf{x} \in \mathscr{Q} \cup \mathscr{D}_{s} \mid s(\mathbf{x}, t) \geq 0.3\right\} .D―s(t)={x∈Q∪Ds∣s(x,t)≥0.3}.
Using D ― l D ¯ l bar(D)_(l)\overline{\mathscr{D}}_{l}D―l and D ― s D ¯ s bar(D)_(s)\overline{\mathscr{D}}_{s}D―s in (2.1) and (2.2) instead of D l D l D_(l)\mathscr{D}_{l}Dl and D s D s D_(s)\mathscr{D}_{s}Ds the averaged boundaries are obtained. The only unmodified boundary is the ingate surface S S S\mathscr{S}S. Excepting some situations discussed below, the averaged conditions (2.7) - (2.13) have the same form if the averaged fields, domains and boundaries are used. The term q s q s q_(s)q_{s}qs in the averaged condition (2.11) vanishes. The condition (2.8) takes the form (2.7) and the shrinkage has to be modeled in a global way. The variation rate of the melt volume
(2.17) v ′ ( t ) = − ρ s ρ l d d t ∫ D ― l ∪ D ― s s ¯ ( x , t ) d x (2.17) v ′ ( t ) = − ρ s ρ l d d t ∫ D ¯ l ∪ D ¯ s   s ¯ ( x , t ) d x {:(2.17)v^(')(t)=-(rho_(s))/(rho_(l))((d))/((d)t)int_( bar(D)_(l)uu bar(D)_(s)) bar(s)(x","t)dx:}\begin{equation*} v^{\prime}(t)=-\frac{\rho_{s}}{\rho_{l}} \frac{\mathrm{~d}}{\mathrm{~d} t} \int_{\overline{\mathscr{D}}_{l} \cup \overline{\mathscr{D}}_{s}} \bar{s}(\mathbf{x}, t) \mathrm{d} \mathbf{x} \tag{2.17} \end{equation*}(2.17)v′(t)=−ρsρl d dt∫D―l∪D―ss¯(x,t)dx
is compensated by the displacement of B ― 1 B ¯ 1 bar(B)_(1)\overline{\mathscr{B}}_{1}B―1. Therefore we add to condition (2.10) a velocity normal to B ― l B ¯ l bar(B)_(l)\overline{\mathscr{B}}_{l}B―l, derived from (2.17).

3. THE NUMERICAL ALGORITHM

We use the finite difference technique applied on a staggered mesh rectangular cells C i j k C i j k C_(ijk)C_{i j k}Cijk. These cells are the averaging volumes defined in the previous section. Since in the following we shall refer only to the mean fields, we renounce to the bar above them.
According to VOF technique [3], an additional field f i j k f i j k f_(ijk)f_{i j k}fijk is introduced to describe the amount of metal in the cell C i j k C i j k C_(ijk)C_{i j k}Cijk. If f i j k = 1 f i j k = 1 f_(ijk)=1f_{i j k}=1fijk=1, then the cell is completely
filled and if f i j k = 0 f i j k = 0 f_(ijk)=0f_{i j k}=0fijk=0, the cell is empty. If 0 < f i j k < 1 0 < f i j k < 1 0 < f_(ijk) < 10<f_{i j k}<10<fijk<1, then the metal occupies a domain C i j k ′ ⊂ C i j k C i j k ′ ⊂ C i j k C_(ijk)^(')subC_(ijk)C_{i j k}^{\prime} \subset C_{i j k}Cijk′⊂Cijk and the cell must have an adjacent completely filled cell. The domain C i j k ′ C i j k ′ C_(ijk)^(')C_{i j k}^{\prime}Cijk′ has common faces with the adjacent filled cells and its thickness is constant. The volume of C i j k ′ C i j k ′ C_(ijk)^(')C_{i j k}^{\prime}Cijk′ is the f i j k f i j k f_(ijk)f_{i j k}fijk fraction of the cell volume. The field s i j k ∈ [ 0 , f i j k ] s i j k ∈ 0 , f i j k s_(ijk)in[0,f_(ijk)]s_{i j k} \in\left[0, f_{i j k}\right]sijk∈[0,fijk] represents the fraction of the cell volume occupied by soliditied metal. We also introduce a field o i j k o i j k o_(ijk)o_{i j k}oijk to describe if the material in C i j k C i j k C_(ijk)C_{i j k}Cijk is solid or fluid. If o i j k = 1 o i j k = 1 o_(ijk)=1o_{i j k}=1oijk=1, then the cell C i j k C i j k C_(ijk)C_{i j k}Cijk belongs to the mold or s i j k ≥ 0.3 f i j k s i j k ≥ 0.3 f i j k s_(ijk) >= 0.3f_(ijk)s_{i j k} \geq 0.3 f_{i j k}sijk≥0.3fijk. Otherwise o i j k = 0 o i j k = 0 o_(ijk)=0o_{i j k}=0oijk=0. Now we can define the discrete analogues of the averaged domains (2.16) and D m D m D_(m)\mathscr{D}_{m}Dm
D m = { x ∈ C i j k ∣ o i j k = 1 and f i j k = 0 } D m = x ∈ C i j k ∣ o i j k = 1  and  f i j k = 0 D_(m)={xinC_(ijk)∣o_(ijk)=1" and "f_(ijk)=0}D_{m}=\left\{\mathbf{x} \in C_{i j k} \mid o_{i j k}=1 \text { and } f_{i j k}=0\right\}Dm={x∈Cijk∣oijk=1 and fijk=0}
(3.1) D l = { x ∈ C i j k ∣ o i j k = 0 and f i j k = 1 } ∪ { x ∈ C i j k ′ ∣ o i j k = 0 and 0 < f i j k < 1 } D s = { x ∈ C i j k ∣ o i j k = 1 and f i j k = 1 } ∪ { x ∈ C i j k ′ ∣ o i j k = 1 and 0 < f i j k < 1 } (3.1) D l = x ∈ C i j k ∣ o i j k = 0  and  f i j k = 1 ∪ x ∈ C i j k ′ ∣ o i j k = 0  and  0 < f i j k < 1 D s = x ∈ C i j k ∣ o i j k = 1  and  f i j k = 1 ∪ x ∈ C i j k ′ ∣ o i j k = 1  and  0 < f i j k < 1 {:[(3.1)D_(l)={xinC_(ijk)∣o_(ijk)=0" and "f_(ijk)=1}uu{xinC_(ijk)^(')∣o_(ijk)=0" and "0 < f_(ijk) < 1}],[D_(s)={xinC_(ijk)∣o_(ijk)=1" and "f_(ijk)=1}uu{xinC_(ijk)^(')∣o_(ijk)=1" and "0 < f_(ijk) < 1}]:}\begin{align*} & D_{l}=\left\{\mathbf{x} \in C_{i j k} \mid o_{i j k}=0 \text { and } f_{i j k}=1\right\} \cup\left\{\mathbf{x} \in C_{i j k}^{\prime} \mid o_{i j k}=0 \text { and } 0<f_{i j k}<1\right\} \tag{3.1}\\ & D_{s}=\left\{\mathbf{x} \in C_{i j k} \mid o_{i j k}=1 \text { and } f_{i j k}=1\right\} \cup\left\{\mathbf{x} \in C_{i j k}^{\prime} \mid o_{i j k}=1 \text { and } 0<f_{i j k}<1\right\} \end{align*}(3.1)Dl={x∈Cijk∣oijk=0 and fijk=1}∪{x∈Cijk′∣oijk=0 and 0<fijk<1}Ds={x∈Cijk∣oijk=1 and fijk=1}∪{x∈Cijk′∣oijk=1 and 0<fijk<1}
Using these domains in (2.1) and (2.2) we obtain the corresponding boundaries. We note that the boundaries B n , B e s , B l m B n , B e s , B l m B_(n),B_(es),B_(lm)B_{n}, B_{e s}, B_{l m}Bn,Bes,Blm, and B s m B s m B_(sm)B_{s m}Bsm are formed only by faces of the cells C i j k C i j k C_(ijk)C_{i j k}Cijk, whereas B l B l B_(l)B_{l}Bl and B s B s B_(s)B_{s}Bs are formed by the internal faces of C ′ i j k C ′ i j k C^(')_(ijk)C^{\prime}{ }_{i j k}C′ijk.
The velocity components are defined at the center of the corresponding faces of the cells and the temperature and pressure at the center of the cells. To describe the pressure and velocity fields in all the mesh we have to use the averaged boundary conditions (2.7)-(2.10). The velocity components on the faces in B l s B l s B_(ls)B_{l s}Bls and B l m B l m B_(lm)B_{l m}Blm vanish according to (2.7) and the average of (2.8). The velocity components on the faces normal to B l s B l s B_(ls)B_{l s}Bls or B l m B l m B_(lm)B_{l m}Blm and included into D n D n D_(n)D_{n}Dn or D s D s D_(s)D_{s}Ds are equal to the components on the adjacent cell in D l D l D_(l)D_{l}Dl multiplied by numerical coefficient γ γ gamma\gammaγ [12]. The cells having a face included in the ingate surface have a special treatment. Not only the velocity components corresponding to faces in S S S\mathscr{S}S, but also those corresponding to the opposite faces are given by (2.9). In addition, the pressure of these cells is all the time equal to the atmospheric pressure p 0 p 0 p_(0)p_{0}p0 and their temperature is T p T p T_(p)T_{p}Tp for t ≤ t ∗ t ≤ t ∗ t <= t^(**)t \leq t^{*}t≤t∗, according to (2.13). Condition (2.10) has to be applied to obtain the velocity components on the faces of the cells with f i j k ∈ ( 0 , 1 ) f i j k ∈ ( 0 , 1 ) f_(ijk)in(0,1)f_{i j k} \in(0,1)fijk∈(0,1) which are not common with a completely filled cell. But the estimation of the normal vector n n n\mathbf{n}n is difficult in the general case and condition (2.10) has to be formulated separately for different types of flows.
For temporal integration we use the explicit finite difference approximation of Euler type. In the following we describe how the fields at time step n n nnn are obtained if the fields at the previous step n − 1 n − 1 n-1n-1n−1 are known. To simplify the formulae we consider the bidimensional case, i.e. v = ( u , 0 , w ) v = ( u , 0 , w ) v=(u,0,w)\mathrm{v}=(u, 0, w)v=(u,0,w) and the independent variables are denoted by x x xxx and z z zzz. We also assume that the cells are squares of side Δ l Δ l Delta l\Delta lΔl.
The Navier-Stokes equation (2.4) has been integrated using the "upstream" approximation [6]. The velocity field obtained does not satisfy the continuity equation (2.3). According to SOLA technique [3] the fields v ( n ) v ( n ) v^((n))\mathbf{v}^{(n)}v(n) and p ( n − 1 ) p ( n − 1 ) p^((n-1))p^{(n-1)}p(n−1) are iteratively modified for the cells with f i k = 1 f i k = 1 f_(ik)=1f_{i k}=1fik=1 such that the continuity equation (2.3) should be satisfied below an limit error ε 0 ε 0 epsi_(0)\varepsilon_{0}ε0. We denote by
(3.2) div i k ( n ) = 1 Δ l ( u i + 1 , k ( n ) − u i k ( n ) ) + 1 Δ l ( w i , k + 1 ( n ) − w i k ( n ) ) (3.2) div i k ( n ) = 1 Δ l u i + 1 , k ( n ) − u i k ( n ) + 1 Δ l w i , k + 1 ( n ) − w i k ( n ) {:(3.2)div_(ik)^((n))=(1)/(Delta l)(u_(i+1,k)^((n))-u_(ik)^((n)))+(1)/(Delta l)(w_(i,k+1)^((n))-w_(ik)^((n))):}\begin{equation*} \operatorname{div}_{i k}^{(n)}=\frac{1}{\Delta l}\left(u_{i+1, k}^{(n)}-u_{i k}^{(n)}\right)+\frac{1}{\Delta l}\left(w_{i, k+1}^{(n)}-w_{i k}^{(n)}\right) \tag{3.2} \end{equation*}(3.2)divik(n)=1Δl(ui+1,k(n)−uik(n))+1Δl(wi,k+1(n)−wik(n))
the divergence of the cell C i k C i k C_(ik)C_{i k}Cik. If there is at least a cell for which | div i k ( n ) | > ε 0 div i k ( n ) > ε 0 |div_(ik)^((n))| > epsi_(0)\left|\operatorname{div}_{i k}^{(n)}\right|>\varepsilon_{0}|divik(n)|>ε0, then for each cell which is not adjacent to a rigid obstacle the velocity components and the pressure are changed by the amounts
(3.3) δ u i + 1 , k ( n ) = − δ u i k ( n ) = δ w i , k + 1 ( n ) = − − δ w i k ( n ) = − 1 4 Δ l ⋅ div i k ( n ) δ p i k ( n ) = − 1 4 ρ l 4 Δ t Δ l 2 div i k ( n ) (3.3) δ u i + 1 , k ( n ) = − δ u i k ( n ) = δ w i , k + 1 ( n ) = − − δ w i k ( n ) = − 1 4 Δ l ⋅ div i k ( n ) δ p i k ( n ) = − 1 4 ρ l 4 Δ t Δ l 2 div i k ( n ) {:[(3.3)deltau_(i+1,k)^((n))=-deltau_(ik)^((n))=deltaw_(i,k+1)^((n))=--deltaw_(ik)^((n))=-(1)/(4)Delta l*div_(ik)^((n))],[deltap_(ik)^((n))=-(1)/(4)(rho_(l))/(4Delta t)Deltal^(2)div_(ik)^((n))]:}\begin{gather*} \delta u_{i+1, k}^{(n)}=-\delta u_{i k}^{(n)}=\delta w_{i, k+1}^{(n)}=--\delta w_{i k}^{(n)}=-\frac{1}{4} \Delta l \cdot \operatorname{div}_{i k}^{(n)} \tag{3.3}\\ \delta p_{i k}^{(n)}=-\frac{1}{4} \frac{\rho_{l}}{4 \Delta t} \Delta l^{2} \operatorname{div}_{i k}^{(n)} \end{gather*}(3.3)δui+1,k(n)=−δuik(n)=δwi,k+1(n)=−−δwik(n)=−14Δl⋅divik(n)δpik(n)=−14ρl4ΔtΔl2divik(n)
If C i k C i k C_(ik)C_{i k}Cik is adjacent to a rigid obstacle on which the velocity vanishes, then the variation of the corresponding velocity component is assigned to the opposite face of the cell. The iterations continue until all the cells with f i k = 1 f i k = 1 f_(ik)=1f_{i k}=1fik=1 have the divergence smaller than ε 0 ε 0 epsi_(0)\varepsilon_{0}ε0 or the number of iterations is greater than a given maximum value. In this way we obtain the solution of the Poisson equation for the pressure field. Summing up (3.2) on the cells with f i k = 1 f i k = 1 f_(ik)=1f_{i k}=1fik=1, we obtain
(3.4) Δ l 2 ∑ f = 1 div i k ( n ) = Δ l ( ∑ S + ∑ B l ) ( u i + 1 , k ( n ) − u i k ( n ) + w i , k + 1 ( n ) − w i k ( n ) ) (3.4) Δ l 2 ∑ f = 1   div i k ( n ) = Δ l ∑ S   + ∑ B l   u i + 1 , k ( n ) − u i k ( n ) + w i , k + 1 ( n ) − w i k ( n ) {:(3.4)Deltal^(2)sum_(f=1)div_(ik)^((n))=Delta l(sum_(S)+sum_(B_(l)))(u_(i+1,k)^((n))-u_(ik)^((n))+w_(i,k+1)^((n))-w_(ik)^((n))):}\begin{equation*} \Delta l^{2} \sum_{f=1} \operatorname{div}_{i k}^{(n)}=\Delta l\left(\sum_{\mathscr{S}}+\sum_{B_{l}}\right)\left(u_{i+1, k}^{(n)}-u_{i k}^{(n)}+w_{i, k+1}^{(n)}-w_{i k}^{(n)}\right) \tag{3.4} \end{equation*}(3.4)Δl2∑f=1divik(n)=Δl(∑S+∑Bl)(ui+1,k(n)−uik(n)+wi,k+1(n)−wik(n))
where the two sums in the right hand side refer to the velocity components normal to the cell faces forming the ingate surface S S S\mathscr{S}S and the free boundary B r B r B_(r)B_{r}Br. The terms corresponding to B l s B l s B_(ls)B_{l s}Bls and B l m B l m B_(lm)B_{l m}Blm are zero because of the condition at the solid boundaries. We multiply (3.4) by the time step Δ t Δ t Delta t\Delta tΔt. Then the sum for S S S\mathscr{S}S is the melt flux entering the mold during a time step. The sum for B l B l B_(l)B_{l}Bl represents the melt flux in the partially filled cells ( 0 < f < 1 ) ( 0 < f < 1 ) (0 < f < 1)(0<f<1)(0<f<1), i.e. it gives the displacement of the free boundary. Normally these two fluxes should be equal. But the divergences vanish only with the approximation ε 0 ε 0 epsi_(0)\varepsilon_{0}ε0. Therefore the left hand side of (3.4) is nonzero and implies a loss of mass. To eliminate this error, the left hand term is distributed in the sum for B l B l B_(l)B_{l}Bl by changing the velocity on the cells with f i k ( n − 1 ) < 1 f i k ( n − 1 ) < 1 f_(ik)^((n-1)) < 1f_{i k}^{(n-1)}<1fik(n−1)<1.
Using the velocity field v ( n ) v ( n ) v^((n))\mathbf{v}^{(n)}v(n) we compute the melt flux in the cells with f i k ( n − 1 ) < 1 f i k ( n − 1 ) < 1 f_(ik)^((n-1)) < 1f_{i k}^{(n-1)}<1fik(n−1)<1 and then the new values of f i k ( n ) f i k ( n ) f_(ik)^((n))f_{i k}^{(n)}fik(n). It is possible that f i k ( n ) < 0 f i k ( n ) < 0 f_(ik)^((n)) < 0f_{i k}^{(n)}<0fik(n)<0 or f i k ( n ) > 1 f i k ( n ) > 1 f_(ik)^((n)) > 1f_{i k}^{(n)}>1fik(n)>1. This means that the free boundary B l B l B_(l)B_{l}Bl is displaced from that cell and the field f i k ( ′ ′ ) f i k ( ′ ′ ) f_(ik)^((''))f_{i k}^{(\prime \prime)}fik(′′) must be modified. First the fluid in excess ( f − 1 ) ( f − 1 ) (f-1)(f-1)(f−1) is distributed from the cells with f > 1 f > 1 f > 1f>1f>1 to the adjacent cells with f < 1 f < 1 f < 1f<1f<1. Then the missing fluid in the cells with f < 0 f < 0 f < 0f<0f<0 is extracted from the completely filled adjacent cells. Finally the fluid in cells with f < 1 f < 1 f < 1f<1f<1 which have not adjacent completely filled cells is redistributed. This displacement technique of the free boundary conserves the mass but not the momentum and the energy. This is the reason why it is applicable only to the flows with streamlines almost normal to the free boundary, like the filling of the mold cavity.
The numerical computation of the temperature field is based on the equality of the heat fluxes on the cell faces [5]. According to (2.14) a new field s i k ( n ) s i k ( n ) s_(ik)^((n))s_{i k}^{(n)}sik(n) corresponds to T i k ( n ) T i k ( n ) T_(ik)^((n))T_{i k}^{(n)}Tik(n) and the volume of the liquid metal is diminished by (see (2.17)):
Δ v = ρ s ρ l Δ l 2 ∑ i , k ( s i k ( n − 1 ) − s i k ( n ) ) Δ v = ρ s ρ l Δ l 2 ∑ i , k   s i k ( n − 1 ) − s i k ( n ) Delta v=(rho_(s))/(rho_(l))Deltal^(2)sum_(i,k)(s_(ik)^((n-1))-s_(ik)^((n)))\Delta v=\frac{\rho_{s}}{\rho_{l}} \Delta l^{2} \sum_{i, k}\left(s_{i k}^{(n-1)}-s_{i k}^{(n)}\right)Δv=ρsρlΔl2∑i,k(sik(n−1)−sik(n))
This volume is eliminated from the cells with f i k ( n ) < 1 f i k ( n ) < 1 f_(ik)^((n)) < 1f_{i k}^{(n)}<1fik(n)<1 and s i k ( n ) < 0.3 f i k ( n ) s i k ( n ) < 0.3 f i k ( n ) s_(ik)^((n)) < 0.3f_(ik)^((n))s_{i k}^{(n)}<0.3 f_{i k}^{(n)}sik(n)<0.3fik(n). From each cell of this type a volume proportional to the mass in that cell ρ ( f i k ( n ) − s i k ( n ) ) + ρ s s i k ( n ) ρ f i k ( n ) − s i k ( n ) + ρ s s i k ( n ) rho(f_(ik)^((n))-s_(ik)^((n)))+rho_(s)s_(ik)^((n))\rho\left(f_{i k}^{(n)}-s_{i k}^{(n)}\right)+\rho_{s} s_{i k}^{(n)}ρ(fik(n)−sik(n))+ρssik(n) is eliminated.

4. A SIMULATION OF AN INGOT CASTING

We applied the numerical model described in the previous sections to steel ingots casting. We used a mesh of 7 × 7 × 20 7 × 7 × 20 7xx7xx207 \times 7 \times 207×7×20 cubical cells with side of 0.112 m . The ingot mold has the wall thickness equal to one cell. The domain D m D m D_(m)D_{m}Dm is formed by the cells in the mesh boundary without the cells of the superior side of the mesh and without the melt ingate cell C 4 , 4 , 1 C 4 , 4 , 1 C_(4,4,1)C_{4,4,1}C4,4,1 in the center of the inferior wall. The inferior face of this cell represents the ingate surface S S S\mathscr{S}S. The upwards directed pouring velocity was constant U = 0.098 m / s U = 0.098 m / s U=0.098m//sU=0.098 \mathrm{~m} / \mathrm{s}U=0.098 m/s, so that 17 rows of cells were filled in t ∗ = 9 t ∗ = 9 t^(**)=9t^{*}=9t∗=9 min. The free boundary B l B l B_(l)B_{l}Bl of the melt was quasihorizontal. Then condition (2.10) was verified taking the pressure and the velocity components continuous on B l [ 12 ] B l [ 12 ] B_(l)[12]B_{l}[12]Bl[12].
A vertical, upwards directed jet is a very unsteady flow, therefore the first row of cells was filled with melt at rest and the temperature varied only because of thermal conduction. Then the whole numerical procedure was used and a flow with a constant pattern during the process of filling occurred. It is a vortex comprising the whole volume of melt, formed by an ascending vertical stream in the middle of the ingot mold above the ingate and by four descending streams at the corners of the ingot mold.
The results obtained allow an effective optimization of the ingot casting technology. Usually the melt in the ingot mold is considered at rest, but our simulation shows that the temperature distribution is essentially affected by the melt flow. Using the numerical simulation, the optimum values of the technological parameters (the pouring velocity and temperature, the initial mold temperature, the variation of the pouring velocity, etc.) can be determined such that the ingot quality should be improved (e.g. the volume of the shrinkage fault). The use of the thermoreactive powder can also be improved such that a maximum amount from the released heat should remain at the ingot top. A detailed description of the simulation results and of their technological importance will be made in other articles.

REFERENCES

  1. RA. Brown, S.H. Davies (editors), Free Boundaries in Viscous Flows, (1991) Springer Verlag, New York.
  2. D.L. Dwoyer, M.Y. Hussaini, R.G. Voigt (editors), Theoretical Approaches to Turbulence, (1985) Springer Verlag, New York.
  3. C.W. Hirt, B.D. Nichols, Volume of stuid (VOF) method for the dynamics of free boundaries, J. Comp. Phys., 39 (1981) 201-225.
  4. K.H. Hoffinamn, J. Sprekels (editors), Free Boundary Problems: Theory and Applications, (1990) Longman, New York.
  5. H.J. Lin, W S. Hwang, Combined fluid flow and heat transfer analysis for the filling of castings, AFS Transactions, 96 (1988) 447-458.
  6. G.I. Marciuk, Methods of Numerical Analysis, (1983) Ed. Academiei, Bucuresti (in Romanian).
  7. Y. Ohtsuka, T. Ono, K. Mizuno, N. Matsubara, Computer simulation system of molten metal flow in diecasting, 57th World Foundry Congress, Osaka, 23-28 September 1990.
  8. T. Piwonka, V. Voller, L. Katgennan (editors), Proceedings of the Sixth International Conference on Modeling of Casting and Welding Processes, held in Palm Coast, Florida, March 21-26, 1993, A Publication of The Minerals, Metals and Materials Society, Printed in USA.
  9. M. Rappaz, M. Özgü, K. Mahin (editors), Proceedings of the Fifth International Conference on Modeling of Casting and Welding Processes, held in Davos, Switzerland, September 16-21,1990, A Publication of The Minerals, Metals and Materials Society, Printed in USA.
  10. V. Soporan, V. Constantinescu, M. Crişan, Solidification of alloys-theoretical preliminaries, (1995) Editura Transilvania, Cluj-Napoca (in Romanian).
  11. V. Soporan, V. Constantinescu, Modeling on a macrostructural level of the alloy solidification process, (1995) Dacia, Cluj-Napoca (in Romanian).
  12. R.A. Stoehr, C. Wang, Coupled heat transfer and fluid flow in the filling of castings, AFS Transactions, 96 (1988) 733-741.
Received 27. II. 1996
"T. Popoviciu" Institute of Numerical Analysis
Republicii 37, PO Box 68
3400 Cluj-Napoca, Romania
Center of Weather Forecasting
Traian Vuia 149
3400 Cluj-Napoca, Romania
Technical University of Cluj-Napoca
C. Daicoviciu 15
3400 Cluj-Napoca, Romania
1996

Related Posts