Buscar

CFD Lecture 08 Numerical Methods of Solving Equilibrium Problems(2015)

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes
Você viu 3, do total de 49 páginas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes
Você viu 6, do total de 49 páginas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes
Você viu 9, do total de 49 páginas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você também pode ser Premium ajudando estudantes

Prévia do material em texto

Numerical Methods for Elliptic PDE 1 
 
 
 
 
Numerical Methods for Solving Equilibrium Problem 
(Elliptic PDE) 
 
Numerical Methods for Elliptic PDE 2 
Consider 2-D Poisson equation of the form 
 
 
2 2
2 2
f f S
x y
∂ ∂
+ =
∂ ∂
, 0 , 0x a y b≤ ≤ ≤ ≤ (1) 
 
 
 xi-1/2 xi+1/2 
yj+1/2 
yj-1/2 
 
Numerical Methods for Elliptic PDE 3 
Using the 2nd-order central difference scheme, Eq. (1) can be discretized over a 
uniform mesh as 
 
 , -1 , , 1 -1, , 1, ,2 2
( 2 ) ( 2 )i j i j i j i j i j i j
i j
f f f f f f
S
y x
+ +− + − ++ =
∆ ∆
 (2a) 
 
In compact form, 
 
 2 2 2 21, , 1 , 1, , 1 ,2(1 )i j i j i j i j i j i jf f f f f S xβ β β− − + ++ − + + + = ∆ (2b) 
i = 2, …, (I-1), j = 2, …., (J-1) 
 
where /x yβ = ∆ ∆ with / ( 1)x a I∆ = − and /( 1)y b J∆ = − . 
 
Numerical Methods for Elliptic PDE 4 
Alternatively, we have 
 
 , 1, , , 1 , , , 1, , , 1 ,
W S P E N
i j i j i j i j i j i j i j i j i j i j i ja f a f a f a f a f b− − + ++ + + + = (2c) 
 
where 
 
2
, ,i j i jb S x= ∆ , , 1
W
i ja = , 
2
,
E
i ja β= , , 1
S
i ja = , 
2
,
N
i ja β= , and , , , , ,( )
P W E S N
i j i j i j i j i ja a a a a= − + + + . 
Numerical Methods for Elliptic PDE 5 
Methods of Solving System of Linear Algebraic Equation
 
Systems of linear algebraic equations can be expressed as 
 
 
1
N
ij j i
j
a f b
=
=∑ ( 1, , )i N=  or [ ]{ } { }A f B= (1) 
 
where [ ]A is the coefficient matrix 
 { }B is the non-homogeneous vector 
 
 
11 12 1
21 22 2
1 2
[ ]
N
N
N N NN
a a a
a a a
A
a a a
 
 
 =
 
 
 


   

, 
1
2{ }
N
f
f
f
f
 
 
 =  
 
  

, and 
1
2{ }
N
b
b
B
b
 
 
 =  
 
  

 
Numerical Methods for Elliptic PDE 6 
There are two basic approaches for solving system of linear algebraic equations: 
 
(a) Direct methods: 
 
 Cramer’s rule (matrix 
inversion) 
 Gauss elimination 
 LU decomposition 
(b) Iterative methods: 
 
 Jacobi iteration 
 Gauss-Seidel iteration 
 Relaxation 
 Conjugate Gradient Methods 
user
註解
沒有用
Numerical Methods for Elliptic PDE 7 
Direct Methods 
 
[LU Decomposition] 
 
The coefficient matrix [A] is decomposed or factorized into 
 
 [ ] [ ][ ]A L U= (1) 
where 
 
11
21 22
1 2
0 0
0
[ ]
N N NN
l
l l
L
l l l
 
 
 =
 
 
 


   

 and 
11 12 1
22 20[ ]
0 0
N
N
NN
u u u
u u
U
u
 
 
 =
 
 
 


   

 (2) 
 (Lower triangular matrix) (Upper triangular matrix) 
Numerical Methods for Elliptic PDE 8 
Now 
 [ ]{ } { }[ ][ ]{ }A f L U f B= = and { }-1 -1
