# Application of Computational Methods

(Difference between revisions)
 Revision as of 05:22, 7 April 2010 (view source)← Older edit Revision as of 05:23, 7 April 2010 (view source)Newer edit → Line 20: Line 20: $\rho \frac{DV}{Dt}=-\nabla p+\nabla \cdot \left( \mu \nabla V \right)+\rho \mathbf{g}+\mathbf{F}$ $\rho \frac{DV}{Dt}=-\nabla p+\nabla \cdot \left( \mu \nabla V \right)+\rho \mathbf{g}+\mathbf{F}$ - |{{EquationRef|(1)}} + |{{EquationRef|(2)}} |} |} Since the flow field is assumed to be independent of the temperature field, a steady-state energy equation is solved in the fluid region once the flow field has been resolved. Since the flow field is assumed to be independent of the temperature field, a steady-state energy equation is solved in the fluid region once the flow field has been resolved. Line 29: Line 29: $\nabla \cdot \left( V\rho h \right)=\nabla \cdot \left( k\nabla T \right)$ $\nabla \cdot \left( V\rho h \right)=\nabla \cdot \left( k\nabla T \right)$ - |{{EquationRef|(1)}} + |{{EquationRef|(3)}} |} |} In the solid region, only conduction is solved. In the solid region, only conduction is solved. Line 38: Line 38: $\nabla ^{2}T=0$ $\nabla ^{2}T=0$ - |{{EquationRef|(1)}} + |{{EquationRef|(4)}} |} |} The free surface is tracked by the Volume of Fluid (VOF) method, where the volume fraction of the fluid, ε, is tracked through each computational cell.  The VOF equation is: The free surface is tracked by the Volume of Fluid (VOF) method, where the volume fraction of the fluid, ε, is tracked through each computational cell.  The VOF equation is: Line 47: Line 47: $\frac{\partial \varepsilon }{\partial t}+V\cdot \nabla \varepsilon =0$ $\frac{\partial \varepsilon }{\partial t}+V\cdot \nabla \varepsilon =0$ - |{{EquationRef|(1)}} + |{{EquationRef|(5)}} |} |} The interface between fluids is represented by a piecewise linear approach, similar to the work of Youngs (1982), in order to greatly limit numerical diffusion of the interface.  Surface tension effects are modeled in the numerical simulations.  The surface tension forces are represented by the F term in eq. (4.327).  A continuum surface force method proposed by Brackbill et al. (1992) is used to model surface tension. The interface between fluids is represented by a piecewise linear approach, similar to the work of Youngs (1982), in order to greatly limit numerical diffusion of the interface.  Surface tension effects are modeled in the numerical simulations.  The surface tension forces are represented by the F term in eq. (4.327).  A continuum surface force method proposed by Brackbill et al. (1992) is used to model surface tension. Line 56: Line 56: $F=\sigma \frac{\rho K\nabla \varepsilon }{0.5\left( \rho _{\ell }+\rho _{v} \right)}$ $F=\sigma \frac{\rho K\nabla \varepsilon }{0.5\left( \rho _{\ell }+\rho _{v} \right)}$ - |{{EquationRef|(1)}} + |{{EquationRef|(6)}} |} |} The curvature, K, is defined as: The curvature, K, is defined as: Line 65: Line 65: $K=-\frac{\nabla ^{2}\varepsilon }{\left| \nabla \varepsilon \right|}$ $K=-\frac{\nabla ^{2}\varepsilon }{\left| \nabla \varepsilon \right|}$ - |{{EquationRef|(1)}} + |{{EquationRef|(7)}} |} |} The fluid properties are calculated by the volume weighted average: The fluid properties are calculated by the volume weighted average: Line 74: Line 74: $\varphi =\varepsilon \varphi _{\ell }+\left( 1-\varepsilon \right)\varphi _{v}$ $\varphi =\varepsilon \varphi _{\ell }+\left( 1-\varepsilon \right)\varphi _{v}$ - |{{EquationRef|(1)}} + |{{EquationRef|(8)}} |} |} The general fluid property, The general fluid property, Line 85: Line 85: $h=\frac{\left[ \varepsilon \rho _{\ell }c_{\ell }+\left( 1-\varepsilon \right)\rho _{v}c_{v} \right]\left( T-T_{ref} \right)}{\rho }$ $h=\frac{\left[ \varepsilon \rho _{\ell }c_{\ell }+\left( 1-\varepsilon \right)\rho _{v}c_{v} \right]\left( T-T_{ref} \right)}{\rho }$ - |{{EquationRef|(1)}} + |{{EquationRef|(9)}} |} |} The boundary conditions are as follows: The boundary conditions are as follows: Line 97: Line 97: $\left( 0 - |{{EquationRef|(1)}} + |{{EquationRef|(10)}} |} |} Line 105: Line 105: [itex]\left( -d_{d} - |{{EquationRef|(1)}} + |{{EquationRef|(11)}} |} |} At the collar At the collar Line 116: Line 116: [itex]u=v=0,\text{ }w=\omega r,\text{ }\frac{\partial T}{\partial x}=0$ $u=v=0,\text{ }w=\omega r,\text{ }\frac{\partial T}{\partial x}=0$ - |{{EquationRef|(1)}} + |{{EquationRef|(12)}} |} |} At the disk surface At the disk surface Line 128: Line 128: $u=v=0,\text{ }w=\omega r,\text{ }\left. k\frac{\partial T}{\partial x} \right|_{\ell }=\left. k\frac{\partial T}{\partial x} \right|_{S}$ $u=v=0,\text{ }w=\omega r,\text{ }\left. k\frac{\partial T}{\partial x} \right|_{\ell }=\left. k\frac{\partial T}{\partial x} \right|_{S}$ - |{{EquationRef|(1)}} + |{{EquationRef|(13)}} |} |} At the disk bottom At the disk bottom Line 139: Line 139: $\left( r_{in} - |{{EquationRef|(1)}} + |{{EquationRef|(14)}} |} |} Line 147: Line 147: [itex]\left( r_{hin} - |{{EquationRef|(1)}} + |{{EquationRef|(15)}} |} |} Line 155: Line 155: [itex]\left( r_{hout} - |{{EquationRef|(1)}} + |{{EquationRef|(16)}} |} |} At the outer boundaries At the outer boundaries Line 165: Line 165: [itex]p=p_{ref}$ $p=p_{ref}$ - |{{EquationRef|(1)}} + |{{EquationRef|(17)}} |} |} {| class="wikitable" border="0" {| class="wikitable" border="0" Line 172: Line 172: if $\left( V\cdot n>0 \right)$ then $\frac{\partial T}{\partial n}=0$, else $T=T_{ref}$ if $\left( V\cdot n>0 \right)$ then $\frac{\partial T}{\partial n}=0$, else $T=T_{ref}$ - |{{EquationRef|(1)}} + |{{EquationRef|(18)}} |} |} At the liquid-vapor interface At the liquid-vapor interface Line 181: Line 181: $p_{\ell }-p_{v}-\left. \mu _{\ell }\frac{\partial V\cdot n}{\partial n} \right|_{\ell }+\left. \mu _{v}\frac{\partial V\cdot n}{\partial n} \right|_{v}=\sigma K$ $p_{\ell }-p_{v}-\left. \mu _{\ell }\frac{\partial V\cdot n}{\partial n} \right|_{\ell }+\left. \mu _{v}\frac{\partial V\cdot n}{\partial n} \right|_{v}=\sigma K$ - |{{EquationRef|(1)}} + |{{EquationRef|(19)}} |} |} Line 189: Line 189: $\left. \mu _{\ell }\frac{\partial V\cdot t}{\partial n} \right|_{\ell }=\left. \mu _{v}\frac{\partial V\cdot t}{\partial n} \right|_{v}$ $\left. \mu _{\ell }\frac{\partial V\cdot t}{\partial n} \right|_{\ell }=\left. \mu _{v}\frac{\partial V\cdot t}{\partial n} \right|_{v}$ - |{{EquationRef|(1)}} + |{{EquationRef|(20)}} |} |} heating: heating: Line 198: Line 198: $\left. k_{\ell }\frac{\partial T}{\partial n} \right|_{\ell }=\left. k_{v}\frac{\partial T}{\partial n} \right|_{v}$ $\left. k_{\ell }\frac{\partial T}{\partial n} \right|_{\ell }=\left. k_{v}\frac{\partial T}{\partial n} \right|_{v}$ - |{{EquationRef|(1)}} + |{{EquationRef|(21)}} |} |} evaporation: evaporation: Line 207: Line 207: $T=T_{sat}$ $T=T_{sat}$ - |{{EquationRef|(1)}} + |{{EquationRef|(22)}} |} |} Rice et al. (2005) solved this problem numerically with the result presented below. Rice et al. (2005) solved this problem numerically with the result presented below.

