Phase averaging model for solidification of binary solution

From Thermal-FluidsPedia

Jump to: navigation, search

Since both liquid and solid are present in the mushy zone, the governing equations in the mushy zone can be obtained by performing volume averaging to each phase to yield a volume-averaging model. The volume-averaged multi-fluid modelis also applicable to alloy solidification. Considering the intrinsic phase-averaged density {{\left\langle {{\rho }_{k}} \right\rangle }^{k}} is the same as the density ρk, assuming that the liquid and solid phases are incompressible, and there are no internal heat generations in either phases, and neglecting viscous dissipation, the volume averaged governing equations become

\frac{\partial }{\partial t}\left( {{\varepsilon }_{k}}{{\rho }_{k}} \right)+\nabla \cdot \left( {{\varepsilon }_{k}}{{\rho }_{k}}{{\left\langle {{\mathbf{V}}_{k}} \right\rangle }^{k}} \right)=\sum\limits_{j=1(j\ne k)}^{\Pi }{{{{{\dot{m}}'''}}_{jk}}}\qquad\qquad(1)
  & \frac{\partial }{\partial t}\left( {{\varepsilon }_{k}}{{\rho }_{k}}{{\left\langle {{\mathbf{V}}_{k}} \right\rangle }^{k}} \right)+\nabla \cdot \left( {{\varepsilon }_{k}}{{\rho }_{k}}{{\left\langle {{\mathbf{V}}_{k}}{{\mathbf{V}}_{k}} \right\rangle }^{k}} \right) \\ 
 & =\nabla \cdot \left( {{\varepsilon }_{k}}{{\left\langle {{\mathbf{\tau }}_{k}} \right\rangle }^{k}} \right)+{{\varepsilon }_{k}}{{\rho }_{k}}{{\mathbf{X}}_{k}}+\sum\limits_{j=1(j\ne k)}^{\Pi }{\left( \left\langle {{\mathbf{F}}_{jk}} \right\rangle +\left\langle {{{{\dot{m}}'''}}_{jk}} \right\rangle {{\left\langle {{\mathbf{V}}_{k,I}} \right\rangle }^{k}} \right)} \\ 
\frac{\partial }{\partial t}\left( {{\varepsilon }_{k}}{{\rho }_{k}}{{\left\langle {{h}_{k}} \right\rangle }^{k}} \right)+\nabla \cdot \left( {{\varepsilon }_{k}}{{\rho }_{k}}{{\left\langle {{\mathbf{V}}_{k}}{{h}_{k}} \right\rangle }^{k}} \right)=-\nabla \cdot \left\langle {{{\mathbf{{q}''}}}_{k}} \right\rangle +\sum\limits_{j=1(j\ne k)}^{\Pi }{\left[ \left\langle {{{{q}'''}}_{jk}} \right\rangle +{{{{\dot{m}}'''}}_{jk}}{{\left\langle {{h}_{k,I}} \right\rangle }^{k}} \right]}\qquad\qquad(3)
\frac{\partial }{\partial t}\left( {{\varepsilon }_{k}}{{\rho }_{k}}{{\left\langle {{\omega }_{k,i}} \right\rangle }^{k}} \right)+\nabla \cdot \left( {{\varepsilon }_{k}}{{\rho }_{k}}{{\left\langle {{\omega }_{k,i}} \right\rangle }^{k}}{{\left\langle {{\mathbf{V}}_{k}} \right\rangle }^{k}} \right)=-\nabla \cdot \left\langle {{{\mathbf{{\dot{m}}''}}}_{k,i}} \right\rangle +\left\langle {{{{\dot{m}}'''}}_{k,i}} \right\rangle +\sum\limits_{j=1(j\ne k)}^{\Pi }{{{{{\dot{m}}'''}}_{jk,i}}}\qquad\qquad(4)

