LINEBURG


<< . .

 13
( 23)



. . >>


including a computer graphics simulation of water, and ¬nally a chapter on
fully two-phase incompressible ¬‚ow. We wrap up this incompressible ¬‚ow
material with a chapter on incompressible ¬‚ames, including applications in
Part IV computer graphics, and a chapter on techniques for coupling a compressible
¬‚uid to an incompressible ¬‚uid. Finally, we turn our attention to heat ¬‚ow

Computational Physics with a chapter discussing the heat equation and a chapter discussing level
set techniques for solving Stefan problems.




While the ¬eld of computational ¬‚uid dynamics is quite broad, a large por-
tion of it is dedicated to computations of compressible ¬‚ow, incompressible
¬‚ow, and heat ¬‚ow. In fact, these three classes of problems can be thought
of as the basic model problems for hyperbolic, elliptic, and parabolic par-
tial di¬erential equations. Volumes have been ¬lled with both broad and
detailed work dedicated to these important ¬‚ow ¬elds; see, for example,
[86], [87], [7], [8], and the references therein. One might assume that there
is little left to add to the understanding of these problems. In fact, many
papers are now concerned with smaller details, e.g., the number of grid
points in a shock wave or contact discontinuity. Other papers are devoted
to rarely occurring pathologies, e.g., slow-moving shock waves and shock
overheating at solid wall boundaries.
One area where signi¬cant new ideas are still needed is multicomponent
¬‚ow, for example, multicomponent compressible ¬‚ow where di¬erent ¬‚u-
ids have di¬erent equations of state, multicomponent incompressible ¬‚ow
where di¬erent ¬‚uids have di¬erent densities and viscosities with surface
tension forces at the interface, and Stefan-type problems where the individ-
ual materials have di¬erent thermal conductivities. These problems have
interfaces separating the di¬erent materials, and special numerical tech-
niques are required to treat the interface. The most commonly used are
front tracking, volume of ¬‚uid, and level set methods. In the next three
chapters we discuss the use of level set methods for the canonical equations
of computational ¬‚uid dynamics.
The ¬rst chapter discusses basic one-phase compressible ¬‚ow, and the
subsequent chapter shows how the level set method can be used for two-
phase compressible ¬‚ow problems where the equations of state di¬er across
the interface. Then, the use of level set techniques for de¬‚agration and det-
onation discontinuities is discussed in the third chapter. The fourth chapter
introduces techniques for coupling an Eulerian grid to a Lagrangian grid.
This is useful, for example, in compressible solid/¬‚uid structure problems.
After this, we turn our attention to incompressible ¬‚ow with a chapter
focused on the basic one-phase equations including a computer graphics
simulation of smoke, a chapter on level set techniques for free surface ¬‚ows
150 14. Hyperbolic Conservation Laws and Compressible Flow

system of conservation laws. These also form the basis for their numerical
modeling.
A conserved quantity, such as mass, can be transported by convective

14 or di¬usive ¬‚uxes. The distinction is that di¬usive ¬‚uxes are driven by
gradients, while convective ¬‚uxes persist even in the absence of gradi-
ents. For most ¬‚ows where compressibility is important, e.g., ¬‚ows with
Hyperbolic Conservation Laws and shock waves, one needs to model only the convective transport and can ig-
nore di¬usion (mass di¬usion, viscosity and thermal conductivity) as well
Compressible Flow as the source terms (such as chemical reactions, atomic excitations, and
ionization processes). Moreover, convective transport requires specialized
numerical treatment, while di¬usive and reactive e¬ects can be treated with
standard numerical methods, such as simple central di¬erencing, that are
independent of those for the convective terms. Sti¬ reactions, however, can
present numerical di¬culties; see, for example, Colella et al. [50]. Conserva-
tion laws with only convective ¬‚uxes are known as hyperbolic conservation
laws. A vast array of physical phenomena are modeled by such systems,
e.g., explosives and high-speed aircraft.
The important physical phenomena exhibited by hyperbolic conserva-
tion laws are bulk convection, waves, contact discontinuities, shocks, and
rarefactions. We brie¬‚y describe the physical features and mathematical
model equations for each e¬ect, and most importantly note the implica-
We begin this chapter by addressing general systems of hyperbolic conserva- tions they have on the design of numerical methods. For more details on
tion laws including numerical techniques for computing accurate solutions numerical methods for conservation laws, see, e.g., LeVeque [105] and Toro
to them. Then we discuss the equations for one-phase compressible ¬‚ow as [164].
an example of a system of hyperbolic conservation laws.