## Revision as of 05:23, 7 April 2010

The computational methods presented in the previous sections provide the solution techniques for solving the governing equations and boundary conditions. However, it is equally important to mathematically set up the physical problem before applying the numerical techniques. In many practical problems for external flow one needs to deal with conjugate effects, free surface, and non-conventional geometry including multiple domain regions with non-uniform boundary conditions. In order to show these effects, and the approach to deal with them, an application is presented below that includes all these non-conventional effects.

Figure 1: Controlled liquid impinging jet onto a rotating disk (Rice et. al. 2005).
Figure 2: Computational domain for a controlled liquid impinging jet onto a rotating disk (Rice et. al. 2005).

The physical problem that will be presented in this section is heat and mass transfer from a controlled liquid impinging jet over a rotating disk, Figure 1, including the conjugate wall effect for both heating with and without evaporation from the free surface. First we need to develop the conservative governing equations and boundary conditions for fluid mechanics, heat transfer, and volume of fluid (VOF) for a controlled liquid impinging jet on a rotating disk with a free surface. Figure 2 shows a computational domain for the three regions: disk, liquid, and vapor, with an extended domain to include the edge effects. The basic assumptions are two dimensional, constant properties, incompressible, and laminar flow. For evaporation, the effect of mass transfer is neglected. It is also assumed that there is constant heat flux along the heated section.

