Transcript Chapter 4: Theoretical Formulation

Chapter 4
Theoretical
Background
4.1 Convergence
QUESTIONS?
 Under what conditions the numerical solution
will coincide with the exact solution?
 What guarantee the numerical solution will be
close to the exact solution of the PDE?
 Consistency – discretization solution process
can be reversed, through a Taylor series
expansion, to recover the original PDEs
 Stability – solution algorithm must be stable
 Convergence – solution of the discretization
equations (approximate solution) approaches
the exact solution for the original PDE when
x  0 and t  0 (grid refinement)
Numerical Convergence




Convergence -- Numerical solution approaches
the exact solution of PDE for each value of the
independent variable
Solution error -- truncation and round-off errors
Numerical Convergence -- Numerical solution
converges to a unique solution with grid
refinement (very expensive to prove)
The converged solution may not be the exact
solution unless the numerical scheme is also
consistent
n
Tjkl
 T ( x j , yk , zl , t n ) as Δx,Δy,Δz,Δt  0
Numerical Convergence

FTCS Scheme
tmax = 5000
Grid refinement is a very expensive process
(especially for unsteady 3D problems)
Numerical Convergence

FTCS
Scheme
Lax Equivalence Theorem

The necessary and sufficient condition for convergence
of a properly posed linear initial value problem
-- Satisfies the consistency condition
-- The algorithm is stable

Applicable to any discretization procedure (not only
finite-differences) that leads to nodal unknowns

For nonlinear boundary value problems or mixed initial
/boundary value problems, Lax equivalence theorem
cannot be rigorously applied. It may be considered as a
necessary, but not always sufficient condition
Consistency + Stability = Convergence
Convergence

Discretization – Replace (approximate)
the PDE by algebraic equations
L(T )  0


La (T )  0
Consistency – Recover (x, t  0)
PDE from algebraic equations
La (T )  0
 L(T )  0 ?
Convergence
L(T )  0
Partial differential
equations (PDE)
Exact
solution
Discretization
Consistency
La (T )  0
System of algebraic
equations
Stability
T T
n
T  T ( x, y, z , t )
Convergence
n
jkl
T jkl
 T ( x j , y k , z l , t n ) x , y , z , t  0
Approximate solution
Exact solution
4.2 Consistency




The discretization (algebraic) equation is
consistent with the original PDE if the two
are equivalent at each grid point as x, y,
z, t 0
But the exact or converged solutions are
unknown!
Consistency is a necessary, but not
sufficient condition
The algorithm must also be STABLE to
achieve CONVERGENCE
Given La (T ), can we recoverL(T )  0?
Numerical Consistency



Use Taylor series expansion and examine the remainder
Substitute the exact solution into the discretization
equation, and compare with the original PDE
The numerical solution  satisfies the discretization
equation exactly (assuming no round-off error), but not
the original PDE
 L(T )  0

 La (T )  0
PDE
Alge braic
Truncation error analysis
 L(T )  0
but 
 La (T )  0
La (T )  L(T )  E (T )  E (T )
Modified equation approach La (T )  L(T )  E (T )  0
4.2.1 FTCS Scheme
T
 2T
L(T ) 

0
Parabolic
2
t
x
T jn  1  T jn
T jn 1  2T jn  T jn 1
La (T ) 

0
2
t
x
t
n 1
n
n
n
T j  sT j  1  (1  2 s )T j  sT j  1 ; s 
x 2
n+1
1
n

s
j1
12s
s
j+1
j
Include both marching and equilibrium behaviors
FTCS Scheme
 Taylor
T jn  1
series expansion
 T
n
 T j  t  
 t


t 
 
2!
j
n
2
 T 
x 
 
T jn 1  T jn  x  
2!
 xj
n
2
  2T
 2
t
 T

2

x

2
3


t 
 
3!
j
n
  3T
 3
t

x 
 
3!
j
n
3
n

  
j
 T 
x 



3 
4!

x

j
3
n
4

1  n 1
  t  n
n
n
n
La (T ) 
 T j  T j   2  T j  1  2T j  T j  1 
t
 x 

n
n
n
2
3
2
3





 T 
t    T   t    T   
1 






t


 t 
  t2 
  t3 
t
2
!
3
!

j

j

j


2





t

x


 2 2  
 x   2!
 T

2

x

2


x 
 
4!
j
n
4

 T

4

x

4


x 
 
6!
j
n
6
 T

6

x

6
 T

4

x

4
n

  
j
 t
