<< . .

. 87
( : 95)



. . >>

 x0 + t
x0 + t(1 ’ x0 ) 0 ¤ x0 ¤ 1,
x(t) = x0 + tu0 (x0 ) =

x0 ≥ 1.
x0

Notice that the characteristic lines do not intersect only if t < 1 (see Figure 13.4,

right).



13.6 Systems of Linear Hyperbolic Equations
Consider the linear hyperbolic systems of the form

‚u ‚u
= 0 x ∈ R, t > 0, (13.30)
+A
‚t ‚x
where u : R — [0, ∞) ’ Rp and A ∈ Rp—p is a matrix with constant
coe¬cients.
This system is said hyperbolic if A is diagonalizable and has real eigen-
values, that is, if there exists a nonsingular matrix T ∈ Rp—p such that

A = TΛT’1 ,
600 13. Parabolic and Hyperbolic Initial Boundary Value Problems

where Λ = diag(»1 , ..., »p ) is the diagonal matrix of the real eigenvalues of
A, while T = (ω 1 , ω 2 , . . . , ω p ) is the matrix whose column vectors are the
right eigenvectors of A (see Section 1.7). The system is said to be strictly
hyperbolic if it is hyperbolic with distinct eigenvalues. Thus

Aω k = »k ω k , k = 1, . . . , p.

Introducing the characteristic variables w = T’1 u, system (13.30) becomes
‚w ‚w
+Λ = 0.
‚t ‚x
This is a system of p independent scalar equations of the form
‚wk ‚wk
+ »k = 0, k = 1, . . . , p.
‚t ‚x
Proceeding as in Section 13.5 we obtain wk (x, t) = wk (x ’ »k t, 0), and thus
the solution u = Tw of problem (13.30) can be written as
p
wk (x ’ »k t, 0)ω k .
u(x, t) =
k=1

The curve (xk (t), t) in the plane (x, t) that satis¬es xk (t) = »k is the k-th
characteristic curve and wk is constant along it. A strictly hyperbolic sys-
tem enjoys the property that p distinct characteristic curves pass through
any point of the plane (x, t), for any ¬xed x and t. Then u(x, t) depends
only on the initial datum at the points x ’ »k t. For this reason, the set of p
points that form the feet of the characteristics issuing from the point (x, t)

D(t, x) = x ∈ R : x = x ’ »k t , k = 1, ..., p , (13.31)

is called the domain of dependence of the solution u(x, t).
If (13.30) is set on a bounded interval (±, β) instead of on the whole real
line, the in¬‚ow point for each characteristic variable wk is determined by the
sign of »k . Correspondingly, the number of positive eigenvalues determines
the number of boundary conditions that can be assigned at x = ±, whereas
at x = β it is admissible to assign a number of conditions which equals the
number of negative eigenvalues. An example is discussed in Section 13.6.1.

Remark 13.2 (The nonlinear case) Let us consider the following non-
linear system of ¬rst-order equations
‚u ‚
+ g(u) = 0 (13.32)
‚t ‚x
where g = (g1 , . . . , gp )T is called the ¬‚ux function. The system is hyperbolic
if the jacobian matrix A(u) whose elements are aij = ‚gi (u)/‚uj , i, j =
1, . . . , p, is diagonalizable and has p real eigenvalues.
13.6 Systems of Linear Hyperbolic Equations 601

13.6.1 The Wave Equation
Consider the second-order hyperbolic equation

‚2u ‚2u
’ γ2 2 = f x ∈ (±, β), t > 0, (13.33)
‚t2 ‚x
with initial data
‚u
(x, 0) = v0 (x), x ∈ (±, β),
u(x, 0) = u0 (x) and
‚t
and boundary data

u(±, t) = 0 and u(β, t) = 0, t > 0. (13.34)

In this case, u may represent the transverse displacement of an elastic
vibrating string of length β ’±, ¬xed at the endpoints, and γ is a coe¬cient
depending on the speci¬c mass of the string and on its tension. The string
is subject to a vertical force of density f .
The functions u0 (x) and v0 (x) denote respectively the initial displace-
ment and the initial velocity of the string.
The change of variables

‚u ‚u
ω1 = , ω2 = ,
‚x ‚t
transforms (13.33) into the following ¬rst-order system

‚ω ‚ω
= 0 x ∈ (±, β), t > 0, (13.35)
+A
‚t ‚x
where

’1
ω1 0
ω= , A= ,
’γ 2
ω2 0