14.1.1 Bulk Convection and Waves
14.1 Hyperbolic Conservation Laws Bulk convection is simply the bulk movement of matter, carrying it from
one spot to another, like water streaming from a hose. Waves are small-
A continuum physical system is described by the laws of conservation of amplitude smooth disturbances that are transmitted through the system
mass, momentum, and energy. That is, for each conserved quantity, the rate without any bulk transport like ripples on a water surface or sound waves
of change of the total amount in some region is given by its ¬‚ux (convective through air. Whereas convective transport occurs at the gross velocity
or di¬usive) through the region boundary, plus whatever internal sources of the material, waves propagate at the “speed of sound” in the system
exist. The integral form of this conservation law is (relative to the bulk convective motion of the system). Waves interact by
superposition, so that they can either cancel out (interfere) or enhance each
d
f (u) · dA =
u dV + s(u) dV (14.1) other.
dt „¦ ‚„¦ „¦
The simplest model equation that describes bulk convective transport is
where u is the density of the conserved quantity, f (u) is the ¬‚ux, and s(u) the linear convection equation
is the source rate. By taking „¦ to be an in¬nitesimal volume and applying
ut + v · u = 0, (14.3)
the divergence theorem, we get the di¬erential form of the conservation
law, where v is a constant equal to the convection velocity. The solution to this
is simply that u translates at the constant speed v. This same equation can
‚u
+ · f (u) = s(u), (14.2) also be taken as a simple model of wave motion if u is a sine wave and v is
‚t
interpreted as the speed of sound. The linear convection equation is also an
which is the basis for the numerical modeling of continuum systems. A important model for understanding smooth transport in any conservation
physical system can be described by a system of such equations, i.e., a
14.1. Hyperbolic Conservation Laws 151 152 14. Hyperbolic Conservation Laws and Compressible Flow

