LINEBURG

 << Пред. стр. страница 20(всего 23)ОГЛАВЛЕНИЕ След. стр. >>
velocity п¬Ѓeld is continuous across the interface. Using this, the jump
(2Вµux )x + (Вµ(uy + vx ))y + (Вµ(uz + wx ))z
px
condition
= , (21.15)
ПЃ ПЃ
[p] в€’ 2 [Вµ] u В· N , v В· N , w В· N В· N = ПѓОє (21.10)
(Вµ(uy + vx ))x + (2Вµvy )y + (Вµ(vz + wy ))z
py
= , (21.16)
can be written for the pressure. Notice that this reduces to [p] = ПѓОє when ПЃ ПЃ
the viscosity is continuous across the interface. Further derivations lead to
(Вµ(uz + wx ))x + (Вµ(vz + wy ))y + (2Вµwz )z
pz
пЈ« пЈ¶ = , (21.17)
[Вµux ] [Вµuy ] [Вµuz ] ПЃ ПЃ
пЈ­ [Вµvx ] [Вµvy ] [Вµvz ] пЈё
in expanded form.
[Вµwx ] [Вµwy ] [Вµwz ]
The two-phase incompressible п¬‚ow equations are discretized in the
пЈ«пЈ« пЈ¶пЈ« пЈ¶T пЈ« пЈ¶
0 0
u same manner as the equations for one-phase incompressible п¬‚ow. First,
пЈ¬
= [Вµ] пЈ­пЈ­ v пЈё пЈ­ T1 пЈё пЈ­ T1 пЈё an intermediate velocity п¬Ѓeld V is computed using equation (18.19).
w Then equation (18.21) is solved to п¬Ѓnd the pressure, which is used in
T2 T2
пЈ« пЈ¶ equation (18.20) to make the velocity п¬Ѓeld divergence free. Due to the dis-
u
continuous nature of several quantities across the interface, special care
+ N N пЈ­ v пЈё NT N
T
(21.11)
is needed in discretizing the viscous terms in equation (18.19) and in
w
discretizing the pressure and density in equations (18.21) and (18.20).
пЈ¶
пЈ« пЈ¶T пЈ« пЈ¶пЈ« пЈ¶T
0 0 u
пЈ­ T1 пЈё пЈ­ T1 пЈё пЈ­ v пЈё N T N пЈ· ,
в€’ пЈё
21.3 Viscous Terms
w
T2 T2

which is useful for discretizing the viscous terms, especially since the right- In the Оґ-function approach the density and the viscosity are numerically
hand side of this equation involves only derivatives that are continuous smeared out, so that they are continuous across the interface. The contin-
across the interface. Notice that all the quantities on the left-hand side of uous viscosity simpliп¬Ѓes the jump conditions, allowing the viscous terms
this equation become continuous across the interface when [Вµ] = 0 forces to be discretized just as they were for one-phase incompressible п¬‚ow. The
the right-hand side to be identically zero. See  for details. only diп¬Ђerence is that ПЃ and Вµ are deп¬Ѓned by equations (21.1) and (21.2)
Since the velocity is continuous across the interface, the material respectively. Averaging of П† is preferred to averaging of other quantities.
derivative, or Lagrangian, acceleration is continuous as well, For example, the viscosity at xi+ 1 ,j,k is deп¬Ѓned as
2

П†i,j,k + П†i+1,j,k
Du Dv Dw
Вµi+ 1 ,j,k = Вµ , (21.18)
= = = 0. (21.12)
2
Dt Dt Dt 2

as opposed to
Since the Navier-Stokes equations (18.2), (18.3), and (18.4) are valid on
Вµ (П†i,j,k ) + Вµ (П†i+1,j,k )
both sides of the interface, these equations do not jump across the interface; Вµi+ 1 ,j,k = , (21.19)
2
i.e., 2