[ ] { }
[ ] [ ][ ]{ } [ ]
I B
L L U f L B
′
=


 
 
⇒ [ ]{ } { }U f B′= ⇒ [ ]{ } { }
{B }
[ ][ ]{ }A f L U f B
′
= =

 ⇒ [ ]{ } { }L B B′ = 
 
It follows that 
 
[ ]{ } { } (Forward substitution)
[ ]{ } { } (Backward substitution)
L B B
U f B
′ =


 ′=
 
 
[LU Decomposition Algorithm] 
 
(a) Forward substitution: using 
[ ]{ } { }L B B′ = to solve for { }B′ . 
 
(b) Backward substitution: solve for 
{ }f by [ ]{ } { }U f B′= . 
Numerical Methods for Elliptic PDE 9 
Tri-Diagonal Matrix Algorithm (TDMA) (Thomas Algorithm) 
 
Consider a system of algebraic equations with tri-diagonal coefficient matrix 
 
1 1 1 1
2 2 2 2 2
-1 -1 -1 -1 -1
0 0 0 0
0 0
0
0 0 0
0 0
0 0
0 0 0 0
i i i i i
N N N N N
N N N N
b c f s
a b c f s
a b c f s
a b c f s
a b f s
     
     
     
     
     =    
    
    
    
         

 
       

      
 

 
or 
Numerical Methods for Elliptic PDE 10 
1 1 1 2 1
-1 1
-1
 (1)
, ( 2, 3, , 1) (2)
 (3)
i i i i i i i
N N N N N
b f c f s
a f b f c f s i N
a f b f s
+
+ =

 + + = = … −


+ =
 
 
We postulate the existence of two vectors { }D and { }E such that 
 
 1i i i if d f e+= + (4) 
 
⇒ -1 -1 -1i i i if d f e= + (5) 
 
Numerical Methods for Elliptic PDE 11 
Substituting Eq. (5) into Eq. (2) and solving for fi gives 
 
 -1 -1 1( )i i i i i i i i ia d f e b f c f s++ + + = ⇒ -1 1 -1( ) ( )i i i i i i i i ib a d f c f s a e++ = − + − 
 
⇒ -11
-1 -1
( )
( ) ( )
i i i i
i i
i i i i i i
i i
c s a ef f
b a d b a d
d e
+
− −
= +
+ +
 
 (6) 
 
From Eqs. (4) and (6) we have 
 
-1
-1
-1
 
( )
( ) 
( )
i
i
i i i
i i i
i
i i i
cd
b a d
s a ee
b a d
− = +


 −
 =
+
 2 i N≤ ≤ (7) 
Numerical Methods for Elliptic PDE 12 
Furthermore, to start the forward substitution of Eq. (7) the coefficients 1d and 
1e can be evaluated from Eq. (1) as 
 
 1 1 1 2 1b f c f s+ = ⇒ 1 11 2 1 2 1
1 1
( ) ( )c sf f d f e
b b
−
= + = + 
⇒ 
1
1
1
1
1
1
( )
( )
cd
b
se
b
 = −

 =

 (8) 
 
Similarly, from Eq. (3) 
 
 -1N N N N Na f b f s+ = ⇒ -1 ( ) ( )N NN N
N N
b sf f
a a
= − + (9) 
Numerical Methods for Elliptic PDE 13 
Further application of Eq. (5) for i = (N-1) and comparison with Eq. (9) gives 
 
 -1 -1 -1 ( ) ( )N NN N N N N
N N
b sf d f e f
a a
= + = − + ⇒ -1 -1[ ( )] [( ) ]N NN N N
N N
b sd f e
a a
+ = − 
 
⇒ -1
-1
( )
( )
N N N
N
N N N
s a ef
b a d
−
=
+
 (16) 
 
Next, the backward substitution can be implemented using Eqs. (5) and (16) as 
 
 1i i i if d f e+= + ( 1), ( 2), , 2, 1i N N= − − … (17) 
 
Numerical Methods for Elliptic PDE 14 
Furthermore, let 
 
-1
-1
i i i i
i i i i
b b a d
s s a e
′ = +


 ′ = −
 (18) 
 
