Baixe CFD I: Métodos Numéricos para Equações Hiperbólicas - Grétar Tryggvason e outras Notas de estudo em PDF para Engenharia Mecânica, somente na Docsity!
Numerical Methods!
for Hyperbolic!
Equations—I!
Grétar Tryggvason!
Spring 2009!
http://users.wpi.edu/~gretar/me612.html!
FTCS and upwind!
Stability in terms of fluxes!
Generalized upwind!
Second order schemes for smooth flow!
Modified Equation!
Conservation!
Computational Fluid Dynamics I
2
f
∂ t
2
− c
2
2
f
∂ x
2
∂u
∂t
− c
2
∂v
∂x
∂v
∂t
∂u
∂x
The wave equation:!
Write as:!
In general:!
∂ u
∂ t
∂ v
∂ t
a
11
a
12
a
21
a
22
∂ u
∂ x
∂ v
∂ x
Most of the issues involved can be addressed by examining:!
∂ f
∂ t
+ U
∂ f
∂ x
Computational Fluid Dynamics I
FTCS and Upwind!
Computational Fluid Dynamics I
We will start by examining the linear advection equation:!
∂ f
∂ t
+ U
∂ f
∂ x
The characteristic for this equation are:!
dx
dt
= U ;
df
dt
Showing that the
initial conditions are
simply advected by a
constant velocity U!
t
f
f
x
Computational Fluid Dynamics I
A simple forward in time, centered in space
discretization yields!
! f
! t
+ U
! f
! x
f
j
n + 1
= f
j
n
Δ t
2 h
U ( f
j + 1
n
− f
j − 1
n
j-1 j j+
n
n+
This scheme is O(Δt, h
2
) accurate, but a stability
analysis shows that the error grows as!
ε
n + 1
ε
n
= 1 − i
U Δ t
2 h
sin kh
Since the amplification
factor has the form 1+ i ()
the absolute value of this
complex number is always
larger than unity and the
method is unconditionally
unstable for this case.!
i
U Δ t
2 h
sin kh
1
ε
n + 1
ε
n
A simple forward in time but “upwind” in space
discretization yields!
∂ f
∂ t
+ U
∂ f
∂ x
f
j
n + 1
= f
j
n
Δ t
h
U ( f
j
n
− f
j − 1
n
j-1 j
n
n+
This scheme
is O(Δt, h)
accurate.!
Another scheme for!
Flow direction!
Computational Fluid Dynamics I
To examine the stability we use the von Neumanʼs method:!
ε
j
n + 1
− ε
j
n
Δ t
U
h
( ε
j
n
− ε
j − 1
n
) = 0
ε
j
n
= ε
n
e
ikx j
ε
n + 1
− ε
n
Δ t
ε
n
h
( 1 − e
− ikh
) = 0
The evolution of the error is governed by:!
Write the error as:!
ε
n + 1
ε
n
U Δ t
h
( 1 − e
− ikh
G =
ε
n + 1
ε
n
= 1 − λ( 1 − e
− ikh
), λ =
U Δ t
h
Amplification factor!
G = 1! " + " e
! ikh
Or:!
Need to find when! G < 1
Computational Fluid Dynamics I
Stable!
Im( G )
1
1
kh
λ
1-λ
G
U Δ t
h
This restriction was first
derived by Courant, Fredrik,
and Levy in 1932, and is
usually called the Courant
condition, or the CFL
condition.!
G = 1 − λ + λ e
− ikh
U Δ t
h
Stability condition: λ<1!
Graphically:!
Computational Fluid Dynamics I
G = 1 − λ + λ e
− ikh
= 1 − λ + λ cos kh − i λ sin kh
G
2
= 1 − λ + λ cos kh
2
2
sin
2
kh
2
+ 2 ( 1 − λ) λ cos kh + λ
2
cos
2
kh + λ
2
sin
2
kh
2
+ 2 ( 1 − λ) λ cos kh + λ
2
= 1 − 2 λ + 2 λ
2
λ cos kh
G
2
= 1 − 2 λ ( 1 − λ) if cos kh = 0
G
2
= 1 if cos kh = 1
G
2
= 1 − λ 4 − 3 λ ( )
if cos kh = − 1
G
2
! 1 if "! 1
Finf the absolute value of the amplification factor!
1!
G
2
λ
1!
cos kh = 0
cos kh = − 1
Computational Fluid Dynamics I
The CFL condition implies that a signal has to travel
less than one grid spacing in one time step!
U Δ t ≤ h
Flow direction!
j-1 j
n
n+
Allowable characteristics!
Not allowable characteristics!
MOVIE FROM MATLAB!
% one-dimensional advection by first order upwind.!
n=80; nstep=100; dt=0.0125; length=2.0;!
h=length/(n-1);y=zeros(n,1);f=zeros(n,1);f(1)=1.0;!
for m=1:nstep,m!
hold off, plot(f); axis([1, n, -0.5, 1.5]);!
pause(0.01);!
y=f;!
for i=2:n-1,!
f(i)=y(i)-(dt/h)*(y(i)-y(i-1)); %upwind!
end;!
end;!
f
j + 1
Consider the following initial conditions:!
f
j − 1
f
j
F
j −1/ 2
=
U
2
( f
j − 1
n
j
n
) = 1.0 F
j +1/ 2
U
( f
j
n
j + 1
n
f
j
n + 1
= f
j
n
Δ t
h
( F
j +1/ 2
n
− F
j −1/ 2
n
So cell j will overflow immediately!!
By considering the fluxes, it is easy to see why the
centered difference approximation is always unstable.!
Computational Fluid Dynamics I
dt=0.25*h!
20 40 60 80
-0.
0
1
Upwind
Although the
upwind
method is
exceptionally
robust, its low
accuracy in
space and
time makes it
unsuitable for
most serious
computations!
Computational Fluid Dynamics I
Numerical Methods!
for Hyperbolic!
Equations—II!
http://users.wpi.edu/~gretar/me612.html!
Grétar Tryggvason!
Spring 2009!
Computational Fluid Dynamics I
FTCS and upwind!
Stability in terms of fluxes!
Generalized upwind!
Second order schemes for smooth flow!
Modified Equation!
Conservation!
Computational Fluid Dynamics I
f
j
n + 1
= f
j
n
Δ t
h
U ( f
j
n
− f
j − 1
n
j-1 j
n
n+
O(Δt, h)
accurate.!
For the linear advection equation:!
Flow direction!
U Δ t
h
U
h
First Order Schemes!
The Upwind Scheme!
∂ f
∂ t
+ U
∂ f
∂ x
Finite Volume point of view:!
x
j 1/
x
j +1/
!
f
j − 1
f
j
x!
f
f
j + 1
f
j
n + 1
= f
j
n
" t
h
( F
j + 1 / 2
n
! F
j! 1 / 2
n
) = f
j
n
" t
h
U ( f
j
n
! f
j! 1
n
F
j + 1 / 2 =
Uf
j
n
F
j − 1 / 2 =
Uf
j − 1
n
First Order Schemes!
Generalized Upwind Scheme (for both U > 0 and U < 0 )!
f
j
n + 1
= f
j
n
−
U Δ t
h
( f
j
n
− f
j − 1
n
), U > 0
f
j
n + 1
= f
j
n
−
U Δ t
h
( f
j + 1
n
− f
j
n
), U < 0
Define:!
U
=
1
2
U + U
, U
−
=
1
2
U − U
The two cases can be combined into a single expression:!
f j
n + 1
= f j
n
−
Δ t
h
U
( f j
n
− f j − 1
n
) + U
−
( f j + 1
n
− f j
n
)
[ ]
First Order Schemes!
Or, substituting!
U , U
f
j
n + 1
= f
j
n
− U
Δ t
2 h
( f
j + 1
n
− f
j − 1
n
) +
U Δ t
2 h
( f
j + 1
n
− 2 f
j
n
j − 1
n
)
central difference + numerical viscosity!
D
num
=
U h
2
Computational Fluid Dynamics I
Other First Order Schemes!
Computational Fluid Dynamics I
Implicit (Backward Euler) Method!
- Unconditionally stable!
- 1st order in time, 2nd order in space!
- Forms a tri-diagonal matrix (Thomas algorithm)!
1
1
1
1
1
−
n
j
n
j
n
j
n
j
f f
h
U
t
f f
n
j
n
j
n
j
n
j
f
t
f
h
U
f
t
f
h
U
−
1
1
1 1
1
j
n
j j
n
j j
n
j j
a f + df + bf = C
−
1
1
1 1
1
First Order Schemes!
Computational Fluid Dynamics I
Lax-Fredrichs method!
- stable for λ < 1!
- 1st order in time, 1nd order in space!
- Conditionally consistent!
f
j
n + 1
1
2
f
j + 1
n
j − 1
n
f
j
n
Δ t
U
2 h
f
j + 1
n
− f
j − 1
n
xx xxx
f
Uh
f
Uh
2
2
λ λ
λ
Error term:!
First Order Schemes!
λ =
U Δ t
h
h
λ
=
h
2
U Δ t
Computational Fluid Dynamics I
Second Order Schemes!
QUICK!
s=1!
f
1
!
s=2!
f
2
!
s=3!
f
3
!
s=4!
f
4
!
s=5/2!
1 − 2 : ( 2 − s ) f
1
2
2 − 3 : ( 3 − s ) f
2
3
1 − 3 :
1
2
( 3 − s ) ( 2 − s ) f 1
[ ] +
1
2
( s − 1 ) ( 3 − s ) f 2
[ ]
s = 5 / 2 : f
5 / 2
=
1
8
( 3 f
3
2
− f
1
)
Interpolate linearly between 1-2 and 2-3!
Then quadratically between 1-2!
The function can now be evaluated at the half point!
Third Order Schemes!
∂ f
i
∂ x
=
1
h
f
i + 1 / 2
− f
i − 1 / 2
QUICK!
Approximate the derivative:!
∂ f
i
∂ x
=
1
h
f
i + 1 / 2
− f
i − 1 / 2
where the value at the half points by
second order upwind interpolation!
! f
! x
"
$
i
%
1
8 h
3 f
i
i & 1
& f
i & 2
[ ] & 3 f
i + 1
i
& f
i & 1
{ [ ]}
Third Order Schemes!
Computational Fluid Dynamics I
The one-step Lax-Wendroff is not easily extended to
non-linear or multi-dimensional problems. The split
version is.!
In the Lax-Wendroff and the MacCormack methods
the spatial and the temporal discretization are not
independent.!
Other methods have been developed where the time
integration is independent of the spatial discretization,
such as the Beam-Warming and various Runge-Kutta
methods!
Second Order Schemes!
Computational Fluid Dynamics I
FTCS!
Unconditionally
Unstable!
Upwind!
Stable for!
Implicit!
Unconditionally!
Stable!
Lax-Friedrichs!
Conditionally
consistent!
Stable for!
t x
f Uf
0
2
1 1
1
=
−
Δ
−
h
f f
U
t
f f
n
j
n
j
n
j
n
j
( )
xx xxx
f
Uh
f
U
t
2
2 2
12
2 6
−Δ − + λ
0
1
1
=
−
Δ
− −
h
f f
U
t
f f
n
j
n
j
n
j
n
j
( )
( ) xxx
xx
f
Uh
f
Uh
2 3 1
6
1
2
2
2
− − +
−
λ λ
λ
λ≤ 1
( ) xx xxx
f
Uh
f
Uh 2
2
1
3
1
2
λ λ
λ
⎟+ −
⎠
⎞
⎜
⎝
⎛
−
( )
0
2
1
1
1
1
1
=
−
Δ
−
−
h
f f
U
t
f f
n
j
n
j
n
j
n
j
xx xxx
f Uh Utf
U t
!
"
$
%
&
( + '
' 2 32
2
3
1
6
1
2
( )
( )
0
2
/ 2
1 1
1 1
1
=
−
Δ
− +
h
f f
U
t
f f f
n
j
n
j
n
j
n
j
n
j
λ≤ 1
Summary!
Computational Fluid Dynamics I
Leap Frog!
Stable for!
Lax-Wendroff I!
Stable for!
Lax-Wendroff II!
Same as LW-I!
Stable for!
MacCormack!
Same as LW-I!
Stable for!
t x
f Uf
0
2 2
1 1
1 1
=
−
Δ
−
h
f f
U
t
f f
n
j
n
j
n
j
n
j
( ) xxx
f
Uh
1
6
2
2
λ −
( )
( ) xxxx
xxx
f
Uh
f
Uh
2
3
2
2
1
8
1
6
!!
!
" "
" "
λ≤ 1
λ≤ 1
λ≤ 1
( )
( )
0
2
2
2
2
22 1 1
1 1
1
=
− +
− Δ
−
Δ
−
h
f f f
U t
h
f f
U
t
f f
n
j
n
j
n
j
n
j
n
j
n
j
n
j
0
/ 2
( )/ 2 1 1
1 / 2
1 / 2
−
Δ
− +
h
f f
U
t
f f f
n
j
n
j
n
j
n
j
n
j
0
1 / 2
1 / 2
12 /
1 / 2
1
=
−
Δ
−
−
h
f f
U
t
f f
n
j
n
j
n
j
n
j λ≤ 1
( )
0
1
=
−
Δ
−
h
f f
U
t
f f
n
j
n
j
n
j
t
j
( ) ( )
0
/ 2 1
1
=
−
Δ
− + −
h
f f
U
t
f f f
t
j
t
j
t
j
n
j
n
j
Summary!
Computational Fluid Dynamics I
Beam-Warming!
Stable for!
QUICK!
Stable for!
ENO!
WENO!
t x
f Uf
( )
( 2 ) 0
2
2
3 4
1 2 2
2
1 2
1
− + =
Δ
−
− +
Δ
−
− −
− −
n
j
n
j
n
j
n
j
n
j
n
j
n
j
n
j
f f f
h
U t
h
f f f
U
t
f f
0 ≤ λ≤ 2
( )( )
( ) ( ) xxxx
xxx
f
Uh
f
Uh
λ λ
λ λ
− − −
− −
1 2
8
1 2
6
2
3
2
Summary!
f j
n + 1
− f j
n
Δ t
U
3 f j
n
n
− f j − 2
n
( )
− 3 f j + 1
n
n
− f j − 1
n
( )
8 h
λ≤ 1
A large number of (conditionally) stable and
accurate methods exists for hyperbolic equations
with smooth solutions!
And Many More!!
Modified Equation!
Looking at the structure of the error terms—by
deriving the Modified Equation—can often lead
to insight into how the approximate solution
behaves!
Computational Fluid Dynamics I
Derive modified equation for upwind difference method:!
1
1
−
n
j
n
j
n
j
n
j
f f
h
U
t
f f
Using Taylor expansion:!
3
3
2 3
2
2
1
t
t
t f
t
f
t
t
f
f f
n
j
n
j
−
3
3
2 3
2
2
1
h
x
h f
x
f
h
x
f
f f
n
j
n
j
Computational Fluid Dynamics I
Substituting!
n
j
n
j
f
t
t
t f
t
f
t
t
f
f
t
3
3
2 3
2
2
3
3
2 3
2
2
h
x
h f
x
f
h
x
f
f f
h
U
n
j
n
j
Therefore,!
tt xx ttt xxx
f
Uh
f
t
f
Uh
f
t
x
f
U
t
f
2 2
It helps the interpretation if all terms are written in! xx xxx
f , f
Computational Fluid Dynamics I
Taking further derivatives (first in time, then in space):!
tt xt ttt xxt tttt xxxt
f
Uh
f
t
f
Uh
f
t
f Uf
2 2
tx xx ttx xxx tttx xxxx
f
U h
f
U t
f
U h
f
U t
Uf U f
2 2 22
2
2
2
u O h
U
u
U
x
f O t
f U
f Uf t
xxt xxx
ttx
ttt
tt xx
Computational Fluid Dynamics I
Similarly, we get!
3
f Uf O th ttt xxx
2
f Uf O th
ttx xxx
f Uf O ( t , h )
xxt xxx
Final form of the modified equation:!
xx xxx
f
Uh
f
Uh
x
f
U
t
f
2
2
λ λ λ
3
, h
2
Δ t , h Δ t
2
, Δ t
3
h
U Δ t
λ =
In finite volume method, equations in conservative forms
are needed in order to satisfy conservation properties.!
As an example, consider a 1-D equation!
∂ f ( x , t )
∂ t
∂ F [ f ( x , t )]
∂ x
where F denotes a general advection/diffusion term, e.g.!
F =
f
2
; F = μ( x )
∂ f
∂ x
∂ F
∂ x
dx
L
∫
F
j + 1 / 2
− F
j − 1 / 2
Δ x
∑
Δ x
X
j- 1/
x
j+ 1/
!
f j + 1
f
j − 1
f
j
x!
In discretized form:!
= + F
j − 1 / 2
− F
j − 3 / 2
+ F
j + 1 / 2
− F
j − 1 / 2
+ F
j + 3 / 2
− F
j + 1 / 2
[ ]
= F
L
− F
0
If F = 0 at the end points of the domain, f is conserved.!
Integrating over the domain L,!
∂ f
∂ t
dx
L
∫
∂ F
∂ x
dx
L
∫
F ( L ) − F (0) = 0
⇒
d
dt
fdx
L
∫
= 0
Computational Fluid Dynamics I
∂ f
∂ t
∂ x
f
2
Examples of Conservative Form!
Discretize using upwind!
∂ x
f
2
2 h
f i
2
− f i − 1
2
( ) Conservative!
∂ f
∂ t
∂ f
∂ x
f
∂ f
∂ x
f i
h
f
i
− f
i − 1
( )
Non-conservative!
versus!
Computational Fluid Dynamics I
∂ F
∂ x
dx
L
∫
≈ F
j + 1 / 2
− F
j − 1 / 2
∑
+ f
i − 1
2
− f
i − 2
2
i
2
− f
i − 1
2
i + 1
2
− f
i
2
∂ F
∂ x
dx
L
∫
≈ F
j + 1 / 2
− F
∑ j − 1 / 2
+ f
i − 1
2
− f
i − 1
f
i
i
2
− f
i
f
i − 1
i + 1
2
− f
i + 1
f
i
∂ x
f
2
i
2 h
f
i
2
− f
i − 1
2
( )
f
∂ f
∂ x
i
f
i
h
f
i
− f
i − 1
( )
Computational Fluid Dynamics I
∂ f
∂ t
∂ x
μ( x )
∂ f
∂ x
Another example (variable diffusion coefficient)!
Discretize!
∂ x
μ( x )
∂ f
∂ x
Conservative!
μ
2
f
∂ x
2
∂ μ
∂ x
∂ f
∂ x
Non-conservative!
versus!
Computational Fluid Dynamics I
∂ f
∂ t
∂ f
∂ x
f
∂ f
∂ t
∂ f
∂ x
∂ t
f
2
∂ x
f
3
d
dt
f
2
dx ∫
∂ x
f
3
dx ∫
f
3
−∞
∞
f ∞ ( )
= f −∞ ( )
d
dt
f
n
n
dx ∫
Similarly, it can
be shown that!
Conservation of higher order quantities!
Consider:!
Integrate:!
Generally it is believed that conservative
schemes are superior to non-conservative ones
and that the more conserved quantities the
discrete approximation preserves, the better the
scheme. While mostly true, there are exceptions
as we will discuss later!!
Numerical Methods!
for Hyperbolic!
Equations—III!
http://users.wpi.edu/~gretar/me612.html!
Grétar Tryggvason!
Spring 2009!
Computational Fluid Dynamics I
Flows with shocks. Solutions to Burgers eqn.
Entropy condition!
Conservation and shock speed!
Advecting a shock with several schemes!
Godunovʼs Theorem!
Artificial Viscosity!
Multidimensional flows!
Computational Fluid Dynamics I
Discontinuous!
Solutions!
Shocks-I!
Computational Fluid Dynamics I
∂ f
∂ t
+ U
∂ f
∂ x
The analytic solution is obtained by characteristics!
dx
dt
= U ;
df
dt
Consider the linear Advection Equation!
t
f
f
x
Discontinuity of solution is allowed!
(Weak Solution)!
L R
R
L
f f
f x x
f x x
f x >
0
0
dx
dt
= U
Computational Fluid Dynamics I
L R
R
L
f f
f x x
f x x
f x >
0
0
x
f
f
t
f
dt
df
f
dt
dx
Inviscid Burgersʼ Equation!
t
x
Characteristics!
0
x
L
f
R
f
= ??
dt
dx
0
x
Shock Wave!
∂ f
∂ t
+ f
∂ f
∂ x
Example!
F =
f
2
C =
F
R
− F
L
f
R
− f
L
f
R
2
− f
L
2
f
R
− f
L
f
R
− f
L
f
R
+ f
L
f
R
− f
L
C =
f
R
+ f
L
Discontinuous!
Solutions!
Shocks-II!
Computational Fluid Dynamics I
Weak solutions to hyperbolic equations may not be unique.!
How can we find a physical solution out of many weak!
solutions?!
In fluid mechanics, the actual physics always includes!
dissipation, i.e. in the form of viscous Burgersʼ equation:!
2
2
x
f
x
f
f
t
f
ε
Therefore, what we are truly seeking is the solution to!
the viscous Burgersʼ equation in the limit of!
ε→ 0
Computational Fluid Dynamics I
Entropy Condition:!
A discontinuity propagating with speed C satisfies the!
entropy condition if!
For a conservation equation!
L R
F ′ f > C > F ′ f
R
R
L
L
f f
Ff F f
C
f f
Ff F f
x
F f
t
f
Version I:!
Version II:!
And some others…!
Entropy Condition!
t
x
Shock Wave!
for!
f L
≥ f ≥ f R
Computational Fluid Dynamics I
Means that the characteristics on both
side must “slope into” the shock!
Start with:!
L R
F ′ f > C > F ′ f
x
F f
t
f
The condition:!
Entropy Condition!
t
x
Shock Wave!
F ' =
∂ F
∂ f
∂ f
∂ t
+ F '
∂ F
∂ x
where!
F ' R
F '
L
Computational Fluid Dynamics I
Similarly, the shock speed is given by!
R
R
L
L
f f
Ff F f
C
f f
Ff F f
Thus!
Entropy Condition!
t
x
Shock Wave!
C =
F
R
− F
L
f
R
− f
L
f
L
≥ f ≥ f
R
for!
Means that the hypothetical shock speed for values
of f between the left and the right state must give
shock speeds that are larger on the left and smaller
on the right.!
In CFD, however, we do not want to solve the viscous!
Burgersʼ equation with extremely small because it!
becomes computationally expensive.!
Nevertheless, in many gas dynamics application, there!
is need to capture discontinuous solution behavior.!
And we want to reproduce the weak solution behavior!
numerically by solving inviscid Burgersʼ equation.!
ε
- Shock speed!
- Entropy condition!
- Wave shape (without smearing or wiggles)!
Examples:!
!Godunov method!
!Approximate Riemann solver!
!Artificial dissipation!
In practice, rather than applying a separate entropy!
condition to the numerical solution procedure, numerical!
methods are developed first and then it can be proved!
that the methods satisfy the entropy condition.!
Entropy Condition!
Computational Fluid Dynamics I
Conservative Schemes are guarantied to
give the correct Shock Speed since they
correspond to a direct application of the
conservation principles.!
Non-conservative schemes may or may not
do so.!
Computational Fluid Dynamics I
In capturing the correct solution behavior for discontinuous!
Initial data, conservative methods are essential (Lax).!
Conservation and shock speed!
Example: inviscid Burgers equation:!
0 or
2
f
t x
f
x
f
f
t
f
Conservative Method!
Computational Fluid Dynamics I
Consider upwind and forward Euler scheme:!
Example: for inviscid Burgers equation!
with discontinuous initial data!
t x
f ff
j − 1
j j + 1
f
j
n + 1
= f
j
n
Δ t
h
f
j
n
f
j
n
− f
j − 1
n
( )
Non-conservative form!
Never moves!!
f
t
1
2
f
2
( )
x
f
j
n + 1
= f
j
n
Δ t
2 h
( f
j
n
2
− ( f
j − 1
n
2
( )
Δ t
2 h
Conservative form!
Conservative Method!