s
 x2


  

j
n



FTCS Scheme
 Taylor
series expansion
n
 T
 T 
n


La (T )  


E
j T 
2 
 x j
 t
2
E nj T  
t   T
 2
2 t
2
n

t   T
 
 3
6 t
j
2
3
n
  x   T
 

4
12

x
j

2
4
n
  x   T
 

6
360

x
j

4
 Consistency
E nj T   0 as x, t  0;  La T   LT 
 Recovers
the original PDE!
 Convergence?
6
n

  
j
Numerical Accuracy



Time- and spatial-derivatives are not independent
We need to know the accuracy of the discretization
equation, not just individual terms in the equation
Use the PDE to relate time- and spatial-derivatives
 T
 2T


 x2
t
 2
   2T
 T
 
 2 
2

t

t

x


 3
  6T 

T
3


  
3
6 
  t
 x 

2
  
2

x

 T

t

2
  
2

x

4
  2T 
2  T
 
   
2 
4

x

x



2m

 T
T
m 
  
m
2m
t

x

m






FTCS Scheme
 Convert
to pure time-derivative
or pure spatial-derivative
2


t

T
 2
E nj T  
2 t
n

t 2   3 T
 
 3
6 t
j
  t  x   T
 2
  
 2 12   t
2
2
t 
1  T

 1    2
2 
6s    t
2

2
n
  t
x
  

360 2
j  6
2
n
4

t 
1
 
1 
2
6
60
s

j
1  T
 s   
6    x4

 x 
2
n

 x 2   4 T
 

4
12

x
j

4
2
n
n

 x 4   6 T
 

6
360

x
j

 T
  3
  t
3
 T
  3
  t
3
n

  
j
n

  
j

 x  2 1    T
 
 s   
6
6
60

x


j

4
n

  
j
6
n
 t
s
2
x

  
j
FTCS Scheme

s = t /x2 is dimensionless
E T  
n
j
1  T
 s   
6    x4

 x 
2
2
4
  t  x   T
 2
  
 2 12   t
2
Time
derivative

s 




s 


2
n
n

 x  2 1    T
 
 s   
60 
60    x 6
j
4
6
  t
x
  

2
6
360

j 
Spatial
derivative
2
4

  
j
n

 
j
Error Cancellation !
4

1


x
1

T


 E nj T  
 s   
6
2 
6    x4
2
n


j
1
 x  2 1    T
 E nj T  
 s   
6
6 
60    x 6
4
 T
  3
  t
3
n
6
n

 x   T
 

6
540

x
j

4
6
n


j
4.2.2 Fully Implicit Scheme
T
 2T
L(T ) 

0
2
t
x
T jn  T jn 1
T jn 1  2T jn  T jn 1
La (T ) 

0
2
t
x
n 1
n 1
n 1
n

sT

(
1

2
s
)
T

sT

T
Level n+1
j 1
j
j 1
j
Level n
 sT jn 1  (1  2 s )T jn  sT jn 1  T jn 1
s
j1
1+2s
1
j
s
n
n-1
j+1
Fully Implicit Scheme
 Taylor
T jn  1
series expansion

t 
 
2!
j
n
 T
 T jn  t  
t
 T 
x 
 
 T jn  x  
2!
 xj
n
T jn 1
2
La (T ) 
2
  2T
 2
t
3


t 
 
3!
j
  2T

2
 x
n
  3T
 3
t
3


x 
 
3!
j
n
n

  
j
4
  3T 

x 

 
3 
4!
 x j
n
  4T

4
 x
n

  
j

1  n
  t  n
n1
n
n
T

T

T

2
T

T

 j 1
 j
j
j
j 1 
2
t 
 x 

 T
1 

 t  
t 
 t




t 
 
2!
j
2
  t   x 
 2

2
 x   2!
n
  2T

2

x

2
 T
 2
t
2

t 
 
3!
j
4


x 
 
4!
j
n

n
  4T

4

x

3
 T
 3
t
3
n

  
j
6


x 
 
6!
j
n
  6T

6

x

n


  

j



Fully Implicit Scheme

Convert to pure time- or spatial-derivative
E nj T   
t   T
 2
2 t
2
n

t   T
 
 3
6 t
j
3
n
n

 x   T
 

4
12

x
j

2
4
  t  x   T
 2
  

 2 12   t
  t
x
  

360 2
j  6
t 
1    2T

  1    2
2 
6s    t

t 2 
1
 
1


2
6
60
s