Eq. (7) becomes then 
 
i
i
i
i
i
i
cd
b
se
b
 = − ′


 ′
 =
′
 2 i N≤ ≤ (19) 
Numerical Methods for Elliptic PDE 15 
Since 
1 1
1
1 1
1 1
1
1 1
 
c cd
b b
s se
b b
 = − = − ′


 ′
 = =
′
 ⇒ 1 1
1
b b
s s
′ =
 ′ =
 (20) 
Eq. (17) can be further expressed as 
 
 1( )i i ii
i
s c ff
b
+′ −=
′
 ( 1), ( 2), , 2, 1i N N= − − … (22) 
 
where 
 
-1
-1
-1
/i i i
i i i i
i i i i
r a b
b b rc
s s r s
′=
 ′ = −
 ′ ′= − (23) 
Numerical Methods for Elliptic PDE 16 
 
Likewise, Eq. (16) becomes 
 
 -1
-1
( )
( )
N N N N
N
N N N N
s a e sf
b a d b
′−
= =
′+
 (24) 
Numerical Methods for Elliptic PDE 17 
[TDMA Algorithm] 
 
(a) Initialize the variables: 1 1b b′ = and 1s s′ = 
(b) Recurrently calculate the following equations in increasing order of 
2,3, , i N= … : 
-1
-1
-1
/
 (Forward substitution)
i i i
i i i i
i i i i
r a b
b b rc
s s r s
′=
 ′ = −
 ′ ′= − 
 
(c) Calculate /N N Nf s b′ ′= 
(d) Calculate the following equations in decreasing order of i: 
 
 1( )i i ii
i
s c ff
b
+′ −=
′
 ( 1), ( 2), , 2, 1i N N= − − … . (Backward substitution)
Numerical Methods for Elliptic PDE 18 
FORTRAN Subroutine of TDMA 
 
 SUBROUTINE TRIDG(IL, IU, A, B, C, S) 
 DIMENSION A(1), B(1), C(1), S(1) 
C….Subroutine to perform solution algorithm for 
the tri-diagonal 
C…..system of algebraicequations 
C….. 
C…..IL = subscript of first equation 
C…..IU = subscript of the last equation 
C…..A = coefficient to the left of the main diagonal 
C…..B = coefficient of the main diagonal 
C…..C = coefficient to right of the main diagonal 
C…..S = Element in the constant vector 
C…… 
C…..Forward substitution to establish upper 
tri-diagonal matrix 
C….. 
 DO 10 I = LL+1, IU 
 R = A(I)/B(I−1) 
 B(I)= B(I)−R*C(I−1) 
 S(I) = S(I)−R*S(I−1) 
10 CONTINUE 
C….. 
C…..Backward substitution 
C…. 
 S(IU) = S(IU)/B(IU) 
 DO 20 I = IU−1, 1, −1 
 S(I) = (S(I)−C(I)*S(I+1))/B(I) 
20 CONTINUE 
RETURN 
 END 
 
 
 
Numerical Methods for Elliptic PDE 19 
Iterative Methods 
[Basic Concepts] 
Consider the system of algebraic equations of 
 
 [ ]{ } { }A f B= (1) 
 
Iteration Algorithm: 
(a) Assume an initial solution vector 0{ }f . 
(b) Generate an improved solution vector { }mf , m = 1, 2, …, based on some 
strategy for reducing the difference between the iterating solution vector 
{ }mf and the true solution vector { }f , i.e., 
 
Numerical Methods for Elliptic PDE 20 
After m iterations we have an appropriate solution vector, { }mf , there is a 
non-zero residual vector { }mR as 
 
 [ ]{ } { } { } {0}m mA f B R− = ≠ (2) 
 
Substrating Eq. (2) from Eq. (1) gives 
 
 [ ]({ } { }) { }m mA f f R− = (3) 
 
Let the convergence error vector { }mε defined as 
 
 { } { } { }m mf fε = − (4) 
 
Eq. (3) becomes then [ ]{ } { }m mA Rε = (5) 
 
 
Numerical Methods for Elliptic PDE 21 
(c) Repeat step (b) till convergence, i.e., 
 