where k can be either (s) or liquid (\ell ), {{\left\langle {{\mathbf{V}}_{k,I}} \right\rangle }^{k}} and {{\left\langle {{h}_{k,I}} \right\rangle }^{k}} are the intrinsic phase-averaged velocity and enthalpy of the kth phase at the interface.

The above equations are generalized and applicable to any two-phase system. In the alloy solidification process, the solid phase can be assumed to be rigid and stationary as in purely columnar growth (Beckermann and Viskanta, 1993), i.e., {{\left\langle {{\mathbf{V}}_{s}} \right\rangle }^{s}}=0. The continuity equations in the solid and liquid phases become

\frac{\partial }{\partial t}\left( {{\varepsilon }_{s}}{{\rho }_{s}} \right)={{{\dot{m}}'''}_{\ell s}}\qquad\qquad(5)
\frac{\partial }{\partial t}\left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }} \right)+\nabla \cdot \left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }} \right)={{{\dot{m}}'''}_{\ell s}}\qquad\qquad(6)


{{\varepsilon }_{s}}+{{\varepsilon }_{\ell }}=1\qquad\qquad(7)
{{{\dot{m}}'''}_{\ell s}}+{{{\dot{m}}'''}_{s\ell }}=0\qquad\qquad(8)

The momentum equation for the solid phase is not needed, and the momentum equation for the liquid phase is

  & \frac{\partial }{\partial t}\left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }} \right)+\nabla \cdot \left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }} \right) \\ 
 & =\nabla \cdot \left( {{\varepsilon }_{\ell }}{{\left\langle {{\mathbf{\tau }}_{\ell }} \right\rangle }^{\ell }} \right)+{{\varepsilon }_{\ell }}{{\rho }_{\ell }}\mathbf{g}+\left\langle {{\mathbf{F}}_{s\ell }} \right\rangle +\left\langle {{{{\dot{m}}'''}}_{s\ell }} \right\rangle {{\left\langle {{\mathbf{V}}_{\ell ,I}} \right\rangle }^{\ell }}-\nabla \cdot \left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{{\mathbf{\hat{V}}}}_{\ell }}{{{\mathbf{\hat{V}}}}_{\ell }} \right\rangle }^{\ell }} \right) \\ 

where {{\left\langle {{\mathbf{V}}_{\ell }}{{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }} in eq. (2) is replaced by

{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }}+\left\langle {{{\mathbf{\hat{V}}}}_{\ell }}{{{\mathbf{\hat{V}}}}_{\ell }} \right\rangle

and {{\mathbf{X}}_{\ell }} is replaced by g, i.e., gravity is the only body force. Since the liquid velocity in the mushy zone is usually laminar because of the characteristic diameter, the product of the velocity deviations, \left\langle {{{\mathbf{\hat{V}}}}_{\ell }}{{{\mathbf{\hat{V}}}}_{\ell }} \right\rangle , is also small and can be neglected. {{\left\langle {{\mathbf{V}}_{\ell ,I}} \right\rangle }^{\ell }} is the intrinsic phase-averaged velocity of the kth phase at the interface. The difference between two adjacent phases results solely from the density difference between the two phases. Since the density change is usually very small during solid-liquid phase change ({{\rho }_{\ell }}\approx {{\rho }_{s}}), the intrinsic phase-averaged velocities for both phases at the solid-liquid liquid interface are approximately equal to each other, i.e., {{\left\langle {{\mathbf{V}}_{\ell ,I}} \right\rangle }^{\ell }}\approx {{\left\langle {{\mathbf{V}}_{s,I}} \right\rangle }^{s}}. The solid phase is rigid and stationary, thus the velocity of the solid phase at the interface, {{\left\langle {{\mathbf{V}}_{s,I}} \right\rangle }^{s}}=0. The fourth term on the right hand side of eq. (9) can therefore be neglected. The momentum equation for the liquid phase becomes

\frac{\partial }{\partial t}\left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }} \right)+\nabla \cdot \left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }} \right)=\nabla \cdot \left( {{\varepsilon }_{\ell }}{{\left\langle {{\mathbf{\tau }}_{\ell }} \right\rangle }^{\ell }} \right)+{{\varepsilon }_{\ell }}{{\rho }_{\ell }}\mathbf{g}+\left\langle {{\mathbf{F}}_{s\ell }} \right\rangle \qquad\qquad(10)

Applying the stress-strain relationship for Newton’s fluid for the first term on the right-hand side and replacing the third term by the Darcy-Brinkman’s relation (see Section 4.6.2), the momentum equation becomes

  & {{\rho }_{\ell }}\left[ \frac{1}{\varepsilon }\frac{\partial {{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }}}{\partial t}+\frac{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }{{{\varepsilon }^{2}}}\nabla \cdot \left\langle {{\mathbf{V}}_{\ell }} \right\rangle  \right]=-\nabla p+ \\ 
 & \frac{{{\mu }_{\ell }}}{\varepsilon }{{\nabla }^{2}}\left\langle {{\mathbf{V}}_{\ell }} \right\rangle +{{\rho }_{\ell }}\mathbf{g}-\frac{{{\mu }_{\ell }}}{\mathbf{K}}\left\langle {{\mathbf{V}}_{\ell }} \right\rangle -\frac{{{\rho }_{\ell }}{{C}_{\ell }}}{{{\mathbf{K}}^{{1}/{2}\;}}}\left| \left\langle {{\mathbf{V}}_{\ell }} \right\rangle  \right|\left\langle {{\mathbf{V}}_{\ell }} \right\rangle  \\ 

where \left\langle {{\mathbf{V}}_{\ell }} \right\rangle ={{\varepsilon }_{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }} is extrinsic phase-averaged velocity of the liquid phase.

The energy equations for the solid and liquid phases are

\frac{\partial }{\partial t}\left( {{\varepsilon }_{s}}{{\rho }_{s}}{{\left\langle {{h}_{s}} \right\rangle }^{s}} \right)=-\nabla \cdot \left\langle {{{\mathbf{{q}''}}}_{s}} \right\rangle +\left\langle {{{{q}'''}}_{s\ell }} \right\rangle +{{{\dot{m}}'''}_{\ell s}}{{\left\langle {{h}_{s,I}} \right\rangle }^{s}}\qquad\qquad(12)
\frac{\partial }{\partial t}\left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{h}_{\ell }} \right\rangle }^{\ell }} \right)+\nabla \cdot \left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }}{{h}_{\ell }} \right\rangle }^{\ell }} \right)=-\nabla \cdot \left\langle {{{\mathbf{{q}''}}}_{\ell }} \right\rangle +\left\langle {{{{q}'''}}_{s\ell }} \right\rangle +{{{\dot{m}}'''}_{s\ell }}{{\left\langle {{h}_{\ell ,I}} \right\rangle }^{\ell }}\qquad\qquad(13)

Adding eqs. (12) and (13) together, considering {{\left\langle {{h}_{\ell }} \right\rangle }^{\ell }}-{{\left\langle {{h}_{s}} \right\rangle }^{s}}={{h}_{s\ell }}, and assuming the solid and liquid phases are in thermal equilibrium, ({{\left\langle {{T}_{s}} \right\rangle }^{s}}={{\left\langle {{T}_{\ell }} \right\rangle }^{\ell }}=T), the energy equation becomes

\left\langle \rho {{c}_{p}} \right\rangle \frac{\partial T}{\partial t}+{{\varepsilon }_{\ell }}{{(\rho {{c}_{p}})}_{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }}\cdot \nabla T=\nabla \cdot \left( \left\langle k \right\rangle \nabla T \right)+{{\rho }_{s}}{{h}_{s\ell }}\frac{\partial {{\varepsilon }_{s}}}{\partial t}\qquad\qquad(14)

where eqs. (5) and (8) are used to determine the mass production rate of the solid and liquid phases per unit volume. The volume-averaged heat capacity and thermal conductivity are defined by

\left\langle \rho {{c}_{p}} \right\rangle ={{\varepsilon }_{\ell }}{{(\rho {{c}_{p}})}_{\ell }}+{{\varepsilon }_{s}}{{(\rho {{c}_{p}})}_{s}}\qquad\qquad(15)
\left\langle k \right\rangle ={{\varepsilon }_{\ell }}{{k}_{\ell }}+{{\varepsilon }_{s}}{{k}_{s}}\qquad\qquad(16)

Assuming there is no source/sink of the ith component in the solid and liquid phases, the conservations of species of the ith component for the solid and liquid phases are

\frac{\partial }{\partial t}\left( {{\varepsilon }_{s}}{{\rho }_{s}}{{\left\langle {{\omega }_{s,i}} \right\rangle }^{s}} \right)=-\nabla \cdot \left\langle {{{\mathbf{{\dot{m}}''}}}_{s,i}} \right\rangle +{{{\dot{m}}'''}_{\ell s,i}}\qquad\qquad(17)
\frac{\partial }{\partial t}\left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\omega }_{\ell ,i}} \right\rangle }^{\ell }} \right)+\nabla \cdot \left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\omega }_{\ell ,i}} \right\rangle }^{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }} \right)=-\nabla \cdot \left\langle {{{\mathbf{{\dot{m}}''}}}_{\ell ,i}} \right\rangle +{{{\dot{m}}'''}_{s\ell ,i}}\qquad\qquad(18)

Adding eqs. (17) and (18) together and considering {{{\dot{m}}'''}_{\ell s,i}}+{{{\dot{m}}'''}_{s\ell ,i}}=0, the conservation of species for the ith component becomes

\frac{\partial }{\partial t}\left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\omega }_{\ell ,i}} \right\rangle }^{\ell }}+{{\varepsilon }_{s}}{{\rho }_{s}}{{\left\langle {{\omega }_{s,i}} \right\rangle }^{s}} \right)+\nabla \cdot \left( {{\varepsilon }_{\ell }}{{\rho }_{\ell }}{{\left\langle {{\omega }_{\ell ,i}} \right\rangle }^{\ell }}{{\left\langle {{\mathbf{V}}_{\ell }} \right\rangle }^{\ell }} \right)=-\nabla \cdot \left( \left\langle {{{\mathbf{{\dot{m}}''}}}_{k,i}} \right\rangle +\left\langle {{{\mathbf{{\dot{m}}''}}}_{s,i}} \right\rangle  \right)\qquad\qquad(19)

Rizzo et al. (2003) presented a numerical simulation of solidification and mushy zone deformation of an organic alloy (SNC-acetone) in a two-dimensional rectangular cavity cooled from its left sidewall and all other three walls were kept adiabatic. Succinonitrile (SCN) is a transparent organic substance in the liquid state, and acetone was added to obtain an alloy behavior. The liquidus temperature for the pure SCN is 58.1 °C. The initial concentration of acetone was 11.3 wt% and the initial temperature of the melt alloy was at 42.35 °C. The cold wall temperature was set to 5 °C in numerical simulation. The rectangular cavity was deformed from its initial size of 50×170 mm to a final size of 40×210 mm.


Beckermann, C., and Viskanta, R., 1993, “Mathematical Modeling of Transport Phenomena during Alloy Solidification,” Applied Mechanics Reviews, Vol. 46, pp. 1-27.

Faghri, A., and Zhang, Y., 2006, Transport Phenomena in Multiphase Systems, Elsevier, Burlington, MA

Rizzo, E. M. da S., Santos, R. G. and Beckermann, C., 2003, “Modeling Solidification and Mushy Zone Deformation of Alloys,” Journal of Brazil Society of the Mechanical Science and Engineering, Vol. 25, pp. 180-184.

Further Reading

External Links