j
2


2
2
2
4
n
1  T
  s   
6    x4

 x 
2
2
4
n
n

 x   T
 

6
360

x
j

4

  
j
3
 T
  3
  t

  
j
n

 x  2 1    T
 
 s   
6 
60    x 6
j
4
n

  
j
n
 T
  3
  t
3
6
6
 t
s
 x2
n

  
j
Unconditionally stable; but only second-order
Richardson Scheme

Centered-Time, Centered Space (CTCS)
La (T ) 
T
n 1
j
T
T jn 1  T jn 1
2 t
n1
j
2s
n 1
j1


x
2
 4 sT  2 sT
n
j
n
j 1
; s
1
n+1
n
 2 sT
n
j 1

T jn 1  2T jn  T jn 1
 4s
2s
1
j
Second-order in space and time
Unconditionally unstable!
j+1
0
 t
x 2
DuFort-Frankel Scheme

Replace Tjn  (Tjn1  Tjn1 ) / 2 in Richardson scheme
La (T ) 
T jn 1  T jn 1
2 t
(1  2 s )T
n 1
j
 (1  2 s )T
n 1


n1
j
x
 2 sT
n
j 1
2
 2 sT
n
j 1
; s
1+2s
n+1
n

T jn1  (T jn 1  T jn 1 )  T jn 1
2s
2s
j1
1 2s
j
Second-order in space and time
Unconditionally stable!
j+1
0
t
x 2
Finite Difference Methods
 2T
T
0

L(T ) 
2
t
x
T jn  1  T jn
FTC S :
Δt
T jn  T jn 1
Fu llyIm plicit:
Rich ardsonS ch e m e:
Du Fort  Fran k e l:
α
T jn 1
2Δ t
α
T jn 1
α
0
Δx
 2T jn  T jn 1
2
α
Δt
T jn  1  T jn 1
2Δ t
 T jn 1
T jn 1  2T jn  T jn 1
0
Δx
T jn 1  2T jn  T jn 1
Ti n 1
2
0
Δx
 (T jn 1  T jn 1 )  T jn 1
2
Δx
2
0
Finite Difference Methods
FTCS
Fully-Implicit
s
1  2s
s
1  2s
s
2s
12s
 4s
1
Richardson
1
1  2s
s
2s
2s
1  2s
2s
1  2s
1  2s
1  2s
DuFort-Frankel
Truncation Errors

Time spatial
n
t 
1   2 T
n
Ej 
 1   2
2 
6 s   t
FTC S :
Fu llyIm plicit:
Time spatial
t 
1   2 T
n
Ej 
  1   2
2 
6 s   t
Rich ardsonS ch e m e:

  t  x 2   2 T
   
 2
 j  2 12   t
t 
1   2 T
n
Du Fort  Fran k e l: E j 
 2 s   2
2 
6 s   t
5  poin t form u la:
n

  t  x 2   2 T
   
 2

 j  2 12   t
t 
1   2 T
n
Ej 
 0   2
2 
6 s   t
E nj 

t
1  0   T2
2
t
2
n


j
n
n


j
n


 x 2   2 T
   0 
 2
12   t
j 


i

  t  x 2   2 T
   s 
 2
 j  2 12   t


j
n
n

 t
  T
    0  2
  t
j  2
DuFort  Frankel : s 
2
1
12
n


j
n
4.3 Stability





1.
2.
Stability is concerned with the growth, or
decay, of errors introduced at any stage of
the computation
Round-off errors - machine dependent
Intermediate solution for an iterative scheme
For propagation problems, a given method is
stable if the accumulated round-off errors
are negligible
For equilibrium problems:
Direct inversion -- round off errors only
Iterative methods – round-off and iteration errors
Numerical Stability


T  numerical solution without round-off errors
T* numerical solution including round-off errors
  T T
*
n 1
n
n
n

 Exact equation: T j  sT j 1  ( 1  2 s )T j  sT j  1

* n 1
* n
* n
* n
Approximat
ion
:
(
T
)

s
(
T
)

(
1

2
s
)
(
T
)

s
(
T
)j 1

j
j 1
j



n 1
j
 s
n
j 1
 (1  2 s)  s
n
j
Error bound -- assume the worst possible
combinations of individual errors
n
j 1
Numerical Stability

FTCS Scheme, s  ½ for stability
s = 0.6
Round-Off Errors
(T*)nj

n
j
T jn
 jn  Tjn  (T*)nj