Flow enters the disk between two circular plates, one being the collar, and the other being the disk. The space between the collar and the disk is δin. After the flow leaves the entrance region between the collar and the disk at rhin, the flow turns from an internal fully developed flow to a free surface flow. The heater provides a constant heat flux at the outside wall, from the bottom of the aluminum disk. The heat is conducted through the disk to the fluid. The basic assumptions of the problem at hand are that the flow field is incompressible, and the fluid is considered to be in the laminar flow regime while all of the fluid properties are constant. When evaporation is being modeled, no mass transfer is modeled and the effects are only considered in the energy equation, because the evaporation rate of the liquid is much less than the inlet mass flow rate of the liquid. The Navier-Stokes equations are solved to compute the fluid flow. The continuity and momentum equations are as follows:

 $\nabla \cdot V=0$ (1)
 $\rho \frac{DV}{Dt}=-\nabla p+\nabla \cdot \left( \mu \nabla V \right)+\rho \mathbf{g}+\mathbf{F}$ (2)

Since the flow field is assumed to be independent of the temperature field, a steady-state energy equation is solved in the fluid region once the flow field has been resolved.

 $\nabla \cdot \left( V\rho h \right)=\nabla \cdot \left( k\nabla T \right)$ (3)

In the solid region, only conduction is solved.

 $\nabla ^{2}T=0$ (4)

The free surface is tracked by the Volume of Fluid (VOF) method, where the volume fraction of the fluid, ε, is tracked through each computational cell. The VOF equation is:

 $\frac{\partial \varepsilon }{\partial t}+V\cdot \nabla \varepsilon =0$ (5)

The interface between fluids is represented by a piecewise linear approach, similar to the work of Youngs (1982), in order to greatly limit numerical diffusion of the interface. Surface tension effects are modeled in the numerical simulations. The surface tension forces are represented by the F term in eq. (4.327). A continuum surface force method proposed by Brackbill et al. (1992) is used to model surface tension.

 $F=\sigma \frac{\rho K\nabla \varepsilon }{0.5\left( \rho _{\ell }+\rho _{v} \right)}$ (6)

The curvature, K, is defined as:

 $K=-\frac{\nabla ^{2}\varepsilon }{\left| \nabla \varepsilon \right|}$ (7)

The fluid properties are calculated by the volume weighted average:

 $\varphi =\varepsilon \varphi _{\ell }+\left( 1-\varepsilon \right)\varphi _{v}$ (8)

The general fluid property, $\varphi$ , represents either density, viscosity or thermal conductivity. The enthalpy is calculated using a mass-weighted average instead of a volume-weighted average.

 $h=\frac{\left[ \varepsilon \rho _{\ell }c_{\ell }+\left( 1-\varepsilon \right)\rho _{v}c_{v} \right]\left( T-T_{ref} \right)}{\rho }$ (9)