{ } {0}mR → ⇒ { } {0}mε → ⇒ { } { }mf f→
user
螢光標示
Numerical Methods for Elliptic PDE 22 
Generalized Iterative Scheme 
 
The system of algebraic equations can be modified by a non-singular 
pre-conditioning matrix [ ]P as 
 
 [ ]{ } { }[ ]( )P A f B= or { }* *[ ] { }A f B= (6) 
 
where *[ ] [ ][ ]A P A= and *{ } [ ]{ }B P B= . 
 
Decomposing the coefficient matrix to 
 
 *[ ] [ ] [ ]A M N= − (7) 
 
 { } *([ ] [ ]) { }M N f B− = (8) 
 
Numerical Methods for Elliptic PDE 23 
An iterative scheme can then expressed as 
 
 { } { }1 *[ ] [ ] { }m mM f N f B+ = + (9) 
 
Alternatively, Eq. (9) can be rearranged as 
 
{ } { } { }1 *
+1
[ ]( ) { } ([ ] [ ])
{ } { }
m m m
m m
M f f B M N f
Rδ
+ − = − −
 
 
 
 +1[ ]{ } { }m mM Rδ = (10) 
 
where { } { }+1 1{ }m m mf fδ += − is the correction vector and is an approximation to 
the convergence vector { }mε . 
 
It follows that +1{ } {0}mδ → ⇒ { } {0}mR → ⇒ { } {0}mε → .
user
註解
讓error逼近0
Numerical Methods for Elliptic PDE 24 
[Jacobi Iteration Method] 
 
(a) Point Iteration 
Let 
 [P] = [I] and [A] = [L]+[D]+[U] (11) 
where 
 
1 0 0
0 1 0
[ ]
0 0 1
I
 
 
 =
 
 
 


   

 21
1 2
0 0 0
0 0
[ ]
0N N
a
L
a a
 
 
 =
 
 
 


   

 
 (Identity matrix) (Lower Triangular matrix) 
Numerical Methods for Elliptic PDE 25 
 
11
22
0 0
0 0
[ ]
0 0 NN
a
a
D
a
 
 
 =
 
 
 


   

 
12 1
2
0
0 0
[ ]
0 0 0
N
N
a a
a
U
 
 
 =
 
 
 


   

 
 (Diagonal matrix) (Upper Triangular matrix) 
 
By setting 
 
 
[ ] [ ]
[ ] ([ ] [ ])
M D
N L U
=


 = − +

 (12) 
 
Numerical Methods for Elliptic PDE 26 
The iteration formula takes then the form of 
 
 { } { }1 *[ ] { } ([ ] [ ])m mD f B L U f+ = − + (Point Jacobi iteration) (13) 
 
For the case of the Poisson equation, Eq. (2c), the Jacobi iteration formula can be 
expressed as 
 
 1, , , 1, , , 1 , 1, , , 1
,
1 [ ( )]m W m S m E m N mi j i j i j i j i j i j i j i j i j i jP
i j
f b a f a f a f a f
a
+
− − + += − + + + (14) 
 
 
Alternatively, Eq. (14) can be rearranged as 
 
 1, , , , 1, , , 1 , , , 1, , , 1
,
1 [ ( )]m m W m S m P m E m N mi j i j i j i j i j i j i j i j i j i j i j i j i jP
i j
f f b a f a f a f a f a f
a
+
− − + += + − + + + + 
Numerical Methods for Elliptic PDE 27 
 ,1, ,
,
m
i jm m
i j i j P
i j
R
f f
a
+ = + (15) 
 
with the residual at the grid (i, j) as 
 
 , , , 1, , , 1 , , , 1, , , 1[ ( )]
m W m S m P m E m N m
i j i j i j i j i j i j i j i j i j i j i j i jR b a f a f a f a f a f− − + += − + + + + (16) 
user
註解
每一個點都能同時進行
Numerical Methods for Elliptic PDE 28 
Point Iteration Algorithm 
(1) Re-order the equations to provide diagonal-dominant coefficient matrix such 
that the magnitude of the diagonal elements are larger than those of other 
elements in the same row, i.e. 
 ii ija a≥ ( )j i≠ 