Neutral stability – round-off error introduced at
each time step may accumulate (although
cancellation often occurs), but never grow in time
Division of small numbers 1/2 may introduce
significant round-off errors
Stability of FTCS Scheme
If there is no round-off error and j = 0 on
all boundaries, then  jn = 0  stable
n

 In practice,
j   max for all j at step n
(max is machine dependent)

 jn 1  s j-n 1  (1  2 s ) jn  s jn 1   s  1  2 s  s max
0  s  1/2   jn 1   s  (1  2 s )  s max   max

n 1
s

1/2


  s  (1  2 s )  s max  ( 4 s  1) max   max

j
Stability limit: s  1/2
Stability - Matrix Method

Eigenvalue of a tridiagonal matrix
b
a


A   






b c
0


a b c

   
0
a b c

a b 
c
 j  b  2 ac cosθ j
j
θj 

J = 4 (j = 1,2)
 2
J = 6 (j=1,2,3,4)
……
3
,
3
 2
 ,
4 4
 2
 ,
5 5
J = 5 (j = 1,2,3)
In general
J 1
, j  1, 2,  , J  2
3
4
3 4
,
,
5 5

2
( J  2 )

,
, ,
J 1 J 1
J 1
,
4.3.1 Matrix Method:
FTCS Scheme
j=1
2
 1n  1
 n 1
 2
 n 1
 3
 
 n 1
 j
 

 Jn11
 n 1
 J
3
0
4
………..
(B.C.)
 s 1n  (1  2 s ) 2n  s 3n
 s 2n  (1  2 s ) 3n  s 4n
J2
 n 1

J1
J
n
 A
n  0 , 1, 2 , 

 s jn 1  (1  2 s ) jn  s jn 1

 s Jn 2  (1  2 s ) Jn 1  s Jn
0
(B.C.)
Eigenvalues
n
n
A   
n 1
n

 
Matrix Method: FTCS Scheme

The magnitude of all eigenvalues are less than 1
 n 1

n
 n 1
 A  A 
2
1
1
  A    
n
n
s
1  2 s

 s



j
1

2
s
s
0

  j  1  4 s sin2 



s
1  2s s
 2( N  1) 
A  




 1  4 s sin2 θ j /2



0
s 1  2s
s  j  1 ,2 ,...,N-2


1

2
s



2 j 
 j  1  1  1  4 s sin    1

 2
 j
 0  4 s sin2 
 2

  2 for all j


1

1
s
s
2
2
2 sin θ/2
Matrix Method:
Fully-Implicit Scheme

Diagonally-dominant for all s
 sT jn11  (1  2 s)T jn1  sT jn11  T jn
s
1  2 s

 s

1

2
s

s
0




s
1  2s  s
C  







0
 s 1  2s
s 


s
1  2 s 

 n 1
C
 n 1

n


1 n

n
 C   A
l e t C   j ;  A   j
 A  C 
1
1
C

1
j
 j  b  2 ac cosθ j  (1  2 s )  2 s cosθ j  1  4 s sin2 (θ j / 2)  1
1
1
Unconditionally
j 


1
for
all
θ
j
 j 1  4 s sin2 (θ j / 2 )
stable
Pick minimum j (negative sign) so that j = 1/ j is maximum
4.3.2 Matrix Method:
General Two-Level Scheme

Linear combination of FTCS and fully-implicit schemes
T jn 1  T jn
t

   Lxx T jn 1  (1   ) Lxx T jn
Lxx T j 

T j 1  2T j  T j  1
x 2

n+1
 = 0, FTCS scheme
 = 1, Fully-Implicit
1 
j 1
n
j
j+1
 = ½, Crank-Nicolson
General Two-Level Scheme

Matrix Method
 s jn11  (1  2 s ) jn 1  sT jn11  s(1   ) jn1  [1  2 s(1   )] jn  s(1   ) jn 1
 n 1
n
 n 1
n
1
A
 B  
 A B 


 s
 ( 1  2 s )

  s

(
1

2
s

)

s

0




 s
( 1  2 s  )  s
A  







0
 s ( 1  2 s )
 s 


 s
(1  2 s )

s( 1   )
 1  2 s( 1   )
 s( 1   )
1  2 s( 1   )
s( 1   )
0


s( 1   )
1  2 s( 1   ) s( 1   )
B  





0








s( 1   ) 1  2 s( 1   )
s( 1   ) 

s( 1   )
1  2 s(1   )
Matrix Method:
General Two-Level Scheme