and the initial conditions are ω1 (x, 0) = u0 (x) and ω2 (x, 0) = v0 (x).
Since the eigenvalues of A are the two distinct real numbers ±γ (rep-
resenting the propagation velocities of the wave) we conclude that system
(13.35) is hyperbolic. Moreover, one boundary condition needs to be pre-
scribed at every end-point, as in (13.34). Notice that, also in this case,
smooth solutions correspond to smooth initial data, while any discontinuity
that is present in the initial data will propagate along the characteristics.

2
‚2u
Remark 13.3 Notice that replacing ‚ 2 by t2 ,
u
by x2 and f by 1, the
‚x2
‚t
wave equation becomes
t 2 ’ γ 2 x2 = 1
602 13. Parabolic and Hyperbolic Initial Boundary Value Problems

which represents an hyperbola in the (x, t) plane. Proceeding analogously
in the case of the heat equation (13.1), we end up with

t ’ νx2 = 1

which represents a parabola in the (x, t) plane. Finally, for the Poisson
2 2
equation (12.90), replacing ‚ u by x2 , ‚ u by y 2 and f by 1, we get
2 ‚y 2
‚x


x2 + y 2 = 1

which represents an ellipse in the (x, y) plane.
Due to the geometric interpretation above, the corresponding di¬erential
operators are classi¬ed as hyperbolic, parabolic and elliptic.



13.7 The Finite Di¬erence Method for Hyperbolic
Equations
Let us discretize the hyperbolic problem (13.26) by space-time ¬nite dif-
ferences. To this aim, the half-plane {(x, t) : ’∞ < x < ∞, t > 0} is dis-
cretized by choosing a spatial grid size ∆x, a temporal step ∆t and the
grid points (xj , tn ) as follows

xj = j∆x, j ∈ Z, tn = n∆t, n ∈ N.

Let us set
» = ∆t/∆x,
and de¬ne xj+1/2 = xj + ∆x/2. We look for discrete solutions un which
j
n
approximate the values u(xj , t ) of the exact solution for any j, n.
Quite often, explicit methods are employed for advancing in time in
hyperbolic initial-value problems, even though they require restrictions on
the value of », unlike what typically happens with implicit methods.
Let us focus our attention on problem (13.26). Any explicit ¬nite-di¬erence
method can be written in the form

un+1 = un ’ »(hn
j+1/2 ’ hj’1/2 ),
n
(13.36)
j
j


j+1/2 = h(uj , uj+1 ) for every j and h(·, ·) is a particular function
where hn n n

that is called the numerical ¬‚ux.


13.7.1 Discretization of the Scalar Equation
We illustrate several instances of explicit methods, and provide the corre-
sponding numerical ¬‚ux.
13.7 The Finite Di¬erence Method for Hyperbolic Equations 603

1. Forward Euler/centred

»
un+1 = un ’ a(un ’ un ) (13.37)
j j+1 j’1
j
2

which can be cast in the form (13.36) by setting

1
hj+1/2 = a(uj+1 + uj ). (13.38)
2

2. Lax-Friedrichs
1n »
(uj+1 + un ) ’ a(un ’ un )
un+1 = (13.39)
j’1 j+1 j’1
j
2 2

which is of the form (13.36) with

1
[a(uj+1 + uj ) ’ »’1 (uj+1 ’ uj )].
hj+1/2 =
2

3. Lax-Wendro¬

»2 2 n
»
’ a(uj+1 ’ uj’1 ) + a (uj+1 ’ 2un + un ) (13.40)
un+1 un n n
= j j j’1
j
2 2

which can be written in the form (13.36) provided that

1
[a(uj+1 + uj ) ’ »a2 (uj+1 ’ uj )].
hj+1/2 =
2

4. Upwind (or forward Euler/uncentred)

» »
un+1 = un ’ a(un ’ un ) + |a|(un ’ 2un + un ) (13.41)
j j+1 j’1 j+1 j j’1
j
2 2

which ¬ts the form (13.36) when the numerical ¬‚ux is de¬ned to be

1
[a(uj+1 + uj ) ’ |a|(uj+1 ’ uj )].
hj+1/2 =
2

The last three methods can be obtained from the forward Euler/centred
method by adding a term proportional to a numerical di¬usion, so that
they can be written in the equivalent form

1 un ’ 2un + un
» j+1 j j’1
’ a(uj+1 ’ uj’1 ) + k
un+1 un n n
= , (13.42)
j
j 2
2 2 (∆x)

where the arti¬cial viscosity k is given for the three cases in Table 13.1.
604 13. Parabolic and Hyperbolic Initial Boundary Value Problems