The boundary conditions are as follows: At the inlet $\left( r=r_{in} \right)$

 $\left( 0 (10)
 $\left( -d_{d} (11)

At the collar $\left( x=\delta _{in},\text{ }r_{in}

 $u=v=0,\text{ }w=\omega r,\text{ }\frac{\partial T}{\partial x}=0$ (12)

At the disk surface $\left( x=0,\text{ }r_{in}

 $u=v=0,\text{ }w=\omega r,\text{ }\left. k\frac{\partial T}{\partial x} \right|_{\ell }=\left. k\frac{\partial T}{\partial x} \right|_{S}$ (13)

At the disk bottom $\left( x=-d_{d} \right)$

 $\left( r_{in} (14)
 $\left( r_{hin} (15)
 $\left( r_{hout} (16)

At the outer boundaries $\left( r=r_{hin},\delta _{in},$\left( r_{hin},$\left( r=r_{out},\delta _{in}

 p = pref (17)
 if $\left( V\cdot n>0 \right)$ then $\frac{\partial T}{\partial n}=0$, else T = Tref (18)

At the liquid-vapor interface

 $p_{\ell }-p_{v}-\left. \mu _{\ell }\frac{\partial V\cdot n}{\partial n} \right|_{\ell }+\left. \mu _{v}\frac{\partial V\cdot n}{\partial n} \right|_{v}=\sigma K$ (19)
 $\left. \mu _{\ell }\frac{\partial V\cdot t}{\partial n} \right|_{\ell }=\left. \mu _{v}\frac{\partial V\cdot t}{\partial n} \right|_{v}$ (20)

heating:

 $\left. k_{\ell }\frac{\partial T}{\partial n} \right|_{\ell }=\left. k_{v}\frac{\partial T}{\partial n} \right|_{v}$ (21)

evaporation:

 T = Tsat (22)

Rice et al. (2005) solved this problem numerically with the result presented below.

Figure 3: Film thickness vs. radial distance at an inlet temperature of 20ºC at δin = 0.254 mm and a liquid flow rate of 7 RPM (Rice et. al. 2005).
Figure 4: Local Nusselt Numbers vs. radial distance for δin = 0.204 mm, inlet temperature of 40ºC and liquid flow rate of 7 RPM (Rice et. al. 2005).

Air was used as the vapor for all of the simulations. The viscosity and thermal conductivity of air were lessened by an order of magnitude for selective cases, and were found not to change the results. The boundary condition at the interface is therefore essentially a zero shear stress for all cases, an insulated boundary for heating cases, and Tsat for evaporation. As was noted before, the flow field is considered to be independent of the thermal field; therefore the flow field is first solved as a time-dependent solution, until a steady state solution has been reached. Once a steady-state solution has been reached, the energy equation is solved for both the fluid and solid regions, as a steady-state solution. The criterion used to determine a steady-state solution was a mass flow rate of liquid leaving the computational domain within 0.05% of the mass flow rate of liquid entering the domain for greater than 0.05 seconds. The flow field was solved using the following procedure:
1. Solve momentum equations
2. Solve continuity (pressure correction) equation, and update pressure and face mass flow rate
3. Repeat steps 1 and 2 until converged, for each time step
A co-located finite volume computational scheme as described in Section 4.8, where both the flow-field variables and the pressure are stored in the cell centers, was used to solve the governing equations. The pressure was discretized in a manner similar to a staggered-grid scheme, while the pressure and velocity were coupled using the SIMPLE algorithm described in the previous section. All of the convective terms in the governing equations were discretized using a second-order upwind scheme. The grid used for the numerical simulations consisted of rectangular-shaped cells, which were produced in four axial layers in the liquid regions, and one axial layer in the solid region. In the first layer of the fluid region, spanning from the disk surface to δin, there were 25 cell rows, with a growth rate of 1.01. The second layer, spanning from δin to 3δin, consisted of 30 cell rows with a growth rate of 1.02. The third layer, spanning from 3δin to 10δin, consisted of 25 cell rows, and had a growth rate of 1.05. In the final layer, spanning from 10δin to 100δin, there were 20 cell rows, with a growth rate of 1.2. In the solid region, there were 15 evenly spaced cell rows. In the radial direction, there were 15 evenly spaced cell columns between rin and rhin, and 75 evenly spaced cell columns between rhin and rout. The grid spaned such a great distance in the axial direction (100δin), because convergence is greatly increased with the added cells, so the computations ran more rapidly.

Typical results for local film thickness δ and Nussult number, $\text{Nu}={q}''r_{in}/[k_{\ell }(T_{w}-T_{in})]$ , as a function of radial distance are presented in Fig. 3 and 4 respectively. The heating boundary condition was used for all cases with an inlet temperature of 20ºC to 40ºC. The evaporation boundary condition was used for all cases with an outlet temperature of 100ºC.