FINITE DIFFERENCE AND PDE Haryo Tomo
Numerical methods: properties Finite differences
- time-dependent PDEs -> robust, simple concept, easy to parallelize, regular grids, explicit method
Finite elements
- static and time-dependent PDEs
Finite volumes
- time-dependent PDEs
-> implicit approach, matrix inversion, well founded, irregular grids, more complex algorithms, engineering problems
-> robust, simple concept, irregular grids, explicit method
Other numerical methods Particle-based methods
- lattice gas methods - molecular dynamics
Boundary element methods
- problems with boundaries (rupture) - based on analytical solutions
Pseudospectral methods
- orthogonal basis functions, special case of FD - spectral accuracy of space derivatives
- granular problems - fluid flow - earthquake simulations -> very heterogeneous problems, nonlinear problems
- only discretization of planes -> good for problems with special boundary conditions (rupture, cracks, etc)
- wave propagation, ground penetrating radar -> regular grids, explicit method, problems with strongly heterogeneous media
What is a finite difference? Common definitions of the derivative of f(x):
f ( x dx) f ( x) x f lim dx0 dx f ( x) f ( x dx) x f lim dx0 dx
f ( x dx) f ( x dx) x f lim dx0 2dx These are all correct definitions in the limit dx->0. But we want dx to remain FINITE
What is a finite difference? The equivalent approximations of the derivatives are:
f ( x dx) f ( x) x f dx
forward difference
f ( x) f ( x dx) x f dx
backward difference
f ( x dx) f ( x dx) x f 2dx
centered difference
The
big question:
How good are the FD approximations?
This leads us to Taylor series....
Taylor Series Taylor series are expansions of a function f(x) for some finite distance dx to f(x+dx) 2 3 4 dx dx dx f ( x dx) f ( x) dxf ' ( x) f '' ( x) f ''' ( x) f '''' ( x) ... 2! 3! 4!
What happens, if we use this expression for
f ( x dx) f ( x) x f dx
?
Taylor Series ... that leads to : 2 3 f ( x dx) f ( x) 1 dx dx ' '' ''' f ( x) f ( x) ... dxf ( x) dx dx 2! 3! f ' ( x) O(dx)
The error of the first derivative using the forward formulation is of order dx. Is this the case for other formulations of the derivative? Let’s check!
Taylor Series ... with the centered formulation we get: 3 f ( x dx / 2) f ( x dx / 2) 1 dx ' ''' f ( x) ... dxf ( x) dx dx 3! f ' ( x) O(dx 2 )
The error of the first derivative using the centered approximation is of order dx2. This is an important results: it DOES matter which formulation we use. The centered scheme is more accurate!
Higher order operators *a |
(2dx) 2 (2dx)3 f ( x 2dx) f (2dx) f ' f ' ' f ''' 2! 3!
*b |
(dx) 2 (dx)3 f ( x dx) f (dx) f ' f ' ' f ''' 2! 3!
*c |
(dx) 2 (dx)3 f ( x dx) f (dx) f ' f ' ' f ''' 2! 3!
*d |
(2dx) 2 (2dx)3 f ( x 2dx) f (2dx) f ' f ' ' f ''' 2! 3! ... again we are looking for the coefficients a,b,c,d with which the function values at x±(2)dx have to be multiplied in order to obtain the interpolated value or the first (or second) derivative! ... Let us add up all these equations like in the previous case ...
Problems: Stability c 2 dt 2 p( x dx) 2 p( x) p( x dx) p (t dt ) 2 dx 2 p (t ) p (t dt ) sdt 2 Stability: Careful analysis using harmonic functions shows that a stable numerical calculation is subject to special conditions (conditional stability). This holds for many numerical problems. (Derivation on the board).
dt c 1 dx
Problems: Dispersion c 2 dt 2 p( x dx) 2 p( x) p( x dx) p (t dt ) 2 dx 2 p (t ) p (t dt ) sdt 2
True velocity
Dispersion: The numerical approximation has artificial dispersion, in other words, the wave speed becomes frequency dependent (Derivation in the board). You have to find a frequency bandwidth where this effect is small. The solution is to use a sufficient number of grid points per wavelength.
Finite Differences - Summary
Conceptually the most simple of the numerical methods and can be learned quite quickly Depending on the physical problem FD methods are conditionally stable (relation between time and space increment) FD methods have difficulties concerning the accurate implementation of boundary conditions (e.g. free surfaces, absorbing boundaries) FD methods are usually explicit and therefore very easy to implement and efficient on parallel computers FD methods work best on regular, rectangular grids
Partial Differential Equations
by Lale Yurttas, Texas A&M University
Part 8
14
PERSAMAAN DIFERENSIAL PARSIAL • Persamaan Umum
2 2 2 a 2 b c 2 d e f g 0 x xy y x y • Menyatakan bagaimana variabel tak bebas Ø berubah
terhadap variabel bebas x,y. Disini a,b,c,d,e,f, dan g mungkin merupakan fungsi dari Ø
Jenis2 PDP • Ditentukan oleh harga b2-4ac < 0, eliptic = 0, parabolic > 0, hyperbolic • Adveksi… • Difusi…. • Gelombang…
JENISJENIS PDP
by Lale Yurttas, Texas A&M University
18
The Laplacian Difference Equations/ 2T 2T 2 0 2 Laplace Equation x y 2T Ti 1, j 2Ti , j Ti 1, j O[(x)2] x 2 x 2 2T Ti , j 1 2Ti , j Ti , j 1 O[(y)2] 2 2 y y Ti 1, j 2Ti , j Ti 1, j Ti , j 1 2Ti , j Ti , j 1 0 2 2 x y x y Laplacian difference Ti 1, j Ti 1, j Ti , j 1 Ti , j 1 4Ti , j 0
equation. by Lale Yurttas, Texas A&M University
Holds for all interior points Chapter 29
by Lale Yurttas, Texas A&M University
Chapter 29
20
• In addition, boundary conditions along the edges must be
specified to obtain a unique solution. • The simplest case is where the temperature at the boundary is set at a fixed value, Dirichlet boundary condition. • A balance for node (1,1) is:
T21 T01 T12 T10 4T11 0 T01 75 T10 0 4T11 T12 T21 0 • Similar equations can be developed for other interior points to
result a set of simultaneous equations. by Lale Yurttas, Texas A&M University
Chapter 29
21
• The result is a set of nine simultaneous equations with nine unknowns: 4T11 T21 T11 4T21 T13 T21
T12 T22
4T31
T11 T21
75 0 T32
4T12
T22
T12
4T22 T22
T31
50 T13
T32 4T32
T12 T22 T32
75 T23 T33
0 50
4T13 T23
175
T13
100
4T23 T33 T23
4T33 150
h
Diffusion
x h 2h D 2 t x
Diffusion Equation
h ud D x
x
h ud ud x x
ud
x
ud h x xt x
x
ud h x xt x t 0 ud h t x
h h D 2 t x
h ud D x
2
Diffusion Equation
Numerical Solution of Diffusion Eq. h
h(i 1) h(i)
h(i 1) x h h(i ) h(i 1) x x 1
h(i 1) h(i) h x x 2
2 h 1 h h 2 x x x 2 x 1
Numerical Calculation of Diffusion Equation
J
t j+
h(i 1, j 1)
j
h(i, j )
1
2 ・・ j-1
t
x 1
2
・・・・・ i-1
i
i+1
・・・・
N
x
h h(i, j 1) h(i, j ) t t h D h(i 1, j ) h(i, j ) h(i, j ) h(i 1, j ) D 2 x x x x h(i 1, j ) 2h(i, j ) h(i 1, j ) D 2 x Unknown 2
h(i, j 1) h(i, j )
Known
Dt 2 {h(i 1, j ) 2h(i, j ) h(i 1, j )} x
Numerical Solution of Diffusion Equation
t
j
j+
J
Boundary Condition (Given)
1
2 ・・ j-1
t
x 1
2
・・・・・ i-1
i
i+1
・・・・
N
x
Initial Condition (Given)
Dt h(i, j 1) h(i, j ) 2 {h(i 1, j ) 2h(i, j ) h(i 1, j )} x do i 2,N-1 t hnew (i ) hold (i ) hold (i 1) 2hold (i ) hold (i 1) 2 D x end do do i 2,N-1 hold (i ) hnew (i ) end do
Contoh d2y dy 2 3 2 y 4 x , 2 dx dx
y(1) 1,
y(2) 6
• Cari solusinya dengan step size 0,2 • Jumlah titik solusi n=((2-1)/0,2)-1=4 • Kita dapatkan 4 persamaan, satu untuk tiap titik yang
dicari.
Penyelesaian dengan Beda Hingga Persamaannya:
yi 1 2 yi yi 1 yi 1 yi 1 2 3 2 yi 4 xi 2 Buat persamaan untuk semua h 2htitik, mulai dari i=1, hingga i=4
x0=a
x1
x2
x3
x4
x5=b
Domain Solusi y5=6
y0=1 x0=1
x5=2
Penyelesaian dg Finite Diff. • Dengan h = 0,2
32,5 yi 1 48 yi 17,5 yi 1 4 xi • Buat persamaan untuk semua titik, mulai dari i=1,
hingga i=4 • Masukkan nilai-nilai x1=1,2 hingga x4 =1,8 dan kondisi batas y0 = 1 dan y5 = 6
2
Kondisi batas • Dirichlet atau fixed boundary,
misal C(0) = Co • Neuman atau natural boundary, misal dC/dx = 0 • Robin/ Cauchy boundary condition, misal dC/dx + C = 0 • Penerapan dalam finite difference dengan menambahkan
imaginary node
penyelesaian persamaan parabolik dengan skema eksplisit
SEHINGGA… persamaan diferensial,
ditulis dalam bentuk metode beda hingga,
laplace ….equation
penyelesaian persamaan eliptik
Stabilitas skema eksplisit
k
k
penyelesaian persamaan parabolik dengan skema implisit