14.1.3 Shock Waves
law. As long as f is smooth and u has no jumps in it, the general scalar
conservation law
A shock is a spatial jump in material properties, like pressure and tem-
perature, that develops spontaneously from smooth distributions and
· f (u) = 0
ut + (14.4)
then persists. The shock jump is self-forming and also self-maintaining.
This is unlike a contact discontinuity, which must be put in the sys-
can be rewritten as
tem initially and will not resharpen itself if it is smeared out by some
ut + f (u) · u = 0, (14.5) other process. Shocks develop through a feedback mechanism in which
strong impulses move faster than weak ones, and thus tend to steepen
where f (u) acts as a convective velocity. That is, locally in smooth parts themselves up into a “step” pro¬le as they travel through the system. Fa-
of the ¬‚ow, a conservation law behaves like bulk convection with velocity miliar examples are the “sonic boom” of a jet aircraft and the “bang”
f (u). This is called the local characteristic velocity of the ¬‚ow. from a gun. These sounds are our perceptions of a sudden jump in air
Bulk convection and waves are important because they imply that sig- pressure.
nals propagate in de¬nite directions at de¬nite speeds. This is in contrast The simplest model equation that describes shock formation is the one-
to di¬usion, which propagates signals in all directions at arbitrarily large dimensional Burgers™ equation
speeds depending on the severity of the driving gradients. Thus we antici-
u2
pate that suitable numerical methods for hyperbolic systems will also have
ut + = 0, (14.6)
directional biases in space, which leads to the idea of upwind di¬erencing 2 x
and a de¬nite relation between the space and time steps (discrete propa-
which looks like the convection equation (14.3) with a nonconstant convec-
gation speed), which will roughly be that the discrete propagation speed
tive speed of u, i.e., ut + uux = 0. Thus larger u values move faster, and
∆x/∆t must be at least as large as the physical propagation speeds (char-
they will overtake smaller values. This ultimately results in the develop-
acteristic speeds) in the problem. The general form of this relation is called
ment of, for example, a right-going shock if the initial data for u constitute
the Courant-Friedrichs-Lewy (CFL) restriction.
any positive, decreasing function.
Wave motion and bulk convection do not create any new sharp features
Shocks move at a speed that is not simply related to the bulk ¬‚ow speed
in the ¬‚ow. The other remaining phenomena are all special because they
or characteristic speed, and they are not immediately evident from exam-
involve discontinuous jumps in the transported quantities. Because smooth
ining the ¬‚ux, in contrast to contacts. Shock speed is controlled by the
features can be accurately represented by a polynomial interpolation, we
di¬erence between in¬‚ux and out¬‚ux of conserved quantity into the region.
expect to be able to develop numerical methods of extremely high accuracy
Speci¬cally, suppose a conserved quantity u with conservation law
for the wave and convective e¬ects. Conversely, since jump functions are
poorly represented by polynomials, we expect little accuracy and perhaps
ut + f (u)x = 0 (14.7)
great di¬culty in numerically approximating the discontinuous phenomena.
has a step function pro¬le with constant values extending both to the left,
uL , and to the right, uR , with a single shock jump transition in between
14.1.2 Contact Discontinuities moving with speed s. Then the integral form of the conservation law (14.1),
applied to any interval containing the shock, gives the relation
A contact discontinuity is a persistent discontinuous jump in mass density
moving by bulk convection through a system. Since there is negligible mass s(uR ’ uL ) = f (uR ) ’ f (uL ), (14.8)
di¬usion, such a jump persists. These jumps usually appear at the point
which is just another statement that the rate at which u appears, s(uR ’
of contact of di¬erent materials; for example, a contact discontinuity sep-
arates oil from water. Contacts move at the local bulk convection speed, uL ), in the interval of interest is given by the di¬erence in ¬‚uxes across
or more generally the characteristic speed, and can be modeled by using the interval. Thus we see that the proper speed of the shock is directly
step-function initial data in the bulk convection equation (14.3). Since con- determined by conservation of u via the ¬‚ux f . This has an important
tacts are simply a bulk convection e¬ect, they retain any perturbations implication for numerical method design; namely, a numerical method will
they receive. Thus we expect contacts to be especially sensitive to numeri- “capture” the correct shock speeds only if it has “conservation form,” i.e.,
cal methods; i.e., any spurious alteration of the contact will tend to persist if the rate of change of u at some node is the di¬erence of ¬‚uxes that are
and accumulate. accurate approximations of the real ¬‚ux f .
14.1. Hyperbolic Conservation Laws 153 154 14. Hyperbolic Conservation Laws and Compressible Flow

