Difference Operators in Numerical Methods: Gershgorin's Theorem & Fourier Analysis, Slides of Mathematical Methods for Numerical Analysis and Optimization

Download Difference Operators in Numerical Methods: Gershgorin's Theorem & Fourier Analysis and more Slides Mathematical Methods for Numerical Analysis and Optimization in PDF only on Docsity!

Numerical Methods for Partial

Differential Equations

Note on AB3 Stability

3 2

z z

z z

ν

z = re^ i^ θ^ ν

Review: Sample Points in Time Space

t n

u u c t x

∂ ∂

∂ ∂

• We chose to sample the ODE at discrete points on the time axis.
• We can also make the same choice for the spatial representation

at each discrete time level.

• We follow this choice with the dilemma of how to reconstruct an

approximation to the derivative of discretely sampled data.

• Assume c is positive

xm (^) − 2 xm (^) − 1 x m xm (^) + 1 xm (^) + 2

x xm (^) − 3 xm (^) + 3

Docsity.com

Review: Choice of Stencil

t n xm (^) − 2 xm (^) − 1 x m xm (^) + 1 xm (^) + 2

x xm (^) − 3 xm (^) + 3

t n xm (^) − 2 xm (^) − 1 x m xm (^) + 1 xm (^) + 2

x xm (^) − 3 xm (^) + 3

t n xm (^) − 2 xm (^) − 1 x m xm (^) + 1 xm (^) + 2

x xm (^) − 3 xm (^) + 3

t n xm (^) − 2 xm (^) − 1 x m xm (^) + 1 xm (^) + 2

x xm (^) − 3 xm (^) + 3

t n xm (^) − 2 xm (^) − 1 x m xm (^) + 1 xm (^) + 2

x xm (^) − 3 xm (^) + 3

And many more combinations Docsity.com

Right Difference Matrix

• Consider a case with 10 points on a periodic interval:
• The semi-discrete scheme for the advection equations is:

m m

du

c u

dt

= δ (^) +

0 0 1 1 2 2

9 9

1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1

u u u u u u

d c dx dx

u u

  (^)  −     (^)  (^) −          (^)  −      (^) −     (^)       (^) −     =  (^) −          (^)  (^) −      (^) −     (^)       (^) −     (^)       −  

Note: I used indexing from zero to match with classical definitions coming up.

Discrete Operator Matrix

• By Gerschgorin’s theorem we conclude that the matrix has all

eigenvalues contained in the circle centered at -c/dx with radius

c/dx

• Thus the eigenvalues are all in the left half plane, with the possible

exception of one or more on the origin  good candidate for time-

stepping. 0

1 1 2 2

9 9

1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1

u u u u u u

d c dx dx

u u

   −      (^) −     (^)       (^) −     (^)  (^) −          (^)  (^) −     =      (^)  −      (^) −     (^)       −     (^)  (^) −           −  

δ (^) +

Docsity.com

Fourier Transform in Matrix Notation

• Introducing the following notation for the Fourier transform matrix and a shift operator:
• We consider the following eigenvalue problem:

ˆ where^ :

if 1

if 1

m M im^ n

M

n nm m nm

m

n

nm m

u u e

u n M u u n M

= − (^) −^ π

= =

 <^ − =   =^ −

∑ F^ F

R

( )

k M

n nk nk k n

k

c c u u u dx

= (^) ∑ R − I =

   (^) or ( )

c u u dx

R − I = λ

 

Right Shift Matrix

• We write down explicitly what the right shift

matrix is:

• Where the box around the addition operator

indicates that this is addition modulo M

• This time the delta indicates the Kronecker Delta:

R nm = δ n + 1, m

cont

( )

( )

1

1 1 2 2 1, 0 0 1 2 2 1

0 2 1 2 2

0 2

1

1

nm M M i^ nj^ i^ mk M M j k j k M i^ nj i^ m^ j M M j i m M M i^ nj^ i^ mj M M j i m M nm

e e M

e e M

e e e M

e

π π

π^ π

π π π

π

δ

δ

−

− − (^) −

= = − (^) − +

=

− (^) −

=

=

=

=

=

∑ ∑

∑

∑

FRF

What we discover after some algebra

(with slight liberties taken with the

notation) is that the Fourier transform

matrix diagonalizes the right shift matrix.

cont

• Finally we return to the and can now find its

eigenspectrum:

• Using
• We find that the matrix is diagonalized by the

Fourier matrix and that

( )

i m

M

nm e nm

π

FRF =

( )

c (^) 1 v v dx

F R − I F =

 

1 for 0,..., 1

i m

c (^) e M m M

dx

π

  = (^)  − (^)  = −  

δ (^) +

Left Difference Matrix

• Consider a case with 10 points on a periodic interval:
• The semi-discrete scheme for the advection equations is:

m m

du

c u

dt

= δ (^) −

0 0 1 1 2 2

9 9

1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1

u u u u u u

d c dx dx

u u

  (^)  −     (^)  (^) −          (^)  −      (^) −     (^)       (^) −     =  (^) −          (^)  (^) −      (^) −     (^)       (^) −     (^)       −  

cont

• Applying Gerschgorin this time we find that the

eigenvalues are all in a disk centered at +c/dx with radius c/dx.

• This is bad news !.
• All the eigenvalues are in the right half plane,

apart from the possiblity of a subset which are at the origin.

• i.e. any non zero mode will be unstable for time

stepping operator (and clearly not physical)

Important Result: Fourier Transform of Left

Shift Matrix

i m

M

nm nm

e

π

FLF =

i.e. the Fourier (similarity) transform of the

left shift matrix is a diagonal matrix

Polynomials of Shift Matrices

• We can express most finite difference operators as polynomials of

shift operators:

• We can apply the Fourier similarity transform to diagonalize this:
• Where is the vector of Fourier coefficents for u.
• Hence each Fourier coefficient satisfies:

0 0

r p j k jk j k

du

a u

dt = =

∑∑ R L

(^1) ( 1 ) ( 1 ) 0 0

r p j k jk j k

dv

a v

dt

− − − = =

F (^) ∑∑ FRF FLF F

v^ ^ = F −^1 u 

( ) ( )

( )

1 1 1 0 0 2 2 2

0 0 0 0

r p j k jk j k r p i^ mj^ i^ mk r p i^ m j^ k m (^) M M M jk m jk m j k j k

dv

a v

dt

dv

a e e v a e v

dt

π π π

− − − = = − −

= = = =

∑∑

∑∑ ∑∑

F FRF FLF F