jπ
 j  b  2 ac cos  j ;  j 
J 1
Eigenvalues
Pick “+” sign for numerator and “” sign for
denominator for “worst possible” conditions
2



(
1

2
s

)

2
s

cos


1

4
s

sin
( j / 2)
j
 A

2


[
1

2
s
(
1


)]

[
2
s
(
1


)]
cos


1

4
s
(
1


)
sin
( j / 2)

j
 B

Stability criterion
1  4 s(1   ) sin2 ( j / 2 )
B
ST 

1
2
A
1  4 s sin ( j / 2 )
1 1
4 s sin2 ( j / 2 )
1  4 s sin ( j / 2 )
2
1  0
4 s sin2 ( j / 2 )
1  4 s sin ( j / 2 )
2
2
Numerical Stability:
General Two-Level Scheme

Stability criterion
0
4 s sin2 ( j / 2 )
1  4 s sin2 ( j / 2 )
 j
s(1  2  ) sin 
 2
2
 1

 2

 2  4 s sin2 ( j / 2 )  2  8 s sin2 ( j / 2 )
1
1




LHS

0

always satisfie d

2
2


2
sin
( j / 2 )
1
1
   s 


2
2( 1  2  )
2( 1  2  )
 = 0, FTCS scheme: s  1/2
 = 1, Fully-Implicit: unconditionally stable
 = ½, Crank-Nicolson: unconditionally (neutrally)
stable, but oscillatory solution may still occur
4.3.3 Matrix Method:
Derivative Boundary Conditions

Modify A and B matrices to account for Neumann BCs
T
x
x 0
T2  T0
 

2x
To
T2
 s on 1  (1  2 s ) 1n 1  s 2n 1  s(1   ) 0n  [1  2 s(1   )] 1n  s(1   ) 2n
but  0n 1   2n 1 and  0n   2n ( 2x is e xact)
 (1  2 s ) 1n 1  2 s 2n 1  [1  2 s(1   )] 1n  2 s(1   ) 2n
 n 1
A
0 
0 
n
 B  C ; C   

 
0 
Error in textbook
Table 4.2 indicates slight
reduction of stability for
Neumann conditions
Numerical Stability

Generalized two-level scheme
Von Neumann Method
 Fourier




stability method
Most commonly used, easy to apply,
straightforward and dependable
Can only be used for linear, initial value
problem (propagation problem) with
constant coefficients
for nonlinear problems with variable
coefficients, the method may still be
applied “locally” to provide necessary, but
not sufficient stability criterion
may also provide heuristic information
about the influence at the boundary
4.3.4 Von Neumann Method:
FTCS Scheme

Expand the error as a Fourier series
 0j   a m e i j
m
m
For linear problems (superposition implied), it
is sufficient to consider just one term
 G  1, stable
n
n i j
 j  G e ; G : comple x
 G  1, unstable

 jn 1  s jn1  (1  2 s ) jn  s jn 1
G n 1 e i j  sG n e i ( j  1 )  (1  2 s )G n e i j  sG n e i ( j  1)
G  1  s(e  i  e i  2 )  1  4 s sin2 (θ/ 2 )
1
G 1 s
2 sin 2 ( /2 )
4.3.5 Von Neumann Method:
General Two-Level Scheme

For linear problems
 jn1   jn
    Lxx jn 1  (1   ) Lxx jn 
t
 jn  G n e i j ; Lxx jn 
 jn1  2 jn   jn1
x 2
G n e i j  i
i

(
e

2

e
)
2
x
  G n 1 e i j  ( 1   ) G n e i j   i
G n 1 e i j  G n e i j
i
 
(
e

2

e
)

2
t

x


G1

2 (cos  1 )
 i
i


 G  ( 1   )


G

(
1


)
(
e

2

e
)

t
x 2
x 2
G  1  2 s(cos  1 )(  G  1   )  4 s sin 2 (  / 2 )(  G  1   )
1  4 s( 1   ) sin 2 (  / 2 )
4 s sin 2 (  / 2 )
G
 1
2
1  4 s sin (  / 2 )
1  4 s sin 2 (  / 2 )
Von Neumann Method




For more complex equations, it may be
necessary to evaluate the amplification factor
G(s,, ) numerically for a range of s,, and 
values
For three time level schemes, need to solve a
quadratic equation in G
For a system of equations of several variables
(i.e., u, v, w, p), need to solve the eigenvalues
for amplification matrix and require |m| 1
This section deals with numerical instability,
but not the physical instability (transition to
turbulence)
Stability and Consistency