14.2 Discrete Conservation Form
The self-sharpening feature of shocks has two implications for numerical
methods. First, it means that even if the initial data are smooth, steep
gradients and jumps will form spontaneously. Thus, our numerical method To ensure that shocks and other steep gradients are captured by the scheme,
must be prepared to deal with shocks even if none are present in the ini- i.e., that they move at the right speed even if they are unresolved, we must
tial data. Second, there is a bene¬cial e¬ect from self-sharpening, because write the equation in a discrete conservation form. That is, a form in which
modest numerical errors introduced near a shock (smearing or small os- the rate of change of conserved quantities is equal to a di¬erence of ¬‚uxes.
cillations) will tend to be eliminated, and will not accumulate. The shock This form guarantees that we discretely conserve the total amount of the
is naturally driven toward its proper shape. Because of this, computing states u (e.g., mass, momentum, and energy) present, in analogy with the
strong shocks is mostly a matter of having a conservative scheme in order integral form given by equation (14.1). More important, this can be shown
to get their speed correct. to imply that steep gradients or jumps in the discrete pro¬les propagate at
the physically correct speeds; see, for example, LeVeque [105].
Usually, conservation form is derived for control volume methods, that
is methods that evolve cell average values in time rather than nodal values.
14.1.4 Rarefaction Waves In this approach, a grid node xi is assumed to be the center of a grid cell
A rarefaction is a discontinuous jump or steep gradient in properties that (xi’1/2 , xi+1/2 ), and we integrate the conservation law (14.7) across this
dissipates as a smooth expansion. A common example is the jump in air control volume to obtain
pressure from outside to inside a balloon, which dissipates as soon as the
ut + f (ui+1/2 ) ’ f (ui’1/2 ) = 0,
¯ (14.9)
balloon is burst and the high-pressure gas inside is allowed to expand.
Such an expansion also occurs when the piston in an engine is rapidly
where u is the integral of u over the cell, and ui±1/2 are the (unknown)
¯
pulled outward from the cylinder. The expansion (density drop) associated
values of u at the cell walls. This has the desired conservation form in that
with a rarefaction propagates outward at the sound speed of the system,
the rate of change of the cell average is a di¬erence of ¬‚uxes. The di¬culty
relative to the underlying bulk convection speed. A rarefaction can be
with this formulation is that it requires transforming between cell averages
modeled by Burgers™ equation (14.6) with initial data that start out as
of u (which are directly evolved in time by the scheme) and cell wall values
a steep increasing step. This step will broaden and smooth out during the
of u (which must be reconstructed) to evaluate the needed ¬‚uxes. While this
evolution.
is manageable in one spatial dimension, in higher-dimensional problems the
A rarefaction tends to smooth out local features, which is generally ben-
series of transformations necessary to convert the cell averages to cell wall
e¬cial for numerical modeling. It tends to diminish numerical errors over
quantities becomes increasingly complicated. The distinction between cell
time and make the solution easier to represent by polynomials, which form
average and midpoint values can be ignored for schemes whose accuracy
the basis for our numerical representation. However, a rarefaction often
is no higher than second order (e.g., TVD schemes), since the cell average
connects to a smooth (e.g., constant) solution region and this results in
and the midpoint value di¬er by only O(∆x2 ).
a “corner,” which is notoriously di¬cult to capture accurately. The main
Shu and Osher [150, 151] proposed a fully conservative ¬nite di¬erence
numerical problem posed by rarefactions is that of initiating the expansion.
scheme on uniform grids that directly evolves nodal values (as opposed
If the initial data is are perfect, symmetrical step, such as u(x) = sign(x),
to the cell average values) forward in time. They de¬ned a numerical ¬‚ux
it may be “stuck” in this form, since the steady-state Burgers™ equation is
function F by the property that the real ¬‚ux divergence is a ¬nite di¬erence
satis¬ed identically (i.e., the ¬‚ux u2 /2 is constant everywhere, and similarly
of numerical ¬‚uxes
in any numerical discretization). However, local analysis can identify this
F (x + ∆x/2) ’ F (x ’ ∆x/2)
stuck expansion, because the characteristic speed u on either side points
f (u)x = (14.10)
away from the jump, suggesting its potential to expand. In order to get ∆x
the initial data unstuck, a small amount of smoothing must be applied to
at every point x. We call F the numerical ¬‚ux, since we require it in our
introduce some intermediate-state values that have a nonconstant ¬‚ux to
numerical scheme, and also to distinguish it from the closely related “phys-
drive expansion. In numerical methods this smoothing applied at a jump
ical ¬‚ux” f (u). It is not obvious that the numerical ¬‚ux function exists, but
where the e¬ective local velocity indicates expansion should occur is called
from relationship (14.10) one can solve for its Taylor expansion to obtain
an “entropy ¬x,” since it allows the system to evolve from the arti¬cial
low entropy initial state to the proper increased entropy state of a free (∆x)2 7(∆x)4
F = f (u) ’ f (u)xxxx ’ · · · ,
f (u)xx + (14.11)
expansion.
24 5760
14.3. ENO for Conservation Laws 155 156 14. Hyperbolic Conservation Laws and Compressible Flow