hdif f
methods k „ (∆t, ∆x)
j+1/2
∆x2
1
’ (uj+1 ’ uj ) O
∆x2
Lax-Friedrichs + ∆t + ∆x
2» ∆t
»a2
’ (uj+1 ’ uj ) O ∆t2 + ∆x2
a2 ∆t2
Lax-Wendro¬
2
|a|
|a|∆x∆t ’ (uj+1 ’ uj ) O(∆t + ∆x)
Upwind
2
TABLE 13.1. Arti¬cial viscosity, arti¬cial ¬‚ux and truncation error for
Lax-Friedrichs, Lax-Wendro¬ and Upwind methods


As a consequence, the numerical ¬‚ux for each scheme can be written equiv-
alently as
hj+1/2 = hF E + hdif f
j+1/2 j+1/2

where hF E is the numerical ¬‚ux of the forward Euler/centred scheme
j+1/2
(which is given in (13.38)) and the arti¬cial di¬usion ¬‚ux hdif f is given
j+1/2
for the three cases in Table 13.1.

An example of an implicit method is the backward Euler/centred scheme

»
a(un+1 ’ un+1 ) = un .
un+1 + (13.43)
j
j j+1 j’1
2
It can still be written in the form (13.36) provided that hn is replaced by
hn+1 . In the example at hand, the numerical ¬‚ux is the same as for the
Forward Euler/centred method, and so is the arti¬cial viscosity.

Finally, we report the following schemes for the approximation of the
second-order wave equation (13.33):

1. Leap-Frog

un+1 ’ 2un + un’1 = (γ»)2 (un ’ 2un + un ) (13.44)
j j+1 j j’1
j j


2. Newmark

un+1 ’ un = ∆tvj + (γ»)2 βwj + ’ β wj ,
n+1 1
n n
j
j 2
(13.45)
(γ»)2
’ θwj + (1 ’ θ)wj
n+1 n+1 n
n
vj vj =
∆t

with wj = uj+1 ’2uj +uj’1 and where the parameters β and θ satisfy
0 ¤ β ¤ 1 , 0 ¤ θ ¤ 1.
2
13.8 Analysis of Finite Di¬erence Methods 605

13.8 Analysis of Finite Di¬erence Methods
Let us analyze the properties of consistency, stability and convergence, as
well as those of dissipation and dispersion, of the ¬nite di¬erence methods
introduced above.


13.8.1 Consistency
As illustrated in Section 11.3, the local truncation error of a numerical
scheme is the residual that is generated by pretending the exact solution
to satisfy the numerical method itself.
Denoting by u the solution of the exact problem (13.26), in the case of
method (13.37) the local truncation error at (xj , tn ) is de¬ned as follows

u(xj , tn+1 ) ’ u(xj , tn ) u(xj+1 , tn ) ’ u(xj’1 , tn )
’a
n
„j = .
∆t 2∆x
The truncation error is

„ (∆t, ∆x) = max|„j |.
n
j,n

When „ (∆t, ∆x) goes to zero as ∆t and ∆x tend to zero independently,
the numerical scheme is said to be consistent.
Moreover, we say that it is of order p in time and of order q in space (for
suitable integers p and q), if, for a su¬ciently smooth solution of the exact
problem, we have
„ (∆t, ∆x) = O(∆tp + ∆xq ).
Using Taylor™s expansion conveniently we can characterize the truncation
error of the methods previously introduced as indicated in Table 13.1. The
Leap-frog and Newmark methods are both second order accurate if ∆t =
∆x, while the forward (or backward) Euler centred method is O(∆t+∆x2 ).
Finally, we say that a numerical scheme is convergent if

max|u(xj , tn ) ’ un | = 0.
lim j
∆t,∆x’0 j,n



13.8.2 Stability
A numerical method for a hyperbolic problem (linear or nonlinear) is said
to be stable if, for any time T , there exist two constants CT > 0 (possibly
depending on T ) and δ0 > 0, such that

¤ CT u0
un ∆, (13.46)


for any n such that n∆t ¤ T and for any ∆t, ∆x such that 0 < ∆t ¤ δ0 ,
0 < ∆x ¤ δ0 . We have denoted by · ∆ a suitable discrete norm, for
606 13. Parabolic and Hyperbolic Initial Boundary Value Problems

instance one of those indicated below
« p
1

= ∆x |vj |p  = sup|vj |.
v p = 1, 2, v (13.47)
∆,p ∆,∞

<< . .

. 87
( : 95)



. . >>