as was done in equation (18.16).
T
В· П„)
p( If one prefers to keep the density and viscosity sharp across the interface,
Vt + V В· в€’ в€’ g = 0.
V+ (21.13)
the sign of the level set function can be used to determine Вµ as Вµв€’ or Вµ+
ПЃ ПЃ
21.3. Viscous Terms 233 234 21. Two-Phase Incompressible Flow

and to determine ПЃ as ПЃв€’ or ПЃ+ . Consider the case where Вµ and ПЃ are both Suppose that П†L в‰¤ 0 and П†M > 0, so that the interface lies between the
independently spatially constant on either side of the interface, allowing associated grid points. Then
the simpliп¬Ѓcation of the viscous terms to |П†L |
Оё= (21.26)
|П†L | + |П†M |
Вµ (uxx + uyy + uzz )
, (21.20)
ПЃ can be used to estimate the interface location by splitting this cell into
two pieces of size Оё x on the left and (1 в€’ Оё) x on the right. At the
interface, we denote the continuous velocity by uI and calculate the jump
Вµ (vxx + vyy + vzz )
, (21.21) as JI = ОёJM + (1 в€’ Оё)JL . Then we discretize the jump condition [Вµux ] = JI
ПЃ
as
uM в€’ uI uI в€’ uL
Вµ (wxx + wyy + wzz ) в€’ Вµв€’
Вµ+ = JI , (21.27)
, (21.22) (1 в€’ Оё) x Оёx
ПЃ
solving for uI to obtain
В· V = 0. Since the velocities are continuous, their п¬Ѓrst
with the aid of
Вµ+ uM Оё + Вµв€’ uL (1 в€’ Оё) в€’ JI Оё(1 в€’ Оё) x
derivatives can be computed directly. However, the jump conditions in uI = , (21.28)
Вµ+ Оё + Вµв€’ (1 в€’ Оё)
equation (21.11) are needed to compute the second derivatives.
The right-hand side of equation (21.11) needs to be computed in order to so that we can write
evaluate the jumps across the interface. First, the continuous velocity п¬Ѓeld is uM в€’ uI uM в€’ u L Вµ JI Оё
Л†
(Вµux )L = Вµ+ =Вµ
Л† + , (21.29)
averaged from the MAC grid to the cell centers. Then central diп¬Ђerencing
Вµв€’
(1 в€’ Оё) x x
is used to compute the п¬Ѓrst derivatives at each cell center. These п¬Ѓrst
where
derivatives are multiplied by the appropriate components of the normal
Вµ+ Вµв€’
and tangent vectors to obtain a numerical estimate for the right-hand side
Вµ=
Л† (21.30)
of equation (21.11), which we denote by the matrix J. Since J is a spatially Вµ+ Оё + Вµв€’ (1 в€’ Оё)
continuous function, spatial averages can be used to deп¬Ѓne J elsewhere. For
deп¬Ѓnes an eп¬Ђective Вµ. Similarly, if П†L > 0 and П†M в‰¤ 0, then
example, Ji+ 1 ,j,k = (Ji,j,k + Ji+1,j,k ) /2.
2
uM в€’ uI uM в€’ uL Вµ JI Оё
Л†
Once J has been computed, the second derivatives are computed using (Вµux )L = Вµв€’ в€’
=Вµ
Л† , (21.31)
Вµ+
(1 в€’ Оё) x x
techniques similar to those developed in . For example, consider the
discretization of Вµuxx at xi+ 1 ,j,k using uM = ui+ 1 ,j,k and its neighbors where
2 2
uL = uiв€’ 1 ,j,k and uR = ui+ 3 ,j,k . We need averaged values of П† and J 11 Вµв€’ Вµ +
2 2
Вµ= в€’
Л† (21.32)
(the appropriate scalar entry of J) at the same three spatial locations as Вµ Оё + Вµ+ (1 в€’ Оё)
the u terms. If П†L , П†M , and П†R are all greater than zero, we deп¬Ѓne
deп¬Ѓnes an eп¬Ђective Вµ.
uM в€’ uL In similar fashion, if П†R > 0 and П†M в‰¤ 0, then
(Вµux )L = Вµ+ (21.23)
x |П†R |
Оё= (21.33)
|П†R | + |П†M |
and
is used to estimate the interface location with (1 в€’ Оё) x on the left and
uR в€’ uM
(Вµux )R = Вµ+ , (21.24) Оё x on the right. Then JI = ОёJM + (1 в€’ Оё)JR is used to discretize the
x
jump condition as
arriving at
uR в€’ uI uI в€’ uM
в€’ Вµв€’
Вµ+ = JI , (21.34)
(Вµux )R в€’ (Вµux )L (1 в€’ Оё) x
Оёx
(Вµuxx )i+ 1 ,j,k = (21.25)
x
2
resulting in
Вµв€’ uM Оё + Вµ+ uR (1 в€’ Оё) в€’ JI Оё(1 в€’ Оё) x
in the standard fashion. A similar discretization holds when all three П†
uI = (21.35)
values are less than or equal to zero. Вµв€’ Оё + Вµ+ (1 в€’ Оё)
21.4. Poisson Equation 235 236 21. Two-Phase Incompressible Flow

and and solve for pI as
uI в€’ uM uR в€’ u M Вµ JI Оё
Л† ОІ + pk+1 Оё + ОІ в€’ pk (1 в€’ Оё) в€’ bОё(1 в€’ Оё) x
(Вµux )R = Вµв€’ в€’
=Вµ
Л† , (21.36) pI = , (21.43)
Вµ+
(1 в€’ Оё) x x ОІ + Оё + ОІ в€’ (1 в€’ Оё)
where so that approximations to the derivatives on the left and right sides of the
в€’+ interface can be written as
ВµВµ
Вµ=
Л† (21.37)
Л†
Вµв€’ Оё + Вµ+ (1 в€’ Оё) pI в€’ p k Л† pk+1 в€’ pk в€’ ОІb(1 в€’ Оё)
ОІв€’ =ОІ (21.44)
ОІ+
Оёx x
deп¬Ѓnes an eп¬Ђective Вµ. If П†R в‰¤ 0 and П†M > 0, then
and
uI в€’ uM uR в€’ u M Вµ JI Оё
Л†
(Вµux )R = Вµ+ =Вµ
Л† + , (21.38) Л†
Вµв€’
(1 в€’ Оё) x x pk+1 в€’ pI Л† pk+1 в€’ pk ОІbОё
ОІ+ =ОІ + , (21.45)
ОІв€’
(1 в€’ Оё) x x
where
Вµ+ Вµв€’ where
Вµ= +
Л† (21.39)
ОІ+ОІв€’
Вµ Оё + Вµв€’ (1 в€’ Оё) Л†
ОІ= (21.46)
ОІ + Оё + ОІ в€’ (1 в€’ Оё)
deп¬Ѓnes an eп¬Ђective Вµ.
deп¬Ѓnes an eп¬Ђective ОІ.
The new equations for the unknowns pk and pk+1 are then
21.4 Poisson Equation (pk+1 в€’a)в€’pk pk в€’pkв€’1
Л† b(1в€’Оё)
в€’ в€’ ОІkв€’ 1
ОІ ОІ+
x x
2
= fk (21.47)
x
Consider solving
and
В· (ОІ p) = f (21.40)
pk+2 в€’pk+1 uk+1 в€’(pk +a)
Л† bОё
в€’ОІ
ОІk+ 3 +
for the pressure p with speciп¬Ѓed jump conditions of [p] = a and [ОІpn ] = b ОІв€’
x x
2
= fk+1 , (21.48)
across the interface. In the context of equation (18.21), ОІ = 1/ПЃ and x
f = ( В· V )/ t. When the Оґ-function method is used, this equation is where we add the a term to correct for the fact the pressure is discontinuous
straightforward to solve, since a and b are set to 0, and equation (21.1) across the interface as well. Note that this correction was not necessary
is used to deп¬Ѓne a smeared-out continuous value for ПЃ. In this case, one in treating the viscosity, since the velocity п¬Ѓeld is continuous across the
obtains a numerically smeared-out pressure proп¬Ѓle that does not include interface. These new equations for pk and pk+1 can be rewritten in standard
surface-tension forces. As discussed above, the forcing function deп¬Ѓned form as
in equation (21.3) can be added to the momentum equations to recover
pk+1 в€’pk pk в€’pkв€’1
Л† в€’ ОІkв€’ 1
ОІ Л† Л†
surface-tension eп¬Ђects. On the other hand, if one wants to model the ОІb(1 в€’ Оё)
ОІa
x x
2
= fk + + (21.49)
surface-tension eп¬Ђects directly, the jump conditions for the pressure cannot ( x)2 ОІ+ x
x
be ignored.
and
First, consider the one-dimensional problem where a standard second-
pk+2 в€’pk+1 pk+1 в€’pk
Л†
order accurate discretization of в€’ОІ
ОІk+ 3 Л† Л†
ОІa ОІbОё
x x
2
= fk+1 в€’ +в€’ (21.50)
pi+1 в€’pi pi в€’piв€’1
в€’ ОІiв€’ 1
ОІi+ 1 ( x)2
x ОІ x
x x
2 2
= fi (21.41)
emphasizing that this discretization yields the standard symmetric linear
x
Л†
system with ОІk+1/2 = ОІ.
can be written for each unknown pi . Suppose that the interface is lo-
More generally, at each grid point i, we write a linear equation of the
cated between xk and xk+1 . As in the treatment of the viscosity term,
form
we discretize the jump condition [ОІpx ] = b, obtaining, for example,
pi+1 в€’pi pi в€’piв€’1
в€’ ОІiв€’ 1
ОІi+ 1
pk+1 в€’ pI pI в€’ p k x x
в€’ ОІв€’
ОІ+ = fi + F L + F R
2 2
= b, (21.42) (21.51)
(1 в€’ Оё) x Оёx x
21.4. Poisson Equation 237 238 21. Two-Phase Incompressible Flow

and assemble the system of linear equations into matrix form. Each ОІk+1/2 are already balanced out on the right-hand side by the appropriate jumps
is evaluated based on the side of the interface that xk and xk+1 lie on, and included in the V term.
Л†
a special ОІ is used when xk and xk+1 lie on opposite sides of the interface. The resulting linear system of equations can still be solved using a PCG
Then if the left arm of the stencil (the line segment connecting xi and gradient with an incomplete Choleski preconditioner, just as in the case
xiв€’1 ) crosses the interface, a nonzero F L is deп¬Ѓned with correction terms of one-phase incompressible п¬‚ow. However, one needs to use caution when
for [p] = a and [ОІpn ] = b. Likewise, if the right arm of the stencil (the plugging the resulting pressure into equation (18.20), since the pressure
line segment connecting xi and xi+1 ) crosses the interface, a nonzero F R is discontinuous across the interface. The pressure derivatives in equa-
is deп¬Ѓned with correction terms for [p] = a and [ОІpn ] = b. tion (18.20) should be computed in exactly the same fashion as they were
The multidimensional approach is treated in a dimension-by-dimension computed in solving equation (18.21); i.e., the correction terms are still
fashion. While the [p] = a jump condition is trivial to apply, some assump- needed.
tions are made in order to obtain a dimension-by-dimension approach for
[ОІpn ] = b. For example, in two spatial dimensions we assume that
[ОІpx ] = [ОІpn ]n1 (21.52)
and
[ОІpy ] = [ОІpn ]n2 , (21.53)
where n1 and n2 are the components of the local unit normal. Although
these equations are not generally true, adding n1 times the п¬Ѓrst equation to
n2 times the second equation leads to the correct [ОІpn ] = b jump condition
and numerically demonstrated convergence to the correct solution. The
errors in this approach can be characterized by adding t1 times the п¬Ѓrst
equation to t2 times the second equation to obtain [ОІpt ] = 0, implying that
the tangential derivative is incorrectly smeared out. The two-dimensional
application of the method consists in writing a linear equation of the form
pi+1,j в€’pi,j pi,j в€’piв€’1,j
в€’ ОІiв€’1/2,j
ОІi+1/2,j x x
x
pi,j+1 в€’pi,j pi,j в€’pi,jв€’1
в€’ ОІi,jв€’1/2
ОІi,j+1/2 y y
+
y
x y
= fi,j + F + F (21.54)
at each grid point, where F x = F L + F R and F y = F B + F T are obtained
by considering each spatial dimension independently using either [ОІpx ] =
[ОІpn ]n1 = bn1 or [ОІpy ] = [ОІpn ]n2 = bn2 , respectively.
Before using the above-described numerical method to solve equa-
tion (18.21), the jump condition given in equation (21.10) needs to be
computed. This can be done with standard central diп¬Ђerencing of the av-
eraged cell-centered velocities, analogous to the way that J was computed
in discretizing the viscous terms. Note that we can set [px /ПЃ] = [py /ПЃ] =
[pz /ПЃ] = 0 in spite of the nonzero jumps in these quantities. Since the full
equations are continuous across the interface, one can take the divergence
of the full equations to derive equation (18.21) without the need for correc-
tion terms. The jumps in the derivatives of the pressure in equation (18.21)
240 22. Low-Speed Flames

to the discontinuity of the density, viscosity, and pressure). Delta-function
formulations smear out this velocity jump, forcing a continuous velocity
п¬Ѓeld across the interface. This can be problematic, since this numerical
smearing adds a compressible character to the п¬‚ow п¬Ѓeld near the inter-
22 face. The divergence-free condition is not exactly satisп¬Ѓed in the separate
subdomains. In addition, diп¬ѓculties arise in computing the interface veloc-
ity, which depends on the local velocity of the unreacted material. Near the
Low-Speed Flames interface, the velocity of the unreacted material contains large O(1) numer-
ical errors where it has been nonphysically forced to be continuous with the
velocity of the reacted material. Partial solutions to these problems were
proposed by Helenbrook et al. , where the authors were able to remove
the numerical smearing of the normal velocity, obtaining a sharp interface
proп¬Ѓle. This method works well as long as п¬‚ame fronts remain well sepa-
rated with moderate curvature; see Helenbrook and Law . This method
cannot treat merging п¬‚ame fronts or individual fronts with relatively high
curvature. These drawbacks were recently overcome by Nguyen et al. 
who extended the work of Kang et al.  to treat this problem.