which shows that the physical and numerical ¬‚ux functions are the same only nodal values of the conserved variables. Their method is faster and eas-
to second-order accuracy in ∆x. Thus, a ¬nite di¬erence discretization can ier to implement than the cell-averaged formulation. In addition, the ¬nite
be based on di¬erence ENO method extends to higher dimensions in a dimension-by-
dimension fashion, so that the one-dimensional method applies unchanged
F (x + ∆x/2) ’ F (x ’ ∆x/2)
ut + =0 (14.12) to higher-dimensional problems. We emphasize that this is not dimensional
∆x
splitting in time, which has accuracy limitations unlike the dimension-by-
to evolve point values of u forward in time using numerical ¬‚ux functions F dimension approach. Shu and Osher also use the method of lines for time
at the cell walls. integration, decoupling the time and space discretizations.
We consider the treatment of a one-dimensional contact discontinuity
to illustrate how the method works. Assuming that the time evolution
takes place exactly, each time step ∆t should rigidly translate the spatial
14.3 ENO for Conservation Laws
pro¬le by the amount v∆t as governed by equation (14.3). Spatially, the
contact is initially represented by a discrete step function, i.e., nodal values
14.3.1 Motivation that are constant at one value uL on nodes x1 , . . . , xJ , and constant at a
Essentially nonoscillatory (ENO) methods were developed to address the di¬erent value uR on all remaining nodes xJ+1 , . . . , xN . To update the
special di¬culties that arise in the numerical solution of systems of non- value ui in time at a given node xi , we ¬rst reconstruct the graph of a
linear conservation laws. Numerical methods for these problems must be function u(x) near xi by interpolating nearby nodal u values, shift that
able to handle steep gradients, e.g., shocks and contact discontinuities, u(x) graph spatially by v∆t (the exact time evolution), and then reevaluate
that may develop spontaneously and then persist in these ¬‚ows. Classical it at the node xi to obtain the updated ui . We require our local interpolant
numerical schemes had a tendency either to produce large spurious oscil- be smooth at the point xi , since in actual practice we are going to use it to
lations near steep gradients or to greatly smear out both these gradients evaluate the derivative term ux there. The simplest symmetric approach to
and the ¬ne details of the ¬‚ow. The primary goal of the ENO e¬ort was smooth interpolation near a node xi is to run a parabola through the nodal
to develop a general-purpose numerical method for systems of conservation data at xi’1 , xi , and xi+1 . This interpolation is an accurate reconstruction
laws that has high accuracy (e.g., third order) in smooth regions and cap- of u(x) in smooth regions, where it works well. However, near the jump
tures the motion of unresolved steep gradients without creating spurious between xJ and xJ+1 the parabola will signi¬cantly overshoot the nodal u
data by an amount comparable to the jump uL ’ uR , and this overshoot
oscillations and without the use of problem-dependent ¬xes or tunable pa-
rameters. The philosophy underlying the ENO methods is simple: When will show up in the nodal values once the shift is performed. Successive
reconstructing a pro¬le for use in a convective ¬‚ux term, one should not time steps will further enhance these spurious oscillations. This approach
use high-order polynomial interpolation across a steep gradient in the data. corresponds to standard central di¬erencing.
Such an interpolant would be highly oscillatory and ultimately corrupt the To avoid the oscillations from parabolic interpolation, we could instead
computed solution. ENO methods use an adaptive polynomial interpola- use a smooth linear interpolation near xi , noting that there are two linear
tion constructed to avoid steep gradients in the data. The polynomial is also interpolants to choose from, namely the line through the data at nodes xi
biased to use data from the direction of information propagation (upwind) and xi’1 , and the line through the data at xi and xi+1 . The direction of in-
for physical consistency and stability. formation propagation determines which should be used. If the convection
The original ENO schemes developed by Harten et al. [81] were based speed v is positive, the data are moving from left to right, and we use xi
on the conservative control volume discretization of the equations, which and xi’1 . This linear interpolation based on upwind nodes will not intro-
yields discrete evolution equations for grid cell averages of the conserved duce any new extrema in u as long as the shift v∆t is less than the width
of the interval ∆x = xi ’ xi’1 , which is exactly the CFL restriction on the
quantities, e.g. mass, momentum, and energy. This formulation has the
disadvantage of requiring complicated transfers between cell averages and time step. The main problem with the linear upwind biased interpolant is
cell center nodal values in the algorithm. In particular, the transfer pro- that it has low accuracy smearing out the jump over more and more nodes.
cess becomes progressively more complicated in one, two, and three spatial If we naively go to higher accuracy by using a higher-order upwind biased
dimensions. The formulation also results in space and time discretizations interpolant, such as running a parabola through xi , xi’1 , and xi’2 to ad-
that are coupled in a way that becomes complicated for higher-order accu- vance ui , we run into the spurious oscillation problem again. In particular,
rate versions. To eliminate these complications, Shu and Osher [150, 151] at nodes xJ+1 and xJ+2 , this upwind parabola will interpolate across the
developed a conservative ¬nite di¬erence form of the ENO method that uses jump and thus have large overshoots. By forcing the parabola to cross a
14.3. ENO for Conservation Laws 157 158 14. Hyperbolic Conservation Laws and Compressible Flow