(2) Re-arrange each equation as of the form 
 ,1, ,
,
m
i jm m
i j i j P
i j
R
f f
a
+ = + (17) 
 with the residual at the grid (i, j) as 
 , , , 1, , , 1 , , , 1, , , 1[ ( )]
m W m S m P m E m N m
i j i j i j i j i j i j i j i j i j i j i j i jR b a f a f a f a f a f− − + += − + + + + (18) 
where the superscript m (= 0, 1, …) denotes the level of iteration. 
Numerical Methods for Elliptic PDE 29 
(3) Assume an initial guess for 0{ }f . 
(4) Calculate new values of 1{ }mf + from the known value { }mf , (m ≥ 0) from Eq. 
(17). 
(5) Repeat Step (4) until a prescribed convergence criterion is satisfied. 
 
Numerical Methods for Elliptic PDE 30 
(b) Line Iteration 
 
To improve iteration convergence, the point Jacobi iteration formula, Eq. (14), 
may be rearranged as line Jacobi iteration formulas: 
 
 1 1 1, 1, , , , 1, , , , 1 , , 1( )
W m P m E m S m N m
i j i j i j i j i j i j i j i j i j i j i ja f a f a f b a f a f
+ + +
− + − ++ + = − + (i = 2,.., (I-1)) (19a) 
or 
 +1 1 +1, , 1 , , , , 1 , , 1, , 1,( )
S m P m N m W m E m
i j i j i j i j i j i j i j i j i j i j i ja f a f a f b a f a f
+
− + − ++ + = − + (j = 2,.., (J-1)) (19b) 
 
Then, a suitable mean, such as TDMA, is used to obtain the solution of 1{ }mf + 
along either ith column, Eq. (19a), or jth row, Eq. (19b). 
Numerical Methods for Elliptic PDE 31 
[Gauss-Seidel Iteration Method] 
 
(a) Point Iteration 
By setting 
 
[ ] ([ ] [ ])
[ ] [ ]
M D L
N U
= +


 = −

 (20) 
 
The iteration formula takes then the form of 
 
 { } { }1 *([ ]+[ ]) { } ([ ])m mD L f B U f+ = − (21) 
(Point Gauss-Seidel iteration) 
user
註解
修正jacob
Numerical Methods for Elliptic PDE 32 
For Eq. (2c), the Gauss-Seidel iteration formula can be expressed as 
 
 1 1 1, , , 1, , , 1 , 1, , , 1
,
1 [ ( )]m W m S m E m N mi j i j i j i j i j i j i j i j i j i jP
i j
f b a f a f a f a f
a
+ + +
− − + += − + + + (22) 
 
Alternatively, 
 
 1 +1 +1, , , , 1, , , 1 , , , 1, , , 1
,
1 [ ( )]m m W m S m P m E m N mi j i j i j i j i j i j i j i j i j i j i j i j i jP
i j
f f b a f a f a f a f a f
a
+
− − + += + − + + + + 
 
 ,1, ,
,
m
i jm m
i j i j P
i j
R
f f
a
+ = + (23) 
with the residual at the grid (i, j) as 
 
 +1 +1, , , 1, , , 1 , , , 1, , , 1[ ( )]
m W m S m P m E m N m
i j i j i j i j i j i j i j i j i j i j i j i jR b a f a f a f a f a f− − + += − + + + + (24) 
user
註解
該點有新值就用新的算沒有就愈舊的
Numerical Methods for Elliptic PDE 33 
(b) Line Iteration 
The line Gauss-Seidel iteration formulas take the forms of: 
 
 1 1 1 +1, 1, , , , 1, , , , 1 , , 1( )
W m P m E m S m N m
i j i j i j i j i j i j i j i j i j i j i ja f a f a f b a f a f
+ + +
− + − ++ + = − + (i = 2,.., (I-1)) (25a) 
 (Row iteration) 
or 
 
 +1 1 +1 +1, , 1 , , , , 1 , , 1, , 1,( )