FTCS scheme -- s  1/2
Fully Implicit scheme -- unconditionally stable
Richardson scheme -- unconditionally unstable
DuFort-Frankel scheme -- unconditionally
stable, but inconsistent!!
DuFort  Frankel :
2






x

T

t
x 2   2 T 
 

n
 2    
E j   st 
 
 2 
12  t  j   x 
12   t  j

2
  2T
T

T
t
2


Consistent wi th
 

;  
2
2 
t
x
t 
 x
2
2
n
n
DuFort-Frankel scheme is consistent with a
hyperbolic wave equation!
4.4 Numerical Accuracy






Convergence, Consistency, and stability: establish
limiting behavior for discretization scheme as x, t  0
Asymptotic rate of convergence: x, t  0
Accuracy: deals with practical approximate solution on a
finite grid
Higher-order scheme may not be more accurate than
lower-order ones if the grid is not fine enough
Higher-order scheme has faster rate of convergence, but
the absolute error for a given x (coarse grid) may still
be larger than low-order schemes
Accuracy is problem-dependent, superiority for a simple
model problem may not necessarily imply the same
superiority for more complex problems
Numerical Accuracy

How to determine accuracy when the exact
solution is not available?
1. Grid-refinement study: very expensive to obtain
grid-independent solutions
2. Comparison with experimental data
3. Comparison with analytic solutions / theory
 log E
Higher-order
Lower-order
coarse
fine
log x
Numerical Accuracy

How to improve numerical accuracy?
(1) Different choices of independent variables –
Cartesian, cylindrical, orthogonal
curvilinear, general curvilinear
(2) Different choices of dependent variables –
vorticity/stream-function or primitive
variables?
(3) Adaptive grid – fine grid in high-gradient
regions
(4) Grid refinement together with Richardson
extrapolation
4.4.1 Richardson Extrapolation


Numerical accuracy may be established through
successive grid refinements
For a sufficiently fine grid, the solution error
reduces like the leading term of truncation error

Richardson extrapolation: cancel the leading term
of truncation error for two numerical solutions with
different grid resolutions to achieve higher-order
accuracy

Assume that the leading term dominate the
truncation error, valid only if x is sufficiently fine
More economical than using higher-order scheme
directly

Richardson Extrapolation

Consider two different grids xa and xb
E
10 5
10  4
10  3
m=2
m=4
10  2
10  1
10

Grid A : Ta ,


Grid B : Tb ,
1
10
2
10
3
10
4
10
5
x
Construct a composite solution
E  O( x )
Tc = a Ta + (1 a)Tb
to eliminate the truncation
n
m
E j  O( x b )
error leading term
n
j
m
a
Richardson Extrapolation

Example: FTCS Scheme (for fixed s)
n
n
2
4
4
6





x a
x a 2 1  T
1  T
 La (Ta )  E nj (Ta ) 
( s  ) 4  
( s  ) 6   
2
6  x  j
6
60  x  j



n
n
2
4
4
6







x


x
1  T
1  T
n
2
b
b
L
(
T
)

E
(
T
)

(
s

)

(
s

) 6   
j
b
 a b
 4
2
6
6
60

 x  j
 x  j

La (Tc )  E nj (Tc )  aE nj (Ta )  (1  a ) E nj (Tb )
n
1  T 
 ( s  ) 4  ax a2  (1  a )x b2
2
6  x  j

4

n

1  T 
 ( s  ) 6  ax a4  (1  a )x b4  
6
60  x  j

6
2


Richardson Extrapolation

FTCS Scheme (for fixed s)
La (Tc ) 
n
1  T 
( s  ) 4  ax a2  (1  a )x b2
2
6  x  j

4

n

s 


s 


1
6
1
6

 2 1   6T 
 ( s  ) 6  ax a4  (1  a )x b4  
6
60  x  j
2


x
b
 ax a2  (1  a )x b2  0  a 
x a2  x b2
 x b4
4
4
 ax a  (1  a )x b  0  a 
x a4  x b4


In typical grid refinement study

s 
1

x b   x a  
2
s 

 x b2
1
1
4
1
 a



T

T

Ta
c
b
2
2
6
3
3
3
x a  x b
 x b4
1
1
16
1
 a



T

T

Ta
c
b
4
4
6
15
15
15
 x a  x b
Richardson Extrapolation

FTCS Scheme (for fixed s)
More accurate than
finer-grid solutions