jump, it no longer re¬‚ects the data on the interval that will be arriving at where the second equality sign comes from
xJ+1 (or xJ+2 ) during the next time step. i
xi+1/2 xj+1/2
The motivation for ENO is that we must use a higher-degree polynomial H(xi+1/2 ) = h(y) dy = h(y) dy (14.18)
interpolant to achieve more accuracy, and it must involve the immediate up- x’1/2 xj’1/2
j=0
wind node to properly represent the propagation of data. But we must also i
avoid polluting this upwind data with spurious oscillations that come from = x f (u(xj )), (14.19)
interpolating across jumps. Thus, the remaining interpolation nodes (after j=0
the ¬rst upwind point) are chosen based on smoothness considerations. In
and the higher divided di¬erences are
particular, this approach will, if at all possible, not run an interpolant across
f (u(xi+1 )) ’ f (u(xi )) 11
a jump in the data. However, very small interpolation overshoots do occur 2
Di+1/2 H = = Di+1/2 f, (14.20)
near extrema in the nodal data, as they must, since any smooth function 2x 2
will slightly overshoot its values as sampled at discrete points near extrema. 12
3
Di H = Di f, (14.21)
This is the sense in which the method is only essentially nonoscillatory. 3
continuing in that manner.
According to the rules of polynomial interpolation, we can take any path
14.3.2 Constructing the Numerical Flux Function
along the divided di¬erence table to construct H, although not all paths
We de¬ne the numerical ¬‚ux function through the relation give good results. ENO reconstruction consists of two important features.
1
First, choose Di H in the upwind direction. Second, choose higher-order
Fi+1/2 ’ Fi’1/2
f (ui )x = (14.13) divided di¬erences by taking the smaller in absolute value of the two pos-
∆x
sible choices. Once we construct H(x), we evaluate H (xi+1/2 ) to get the
as in equation (14.12). To obtain a convenient algorithm for computing this
numerical ¬‚ux Fi+1/2 .
numerical ¬‚ux function, we de¬ne h(x) implicitly through the equation
x+ x/2
1
14.3.3 ENO-Roe Discretization (Third-Order Accurate)
f (u(x)) = h(y) dy, (14.14)
x x’ x/2
For a speci¬c cell wall, located at xi0 +1/2 , we ¬nd the associated numerical
and note that taking a derivative on both sides of this equation yields
¬‚ux function Fi0 +1/2 as follows. First, we de¬ne a characteristic speed
x/2) ’ h(x ’
h(x + x/2) »i0 +1/2 = f (ui0 +1/2 ), where ui0 +1/2 = (ui0 + ui0 +1 )/2 is de¬ned using a
f (u(x))x = , (14.15)
x standard linear average. Then, if »i0 +1/2 > 0, set k = i0 . Otherwise, set
k = i0 + 1. De¬ne
which shows that h is identical to the numerical ¬‚ux function at the cell
walls. That is, Fi±1/2 = h(xi±1/2 ) for all i. We calculate h by ¬nding its 1
Q1 (x) = (Dk H)(x ’ xi0 +1/2 ). (14.22)
primitive 2 2 2
If |Dk’1/2 H| ¤ |Dk+1/2 H|, then c = Dk’1/2 H and k = k ’ 1. Otherwise,
x
2
c = Dk+1/2 H and k = k. De¬ne
H(x) = h(y) dy (14.16)
x’1/2
Q2 (x) = c(x ’ xk’1/2 )(x ’ xk+1/2 ). (14.23)
using polynomial interpolation, and then take a derivative to get h. Note
3 3 3 3
If |Dk H| ¤ |Dk +1 H|, then c = Dk H. Otherwise, c = Dk +1 H. De¬ne
that we do not need the zeroth-order divided di¬erences of H that vanish
with the derivative. Q3 (x) = c (x ’ xk ’ xk ’ xk
’1/2 )(x +1/2 )(x +3/2 ). (14.24)
0
The zeroth order divided di¬erences, Di+1/2 and all higher-order even
Then
divided di¬erences of H exist at the cell walls and have the subscript i±1/2.
1
The ¬rst order divided di¬erences Di and all higher-order odd divided Fi0 +1/2 = H (xi0 +1/2 ) = Q1 (xi0 +1/2 )+Q2 (xi0 +1/2 )+Q3 (xi0 +1/2 ), (14.25)
di¬erences of H exist at the grid points and will have the subscript i. The
which simpli¬es to
¬rst-order divided di¬erences of H are
1
3(i0 ’ k )2 ’ 1 ( x)2 .
Fi0 +1/2 = Dk H + c (2(i0 ’ k) + 1) x+c
H(xi+1/2 ) ’ H(xi’1/2 )
1
Di H = = f (u(xi )), (14.17) (14.26)
x
14.3. ENO for Conservation Laws 159 160 14. Hyperbolic Conservation Laws and Compressible Flow

14.4 Multiple Spatial Dimensions
14.3.4 ENO-LLF Discretization (and the Entropy Fix)
The ENO-Roe discretization can admit entropy-violating expansion shocks
In multiple spatial dimensions, the ENO discretization is applied inde-
near sonic points. That is, at a place where a characteristic velocity changes
pendently using a dimension-by-dimension discretization. For example,
sign (a sonic point) it is possible to have a stationary expansion shock so-
consider a two-dimensional conservation law
lution with a discontinuous jump in value. If this jump were smoothed out
even slightly, it would break up into an expansion fan (i.e., rarefaction) ut + f (u)x + g(u)y = 0 (14.30)
and dissipate, which is the desired physical solution. For a speci¬c cell wall
on a rectangular 2D grid. Here, we sweep through the grid from bottom
xi0 +1/2 , if there are no nearby sonic points, then we use the ENO-Roe dis-
to top performing ENO on 1D horizontal rows of grid points to evaluate
cretization. Otherwise, we add high-order dissipation to our calculation of
the f (u)x term. The g(u)y term is evaluated in a similar manner, sweeping
Fi0 +1/2 to break up any entropy-violating expansion shocks. We call this
through the grid from left to right performing ENO on 1D vertical rows
entropy-¬xed version of the ENO-Roe discretization ENO-Roe ¬x or just
of grid points. Once we have a numerical approximation to each of the

<< . .

 13
( 23)



. . >>

Copyright Design by: Sunlight webdesign