Numerical Solutions to Partial Differential Equations Zhiping Li - - PowerPoint PPT Presentation
Numerical Solutions to Partial Differential Equations Zhiping Li - - PowerPoint PPT Presentation
Numerical Solutions to Partial Differential Equations Zhiping Li LMAM and School of Mathematical Sciences Peking University More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Evolution Equations and
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Evolution Equations and Their FDM
Initial-Boundary Value Problems of Evolution Equations
∂u ∂t (x, t) = L(u(x, t)) + f (x, t), ∀ (x, t) ∈ Ω × (0, tmax], g(u(x, t)) = g0(x, t), ∀ (x, t) ∈ ∂Ω1 × (0, tmax], u(x, 0) = u0(x), ∀ x ∈ Ω, where L(·) is a (linear) differential operator acting on u with respect to x, and is assumed not explicitly depend on the time t.
2 / 30
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Evolution Equations and Their FDM
Initial-Boundary Value Problems of Evolution Equations
Definition An initial value problem is said to be well posed with respect to the norm · of a Banach space X, if it holds
1 for any given initial data u0 ∈ X, i.e. u0 < ∞, there exists
a solution;
2 there is a constant C > 0, such that, if v, w are the solutions
- f the problem with initial data v0, w0 ∈ X respectively, then
v(·, t) − w(·, t) ≤ Cv0(·) − w0(·), ∀ t ∈ [0, tmax]. Remark: Similarly, we can define the well-posedness for the initial-boundary value problems.
3 / 30
Uniformly Well Conditioned Difference Schemes
We consider the difference scheme of the following form B1Um+1 = B0Um + F m.
1 The difference operators B1 and B0 are independent of m. 2 (BαUm)j′ = j∈JΩ bα j′,jUm j ,
∀ j′ ∈ JΩ, α = 0, 1.
3 F m contains information on the inhomogeneous boundary
conditions as well as the source term of the difftl. eqn..
4 B1 = O(τ −1), B1 is invertible, and B−1 1
is uniformly well conditioned, i.e. there exists a constant K > 0, such that B−1
1 ≤ Kτ. 5 Under the above assumptions, we can rewrite the difference
scheme as Um+1 = B−1
1
[B0Um + F m] .
Remark: In vector case of dimension p, Um
j
∈ Rp, and bα
j′,j ∈ Rp×p.
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Consistency, Order of Accuracy and Convergence
Truncation Error of Difference Schemes
Definition Suppose that u is an exact solution to the problem, define T m := B1um+1 − [B0um + F m] , as the truncation error of the scheme for the problem.
Remark: By properly scaling the coefficients of the scheme (so that F m = f m), the definition is consistent with T m := {[B1(∆+t + 1) − B0] − [∂t − L]}um.
5 / 30
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Consistency, Order of Accuracy and Convergence
Consistency of Difference Schemes
Definition The difference scheme is said to be consistent with the problem, if for all sufficiently smooth exact solution u of the problem, the truncation error satisfies T m
j
→ 0, as τ(h) → 0, ∀mτ ≤ tmax, ∀j ∈ JΩ. In particular, the difference scheme is said to be consistent with the problem in the norm · , if τ
m−1
- l=0
T l → 0, as τ(h) → 0, ∀mτ ≤ tmax.
6 / 30
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Consistency, Order of Accuracy and Convergence
Order of Accuracy of Difference Schemes
Definition The difference scheme is said to have order of accuracy p in τ and q in h, if p and q are the largest integers so that |T m
j | ≤ C [τ p + hq] ,
as τ(h) → 0, ∀mτ ≤ tmax, ∀j ∈ JΩ, for all sufficiently smooth solutions to the problem, where C is a constant independent of j, τ and h.
Remark: Similarly, we can define the order of accuracy of a scheme with respect to a norm · , if τ
m−1
- l=0
T l ≤ C [τ p + hq] , as τ(h) → 0, ∀mτ ≤ tmax.
7 / 30
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Consistency, Order of Accuracy and Convergence
Convergence of Difference Schemes
Definition The difference scheme is said to be convergent in the norm · with respect to the problem, if the difference solution Um to the scheme satisfies Um − um → 0, as τ(h) → 0, mτ → t ∈ [0, tmax], for all initial data u0 with which the problem is well posed in the norm · .
Remark: Similarly, we can define the order of convergence of a scheme to be p in τ and q in h with respect to the norm · , if Um − um ≤ C [τ p + hq] , as τ(h) → 0, mτ → t ∈ [0, tmax].
8 / 30
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Stability and Lax Equivalence Theorem
Lax-Richtmyer Stability of Finite Difference Schemes
Definition A finite difference scheme is said to be stable with respect to a norm · and a given refinement path τ(h), if there exists a constant K1 > 0 independent of τ, h, V 0 and W 0, such that V m − W m ≤ K1V 0 − W 0, ∀ mτ ≤ tmax, as long as V m − W m is a solution to the homogeneous difference scheme (meaning F m = 0 for all m) with initial data V 0 − W 0.
9 / 30
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Stability and Lax Equivalence Theorem
Stability and Uniform Well-Posedness of Finite Difference Schemes
The stability of finite difference schemes are closely related to the uniform well-posedness of the corresponding discrete problems. For linear problems, since V m − W m = B−1
1 B0(V m−1 − W m−1),
thus, a scheme is (Lax-Richtmyer) stable if and only if
- B−1
1 B0
m ≤ K1, ∀ mτ ≤ tmax.
Remark: In general, the uniform well-posedness of the scheme with respect to the boundary data and source term can be derived from the Lax-Richtmyer stability and the uniform invertibility of the scheme.
10 / 30
Lax Equivalence Theorem
We always assume below that, for any initial data u0 with which the corresponding problem is well posed, there exists an sequence
- f sufficiently smooth solutions vα, s.t. limα→∞ v0
α − u0 → 0.
Theorem For a uniformly solvable linear finite difference scheme which is consistent with a well-posed linear evolution problem, the stability is a necessary and sufficient condition for its convergence.
1 The original continuous problem must be well-posed. 2 Linear and linear, can’t emphasize more! Not true if nonlinear. 3 The consistency is also a crucial condition. 4 Do not forget uniform solvability of the scheme. 5 ”Stability ⇔ Convergence” (not hold if miss any condition).
Proof of Sufficiency of the Lax Equivalence Theorem (stability ⇒ convergence)
If u ∈ X is a sufficiently smooth solution of the problem, then, by the definition of the truncation error, we have B1(Um+1 − um+1) = B0(Um − um) − T m,
- r equivalently Um+1 − um+1 = (B−1
1 B0)(Um − um) − B−1 1 T m.
Recursively, and by assuming U0 = u0, we obtain Um − um = −
m−1
- l=0
(B−1
1 B0)lB−1 1 T m−l−1.
Thus, by the uniform solvability B−1
1 ≤ Kτ and the stability
(B−1
1 B0)l ≤ K1, we have Um − um ≤ KK1τ m−1
l=0 T l, ∀m > 0.
Therefore, by the definition of the consistency, we have lim
τ(h)→0 Um − um = 0, 0 ≤ mτ ≤ tmax.
Proof of Sufficiency of the Lax Equivalence Theorem (continue)
For a general solution u, let vα be the smooth solution sequence satisfying limα→∞ v0
α − u0 → 0. 1 ∀ε > 0, ∃A > 0, such that v0 α − u0 < ε for all α > A. 2 For fixed β > A, let V m β be the solution of the difference
scheme with V 0
β = v0 β. 3 ∀ε > 0, ∃h(ε) > 0, s.t. V m β − vm β < ε, for all h < h(ε).
Thus, by the stability and the uniform invertibility of the scheme and the well-posedness of the problem that, if h < h(ε), then
Um − um ≤ Um − V m
β + V m β − v m β + v m β − um ≤ (K1 + 1 + C)ε.
Since ε is arbitrary, this implies lim
τ(h)→0 Um − um = 0, 0 ≤ mτ ≤ tmax.
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Stability and Lax Equivalence Theorem
Proof of Necessity of the Lax Equivalence Theorem (convergence ⇒ stability)
1 (B−1 1 B0)m h : bounded linear operators in (X, · ). 2 (B−1 1 B0)m h ≻ (B−1 1 B0) ˆ m ˆ h , either h < ˆ
h, or h = ˆ h and m > ˆ m.
3 By the resonance theorem, if, for each given u0 ∈ X,
Stmax(u0) := sup
{h>0, m≤τ −1tmax}
- (B−1
1 B0)m h u0
- < ∞,
then the sequence {(B−1
1 B0)m h } is uniformly bounded, and
consequently the scheme is stable. Suppose for some u0 ∈ X, ∃τk → 0, mkτk → t ∈ [0, tmax], s.t. lim
k→∞ (B−1 1 B0)mk hk u0 = ∞.
14 / 30
More on Consistency, Stability and Convergence Consistency, Stability and Convergence of FDM Stability and Lax Equivalence Theorem
Proof of Necessity of the Lax Equivalence Theorem (continue)
4 u, U with initial data u0, w and W with initial data w0 = 0.
(B−1
1 B0)m h u0 = Um h −W m h = (Um h −um)+(um−wm)+(wm−W m h ). 5 Convergence ⇒
lim
k→∞
- Umk
hk − umk + wmk − W mk hk
- = 0.
6 Well-posedness ⇒ um − wm ≤ Cu0 − w0 = Cu0. 7 Therefore lim k→∞ (B−1 1 B0)mk hk u0 ≤ Cu0. (A contradiction).
15 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM In What Norm Should the Stability be Analyzed
How to Choose a Norm for Stability Analysis
1 Maximum principle ⇒ · ∞ a good candidate. 2 Discrete maximum principle ⇒ · ∞ stability. 3 Constant-coefficient equation, periodic boundary conditions
⇒ · 2 a good choice.
4 Fourier analysis method ⇒ · 2 stability conditions. 5 Local Fourier analysis ⇒ necessary · 2 stability conditions. 6 More general cases ⇒ · 1, · p, energy norm, e.t.c..
16 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM In What Norm Should the Stability be Analyzed
Small Disturbance on a Scheme Does Not Change Its Stability
Theorem Suppose Um+1 = B−1
1 B0Um is stable in · . Suppose E(τ) is
uniformly bounded in · . Then, the scheme Um+1 =
- B−1
1 B0 + τE(τ)
- Um
is also stable in the norm · .
17 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM In What Norm Should the Stability be Analyzed
Small Disturbance on a Scheme Does Not Change Its Stability
Proof: (A + B)2 = A2 + AB + BA + B2, in general (A + B)m has 2m terms, there are C m
j
terms having j of B and (m − j) of A multiplied together in various order. (B−1
1 B0)k ≤ K1, for all k, and E(τ) ≤ K2.
Therefore, if mτ ≤ tmax, we have
- B−1
1 B0 + τE(τ)
m ≤
m
- j=0
C m
j K j+1 1
(τK2)j = K1(1 + K1K2τ)m ≤ K1etmaxK1K2.
- 18 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Fourier Analysis and von Neumann Condition of L2 Stability
Fourier Analysis in Constant-Coefficient and Periodic Case (Ω = (−1, 1)n)
1 Coefficients of a scheme: bα j′,j = bα (j−j′) (see (4.1.3)). 2 Fourier mode for a one-step scalar scheme: Um j
= λm
k eiπk·jh,
amplification factor: λk =
- j b0
(j−j′)eiπk·(j−j′)h
- j b1
(j−j′)eiπk·(j−j′)h =
ˆ B0(k) ˆ B1(k) .
3 In vector cases, apply the discrete Fourier transform:
ˆ U(k) = 1 ( √ 2N)n
N
- j1=−N+1
· · ·
N
- jn=−N+1
Uje−iπk·j 1
N ,
discrete inverse Fourier transform (i.e. Fourier expansion): Uj = 1 ( √ 2)n
N
- k1=−N+1
· · ·
N
- kn=−N+1
ˆ U(k)eiπk·j 1
N . 19 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Fourier Analysis and von Neumann Condition of L2 Stability
Fourier Analysis in Vector Case (Constant-Coefficient and Periodic)
4 To obtain the characteristic equations, substitute the Fourier
mode Um
j (k) = ˆ
Um(k)eiπk·j 1
N into the homogeneous scheme
- j b1
j′,j ˆ
Um+1(k)eiπk·j 1
N =
j b0 j′,j ˆ
Um(k)eiπk·j 1
N .
5 Denote ˆ
Bα =
j bα j′,jeiπk·(j−j′) 1
N =
j bα (j−j′)eiπk·(j−j′) 1
N .
6 We have ˆ
B1 ˆ Um+1(k) = ˆ B0 ˆ Um(k).
7 Characteristic equations in the space of frequencies:
ˆ Um+1(k) = ˆ B1(k)−1 ˆ B0(k)ˆ Um(k),
20 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Fourier Analysis and von Neumann Condition of L2 Stability
Fourier Analysis in Vector Case (Constant-Coefficient and Periodic)
9 Amplification matrix: G(k) = ˆ
B1(k)−1 ˆ B0(k) ∈ Cp×p,
10 Fourier mode: Um j
= λm
k eiπk·jh ˆ
U(k) with λk ˆ U(k) = G(k)ˆ U(k).
11 If G(k) is a normal matrix, i.e. ¯
G TG = G ¯ G T, then G has p
- rthogonal eigenvectors. Therefore Um2 = |λk|mU02.
12 In the general case, ˆ
Um(k) = [G(k)]m ˆ U0(k), therefore L2 stability ⇔ [G(k)]m 2 ≤ K1, ∀ k, ∀ mτ ≤ tmax.
21 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Fourier Analysis and von Neumann Condition of L2 Stability
The von Neumann Condition of L2 Stability
Theorem A necessary condition for the difference scheme to be L2 stable along a refinement path is that there exists a constant K ≥ 0 such that every eigenvalue λk of the amplification matrix G(k) satisfies |λk| ≤ 1 + Kτ, ∀ τ ≤ tmax, ∀ k.
22 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Fourier Analysis and von Neumann Condition of L2 Stability
The von Neumann Condition of L2 Stability
Proof:
1 L2 stable ⇒ for |λk|m ≤ K1 for m = [tmax/τ]. 2 Let f (x) = xs, (s < 1), since f ′′(x) = s(s − 1)xs−2 < 0, f (x)
is a concave function of x, hence f (x) ≤ f (1) + f ′(1)(x − 1).
3 K1 ≤ 1 ⇒ |λk| ≤ 1. K1 > 1 ⇒ |λk| ≤ K
1 m
1 ≤ 1 + 1 m(K1 − 1). 4 If 2τ ≤ tmax, m ≥ tmax
τ − 1 ≥ tmax 2τ .
5 The theorem holds for K = 2(K1 − 1)/tmax.
- 23 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Fourier Analysis and von Neumann Condition of L2 Stability
In General, the von Neumann Condition Is Not Sufficient for Stability
1 If G(k) is a normal matrix, then the von Neumann condition
is also a sufficient condition for a scheme to be L2 stable.
2 The conclusion is however not valid if G(k) is not normal. 3 In general, we can not expect to obtain a sufficient condition
by the local Fourier analysis (see Exercise 4.3).
4 In practical computations, the constant K can be a problem.
If K is too big, then the von Neumann condition can fail to provide useful information on the grid ratio.
24 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Strong Stability and Practical Stability Conditions
More Practical Conditions Than the von Neumann Condition Is Necessary
1 Initial value problem of ut + aux = cuxx, (c > 0). 2 Explicit central difference scheme:
Um+1
j
− Um
j
τ + a Um
j+1 − Um j−1
2h = c Um
j+1 − 2Um j + Um j−1
h2 .
3 |λk|2 = [1 − 2µ(1 − cos kh)]2 + [ν sin kh]2. 4 kh = π ⇒ µ ≤ 1 2 a necessary condition for the L2 stability. 5 µ ≤ 1 2
⇒ [1 − 2µ(1 − cos kh)]2 ≤ 1, ν2 = a2
c µτ ≤ a2 2c τ
⇒ |λk|2 ≤ 1 + a2
2c τ ⇒ the von Neumann condition is satisfied.
25 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Strong Stability and Practical Stability Conditions
More Practical Conditions Than the von Neumann Condition Is Necessary
6 In fact, µ ≤ 1 2 ⇔ L2 stability. 7 Recall |λk|2 = [1 − 2µ(1 − cos kh)]2 + [ν sin kh]2. 8 Take µ = 1 4, ν ≥ 1 (P´
eclet number ah
c = ν µ ≥ 4), 9 then, for kh = π/2, the amplification factor satisfies |λk| ≥ 5 4,
hence, the corresponding Fourier mode grows dramatically.
10 Yet, all Fourier mode solutions to the original problem decays.
(dissipation speed: e−ck2 ) 26 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Strong Stability and Practical Stability Conditions
A Practical Stability Condition for Difference Schemes Definition Let ˆ u(k, t) be the Fourier mode solutions to the differential equation, denote α = inf{β ≥ 0 : |ˆ u(k, t + τ)| ≤ eβτ|ˆ u(k, t)|, ∀ k}. If all of the amplification factors of the discrete Fourier mode solutions to the difference scheme satisfy |λk| ≤ eατ, ∀ k, then, the difference scheme is said to be practically stable, or to have practical stability. In particular, if α = 0, the corresponding difference scheme is said to be strongly stable, or to have strong stability.
27 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Strong Stability and Practical Stability Conditions
A Practical Stability Condition for Difference Schemes If α > 0, the discrete Fourier mode solutions can grow exponentially fast; however, the growth speed ≤ the fastest growth speed of real Fourier mode solution; the discrete Fourier modes are ”relatively” not in growth.
28 / 30
More on Consistency, Stability and Convergence More on the Stability of FDM Strong Stability and Practical Stability Conditions
A Practical Stability Condition for Difference Schemes If α = 0, the strong stability requires |λk| ≤ 1. In the above example, |λk|2 = [1 − 2µ(1 − cos kh)]2 + [ν sin kh]2. |λk| ≤ 1, ∀k ⇔ µ ≤ 1
2 and τ ≤ 2c a2 (see (3.4.10)).
Since ν2 = aτ
h
2 = a2
c µτ, τ ≤ 2c a2 ⇔ ν2 ≤ 2µ.
Take µ = 1/4, ν ≤ 1/ √ 2 ⇒ |λk|2 ≤ cos2 1
2kh ≤ 1, ∀k.
29 / 30