S m P m Nm W m E m
i j i j i j i j i j i j i j i j i j i j i ja f a f a f b a f a f
+
− + − ++ + = − + (j = 2,.., (J-1)) (25b) 
 (Column iteration) 
 
Numerical Methods for Elliptic PDE 34 
[Relaxation Techniques] 
 
The idea of relaxation technique is to vary the convergence rate of a (point or line) 
iteration scheme by modifying the propagation of the correction 
1 1{ } { } { }m m mf f fδ + += − through the grid system as 
 
 
1
1 * * 1
{ }
{ } { } (1 ){ } { } ({ } { }) { } { }
m
m m m m m m
f
f f f f f f f f
δ
ω ω ω ω δ
+
+ +
=
= + − = + − = +

 (26) 
 
where 
 *{ }f denotes the newly calculated value from the iteration scheme
Numerical Methods for Elliptic PDE 35 
ω is the relaxation factor with 0 < ω < 2, 
 
 
0 1, (Under-relaxation) 
1, (No relaxation)
1 2, (Over-relaxation)
ω
ω
ω
< <
 =
 < <
 
 The relaxation technique may be further classified into 
 
 Stationary relaxation in which ω = constant through the iteration. 
 Non- Stationary relaxation in which ω ≠ constant through the iteration. 
 
 
Numerical Methods for Elliptic PDE 36 
(a) Point Relaxation 
The point relaxation scheme may be expressed as 
 
 { } { }* *[ ] [ ] { }mM f N f B= + and 1 *{ } { } (1 ){ }m mf f fω ω+ = + − 
 
⇒ { } { } { }1 * *[ ] [ ]( (1 ){ })=[ ] { }m m mM f M f f N f Bω ω+ = + − + 
 
{ } { } { }
{ } { }
1 * *
[ ]
* *
[ ] ( [ ] (1 )[ ]) { } ([ ] ([ ] [ ])) { }
([ ] [ ]) { } ((1 )[ ] [ ]) { }
m m m
A
m m
M f N M f B M M N f B
M A f B A N f B
ω ω ω ω ω
ω ω ω ω
+
=
= + − + = − − +
= − + = − + +

 
 
As a result, the point relaxation iteration scheme can be generally expressed as 
 
 { } { }1 *[ ] ((1 )[ ] [ ]) { }m mM f A N f Bω ω+ = − + + (27) 
Numerical Methods for Elliptic PDE 37 
Gauss-Seidel Relaxation Scheme 
 
Incorporating the point relaxation technique, the Gauss-Seidel point iteration 
scheme takes the form of 
 
 * 1 1, , , 1, , , 1 , 1, , , 1
,
1 [ ( )]W m S m E m N mi j i j i j i j i j i j i j i j i j i jP
i j
f b a f a f a f a f
a
+ +
− − + += − + + + (28a) 
and 
 1 *, , ,(1 )
m m
i j i j i jf f fω ω
+ = + − (28b) 
 
 For 1 < ω < 2, the Gauss-Seidel iteration scheme is known as a successive 
over-relaxation (SOR) scheme. 
 
user
螢光標示
user
螢光標示
Numerical Methods for Elliptic PDE 38 
[Ex] 
 
 
 
Numerical Methods for Elliptic PDE 39 
 
 
 
Numerical Methods for Elliptic PDE 40 
 
Numerical Methods for Elliptic PDE 41 
 
Numerical Methods for Elliptic PDE 42 
(b) Line Relaxation 
[Uni-Directional Line Relaxation] 
The Gauss-Seidel line relaxation formula takes first either row or column 
iterations as: 
 
 * * * +1, 1, , , , 1, , , , 1 , , 1( )
W P E S m N m
i j i j i j i j i j i j i j i j i j i j i ja f a f a f b a f a f− + − ++ + = − + (i = 2,.., (I-1)) (29a) 
or 
 * * * +1, , 1 , , , , 1 , , 1, , 1,( )
S P N W m E m
i j i j i j i j i j i j i j i j i j i j i ja f a f a f b a f a f− + − ++ + = − + (j = 2,.., (J-1)) (29b) 
 
