Numerical Solution of Partial Differential Equations Praveen. C - PowerPoint PPT Presentation

Numerical Solution of Partial Differential Equations Praveen. C praveen@math.tifrbng.res.in Tata Institute of Fundamental Research Center for Applicable Mathematics Bangalore 560065 http://math.tifrbng.res.in 31 January 2009 Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 1 / 40

A continuous function u ( x ) = sin(2 πx ) , x ∈ [0 , 1] u ( x ) x Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 2 / 40

Space of continuous functions • C ([0 , 1]) = Space of continuous functions • C ([0 , 1]) has infinitely many elements in it = ⇒ We say it is infinite dimensional • Example: u ( x ) = sin(2 πx ) is an element of C ([0 , 1]) • A computer has finite memory = ⇒ It cannot represent infinite dimensional objects Hence we need to approximate an infinite dimensional object using a finite dimensional space Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 3 / 40

Discrete approximation u ( x ) = sin(2 πx ) , x ∈ [0 , 1] u ( x ) x 1 Sample at discrete set of N points, with spacing h = N − 1 x i = ( i − 1) h, 1 ≤ i ≤ N u i = u ( x i ) Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 4 / 40

Discrete approximation u ( x ) = sin(2 πx ) , x ∈ [0 , 1] u ( x ) x U = [ u 1 , u 2 , . . . , u N ] ∈ R N U provides a finite dimensional approximation to the continuous function Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 5 / 40

Discrete approximation u ( x ) = sin(2 πx ) , x ∈ [0 , 1] u ( x ) x Piecewise linear approximation u h ( x ) N → ∞ = ⇒ h → 0 and u h ( x ) → u ( x ) Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 6 / 40

Ordinary Differential Equation • u = u ( x ) : function of single variable x • Second order ODE − d 2 u d x 2 = sin( x ) , x ∈ (0 , 2 π ) with boundary condition u (0) = u (2 π ) = 0 • Exact solution u ( x ) = sin( x ) This can be called a symbolic solution. Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 7 / 40

Ordinary Differential Equation • If we have a general ODE − d 2 u d x 2 = f ( x ) , x ∈ (0 , 2 π ) there may not exist an explicit, analytical solution • Example: − d 2 u d x 2 = exp( − x 2 ) , x ∈ (0 , 2 π ) = ⇒ Need for numerical solution • Numerical techniques for ODE/PDE 1 Finite Difference Method 2 Finite Volume Method 3 Finite Element Method Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 8 / 40

Finite Difference Method • Consider a general ODE − d 2 u d x 2 = f ( x ) , x ∈ ( a, b ) with boundary conditions u ( a ) = u a , u ( b ) = u b • Instead of finding a function that solves the ODE, find a discrete approximation to the solution • Divide computational domain ( a, b ) into N + 1 intervals each of size h = b − a N + 1 • Computational grid x i = ih, 0 ≤ i ≤ N + 1 Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 9 / 40

Finite Difference Method • Discrete solution x 0 , x 1 , x 2 , . . . , x N , x N +1 u 0 , u 1 , u 2 , . . . , u N , u N +1 u i ≈ u ( x i ) , i = 0 , . . . , N + 1 • From boundary condition u 0 = u a , u N +1 = u b Hence, we need to find U = [ u 1 , u 2 , . . . , u N ] ∈ R N This is a finite dimensional problem Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 10 / 40

Finite Difference Method • We must approximate the ODE at the interior grid points x 1 , x 2 , . . . , x N • Differentiation is limit of finite difference d u ( x + h ) − u ( x ) d xu ( x ) = lim h h → 0 • Finite difference approximation d xu ( x i ) ≈ u i +1 − u i d h Finite Difference Method Approximate differential operators by finite difference operators Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 11 / 40

FDM: First derivative • From Taylor’s formula u i +1 = u i + hu ′ ( x i ) + h 2 2 u ′′ ( ξ ) • Forward difference approximation δ + u i := u i +1 − u i = u ′ ( x i ) + h 2 u ′′ ( ξ ) h � �� � O ( h ) Forward difference is first order accurate • Backward difference approximation δ − u i := u i − u i − 1 = u ′ ( x i ) + O ( h ) h • Central difference approximation δ 0 u i := u i +1 − u i − 1 = u ′ ( x i ) + O ( h 2 ) 2 h Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 12 / 40

FDM: Second derivative • Second order derivative u ′ ( x + h/ 2) − u ′ ( x − h/ 2) u ′′ ( x ) = lim h h → 0 • Finite difference approximation u i +1 − u i − u i − u i − 1 = u i − 1 − 2 u i + u i +1 h h ∆ u i = h 2 h • Second order accurate ∆ u i = u ′′ ( x i ) + O ( h 2 ) Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 13 / 40

Back to the ODE • Replace ODE − u ′′ ( x ) = f ( x ) , x ∈ ( a, b ) with finite difference approximation − ∆ u i = f i , i = 1 , 2 , . . . , N • At i = 1 h 2 u 1 − 1 2 h 2 u 2 = f 1 + 1 h 2 u a • For i = 2 , . . . , N − 1 − 1 h 2 u i − 1 + 2 h 2 u i − 1 h 2 u i +1 = f i • At i = N − 1 h 2 u N − 1 + 2 h 2 u N = f N + 1 h 2 u b Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 14 / 40