22.2 Governing Equations
22.1 Reacting Interfaces We ignore viscous eп¬Ђects and consider the equations for inviscid incom-
pressible п¬‚ow
In Chapter 21 the interface moved with the local п¬‚uid velocity only,
p
and individual п¬‚uid particles did not cross the interface. In this chap- Vt + V В· V+ =0 (22.1)
ПЃ
ter we consider interfaces across which a chemical reaction is converting
one incompressible п¬‚uid into another. The interface moves with the lo- independently for each п¬‚uid. The interface velocity is W = DN , where D is
cal velocity of the unreacted п¬‚uid plus a reaction term that accounts for the normal component of the interface velocity deп¬Ѓned by D = (VN )u + S.
the conversion of one п¬‚uid into the other as material moves across the The вЂњuвЂќ subscript indicates that the normal velocity is calculated using
interface. Consider an interface separating liquid and gas regions where
the velocity of the unreacted material only. This is important to note,
the liquid is actively vaporizing into the gaseous state. Juric and Tryg- since VN is discontinuous across the interface. The п¬‚ame speed is deп¬Ѓned
gvason  developed a front-tracking approach to this problem using a as S = So + ПѓОє where, So and Пѓ are constants and Оє is the local curvature
Оґ-function formulation to treat the interface boundary conditions. Son and of the interface.
Dir  and Welch and Wilson  developed level-set-based and volume- Conservation of mass and momentum imply the standard Rankine-
of-п¬‚uid-based (respectively) approaches to this same problem also using a Hugoniot jump conditions across the interface
Оґ-function formulation.
Another example of reacting interfaces occurs in premixed п¬‚ames. As- [ПЃ(VN в€’ D)] = 0 (22.2)
suming an inп¬Ѓnitely thin п¬‚ame front allows us to treat the п¬‚ame front as a ПЃ(VN в€’ D)2 + p = 0 (22.3)
discontinuity separating two incompressible п¬‚ows. The unreacted material
as well as continuity of the tangential velocities, [VT1 ] = [VT2 ] = 0, as long
undergoes reaction as it crosses the interface, producing a lower-density
as S = 0. Note that S = 0 only in the case of a contact discontinuity (not
(higher-volume) reacted material. Qian et al.  devised a front-tracking
a п¬‚ame). Denoting the mass п¬‚ux in the moving reference frame (speed D)
approach to this problem using a Оґ-function formulation.
by
Typically, the density is discontinuous across the interface. Thus, mate-
rial must instantaneously expand as it crosses the interface, implying that
M = ПЃr ((VN )r в€’ D) = ПЃu ((VN )u в€’ D) (22.4)
the normal velocity is discontinuous across the interface as well (in addition
22.3. Treating the Jump Conditions 241 242 22. Low-Speed Flames

allows us to rewrite equation (22.2) as [M ] = 0. Here the вЂњrвЂќ subscript to obtain
denotes a reacted material quantity. Substitution of D = (VN )u + S into 1 1
uG ur в€’ M в€’
= n1 , (22.12)
equation (22.4) yields u
ПЃr ПЃu
1 1
M = в€’ПЃu S, (22.5) G
vr в€’ M в€’
vu = n2 , (22.13)
ПЃr ПЃu
which is a rather simple quantity for computations.
 << Пред. стр. страница 20(всего 23)ОГЛАВЛЕНИЕ След. стр. >>