Then 
 1 *, , ,(1 )
m m
i j i j i jf f fω ω
+ = + − (30) 
Numerical Methods for Elliptic PDE 43 
[Alternating Direction Line Relaxation] 
This scheme consists of line relaxation alternately along row and column 
directions as 
 
(1) Row Line Relaxation: 
 
 * * * +1, 1, , , , 1, , , , 1 , , 1( )
W P E S m N m
i j i j i j i j i j i j i j i j i j i j i ja f a f a f b a f a f− + − ++ + = − + (i = 2,.., (I-1)) (31a) 
 
 1 *, , ,(1 )
m m
i j i j i jf f fω ω
+ = + − (31b) 
 
(2) Column Line Relaxation: 
 
 * * * +2 +1, , 1 , , , , 1 , , 1, , 1,( )
S P N W m E m
i j i j i j i j i j i j i j i j i j i j i ja f a f a f b a f a f− + − ++ + = − + (j = 2,.., (J-1)) (32b) 
 
 2 * +1, , ,(1 )
m m
i j i j i jf f fω ω
+ = + − (32b) 
Numerical Methods for Elliptic PDE 44 
Convergence Criteria for Iterative Methods 
 
The convergence of an iterative process is achieved when the desired accuracy 
criterion (or criteria) is satisfied, i.e. 
 
+1{ } {0}mδ → ⇒ +1{ } {0}mε → 
 
where 
 
{ } { }+1 1{ }m m mf fδ += − is the correction vector 
 
+1 +1{ } { } { }m mf fε = − is the convergence vector 
 
Numerical Methods for Elliptic PDE 45 
Convergence criteria can be specified in terms of 
 
 Absolute error: 
 
 1 1, , , ,( )
m m m
i j a i j i j i jf f f fδ δ
+ += = − 
 
 Relative error: 
 
 
1 1
, , ,
, 1 1
, ,
( )
m m m
i j i j i j
i j r m m
i j i j
f f f
f
f f
δ
δ
+ +
+ +
−
= = 
 
Numerical Methods for Elliptic PDE 46 
[Absolute Convergence Criteria] 
 
(a) Infinity norm: 
 
 1, , max{ }
m m
a i j i jf f fδ ε
+
∞
= − ≤ 
 
(b) The 1st norm: 
 
 1, ,1{ }
m m
a i j i j
j i
f f fδ ε+= − ≤∑∑ 
 
(c) The 2nd norm: 
 
 1 2 1/ 2, ,2{ } [ ( ) ]
m m
a i j i j
j i
f f fδ ε+= − ≤∑∑ 
 
Numerical Methods for Elliptic PDE 47 
[Relative Convergence Criteria] 
 
(a) Infinity norm: 
 
1
, , max
1
, max
{ }
m m
i j i j
r m
i j
f f
f
f
δ ε
+
∞ +
−
= ≤ 
 
(b) The 1st norm: 
 
1
, ,
11
,
{ }
m m
i j i j
r m
i jj i
f f
f
f
δ ε
+
+
−
= ≤∑∑ 
 
(c) The 2nd norm: 
 
1
, , 2 1/ 2
12
,
{ } [ ( ) ]
m m
i j i j
r m
i jj i
f f
f
f
δ ε
+
+
−
= ≤∑∑ 
Numerical Methods for Elliptic PDE 48 
Comparison of the Methods for Solving Systems of Linear Algebraic Equations 
 
Method 
Approx. 
Max. No. of 
Equations 
Stability Precision Breadth of Applications 
Programming 
Effort Comments 
Cramer’s rule 3 - Affected by round-off error Limited - 
Excessive 
computational 
effort required 
for more than 
three equations 
LU 
Decomposition 100 - 
Affected by 
round-off error General Moderate 
Preferred 
elimination 
method 
Gauss-Seidel 1000 
May not 
converge if 
not 
diagonally 
dominant 
Excellent 
Appropriate 
only for 
diagonally 
dominant 
systems 
Easy 
 
Numerical Methods for Elliptic PDE 49 
 
 
 
END

Outros materiais