FDM for ODE For N = 10 , putting all equations together       f 1 + u a u 1 2 − 1 0 0 0 0 0 0 0 0 h 2 u 2 f 2 − 1 2 − 1 0 0 0 0 0 0 0       u 3 f 3 0 − 1 2 − 1 0 0 0 0 0 0             0 0 − 1 2 − 1 0 0 0 0 0 u 4 f 4       1  0 0 0 − 1 2 − 1 0 0 0 0   u 5   f 5  =       h 2 0 0 0 0 − 1 2 − 1 0 0 0 u 6 f 6             0 0 0 0 0 − 1 2 − 1 0 0 u 7 f 7             0 0 0 0 0 0 − 1 2 − 1 0 u 8 f 8       0 0 0 0 0 0 0 − 1 2 − 1 u 9 f 9 f 10 + u b 0 0 0 0 0 0 0 0 − 1 2 u 10 h 2 or A h U h = b h We have N equations for the N unknowns: [ u 1 , u 2 , . . . , u N ] Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 15 / 40

FDM for ODE • Matrix A h is invertible = ⇒ Solution to discrete problem exists • Efficient solution using Thomas Tri-diagonal algorithm 1 Numerical Exact 0.8 0.6 0.4 0.2 u(x) 0 −0.2 −0.4 −0.6 −0.8 −1 0 1 2 3 4 5 6 x Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 16 / 40

Partial Differential Equations • Problems involving more than one independent variable u ( x, t ) : x is space, t is time u ( x, y ) : x, y denotes two space coordinates u ( x, y, t ) : x, y denotes two space coordinates, t is time = ⇒ Leads to Partial Differential Equation • One space and one time: u ( x, t ) ◮ Hyperbolic equation ∂ 2 u ∂t 2 = c 2 ∂ 2 u ∂x 2 ◮ Elliptic equation ∂t = µ∂ 2 u ∂u ∂x 2 ◮ Parabolic equation ∂x = µ∂ 2 u ∂u ∂t + c∂u ∂x 2 Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 17 / 40

Simplest hyperbolic PDE • Linear, scalar, convection (advection) equation for u ( x, t ) ∂u ∂t + c∂u ∂x = 0 , x ∈ R with initial condition u ( x, 0) = u 0 ( x ) • Exact solution u ( x, t ) = u 0 ( x − ct ) u t = 0 t = ∆ t c ∆ t c > 0 x Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 18 / 40

Hyperbolic PDE Wave A phenomenon in which some recognizable feature propagates with a recognizable speed Hyperbolic PDE A PDE which has wave-like solutions • Waves propagate in specific directions: • Linear, convection equation ◮ c > 0 = ⇒ wave moves to the right ◮ c < 0 = ⇒ wave moves to the left ◮ c is the speed at which waves propagate ◮ Finite speed of propagation ◮ Preserves shape of initial condition Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 19 / 40

Hyperbolic PDE • Scalar, convection equation � ∂ � ∂t + c ∂ u = 0 ∂x contains one wave • Second order wave equation ∂ 2 u ∂t 2 = c 2 ∂ 2 u ∂x 2 ◮ can be factored � ∂ � � ∂ � ∂t + c ∂ ∂t − c ∂ u = 0 ∂x ∂x ◮ contains two waves, with speed + c and − c ◮ In fact, general solution is u ( x, t ) = f ( x − ct ) + g ( x + ct ) Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 20 / 40

Elliptic PDE • Example: Heat equation ∂t = µ∂ 2 u ∂u ∂x 2 , x ∈ R with initial condition u ( x, 0) = u 0 ( x ) t t+ ∆ t • No waves; initial condition is damped or dissipated Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 21 / 40

Parabolic PDE • Convection-diffusion equation ∂x = µ∂ 2 u ∂u ∂t + c∂u ∂x 2 contains convection and diffusion t+ ∆ t t • Damped wave-like solutions Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 22 / 40

FDM for u t + cu x = 0 • Given u ( x, 0) = u 0 ( x ) , find solution for t > 0 : Initial Value Problem • Space-time grid t n ∆ t x ∆ x i • Numerical solution u n i u n i ≈ u ( x i , t n ) Numerical solution computed only at grid points Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 23 / 40

FDM for u t + cu x = 0 • Forward difference in time ∂tu ( x i , t n ) ≈ u n +1 − u n ∂ i i ∆ t ∂ • Three choices for ∂x 1 Backward difference ∂xu ( x i , t n ) ≈ u n i − u n ∂ i − 1 ∆ x 2 Forward difference ∂xu ( x i , t n ) ≈ u n i +1 − u n ∂ i ∆ x 3 Central difference ∂xu ( x i , t n ) ≈ u n i +1 − u n ∂ i − 1 2∆ x Praveen. C (TIFR-CAM) Numerical PDE Jan 31, 2009 24 / 40