scientific computing i
play

Scientific Computing I Module 7: Solving the 1D Heat Equation - PowerPoint PPT Presentation

Lehrstuhl Informatik V Scientific Computing I Module 7: Solving the 1D Heat Equation Miriam Mehl based on Slides by Michael Bader (Winter 09/10) Winter 2011/2012 Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing


  1. Lehrstuhl Informatik V Scientific Computing I Module 7: Solving the 1D Heat Equation Miriam Mehl based on Slides by Michael Bader (Winter 09/10) Winter 2011/2012 Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 1

  2. Lehrstuhl Informatik V Part I: Analytical Solution The Heat Equation in 1D Analytical Solutions Computing Analytical Solutions – Steps Separation of Variables Ensuring the Initial Conditions – Fourier’s Method Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 2

  3. Lehrstuhl Informatik V The Heat Equation in 1D • remember the heat equation without heat sources: T t = κ ∆ T • we examine the 1D case, and set κ = 1 to get: u t = u xx x ∈ ( 0 , 1 ) , t > 0 for • using the following initial and boundary conditions: u ( x , 0 ) = g ( x ) , x ∈ ( 0 , 1 ) u ( 0 , t ) = u ( 1 , t ) = 0 , t > 0 x ? t Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 3

  4. Lehrstuhl Informatik V Computing Analytical Solutions – Steps Steps: 1. try to find some solution of the PDE 2. try to satisfy boundary conditions 3. try to satisfy initial conditions Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 4

  5. Lehrstuhl Informatik V Separation of Variables • assumption: u ( x , t ) = X ( x ) · T ( t ) • insert this assumption into the heat equation ∂ 2 ∂ ∂ t ( X ( x ) · T ( t )) = ∂ x 2 ( X ( x ) · T ( t )) X ( x ) · T t ( t ) = T ( t ) · X xx ( x ) or • divide by X ( x ) T ( t ) , and get: T t ( t ) T ( t ) = X xx ( x ) X ( x ) Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 5

  6. Lehrstuhl Informatik V Transforming the PDE into two ODEs • separation of variables leads to: T t ( t ) T ( t ) = X xx ( x ) X ( x ) = − λ. • thus, we obtain two ODEs: X xx ( x ) + λ X ( x ) = 0 (1) T t ( t ) + λ T ( t ) = 0 . (2) • solve X ( x ) -part: � √ � √ � � X ( x ) = sin λ x or X ( x ) = cos λ x , • solve T ( t ) -part: T ( t ) = e − λ t . Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 6

  7. Lehrstuhl Informatik V Ensuring the Boundary Conditions • We have computed general solutions � √ � √ � e − λ t or u ( x , t ) = cos � e − λ t . u ( x , t ) = sin λ x λ x 1 0.5 0 The boundary conditions of −0.5 our example are −1 1 u ( 0 , t ) = u ( 1 , t ) = 0 for all t . 0.5 0.1 0.08 0.06 0.04 x 0 0.02 0 t • Thus, we can reduce possible solutions to u k ( x , t ) = sin ( k π x ) e − ( k π ) 2 t , k = 1 , 2 , 3 . . . Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 7

  8. Lehrstuhl Informatik V Ensuring the Initial Conditions – Fourier’s Method • The functions u k ( x , t ) = sin ( k π x ) e − ( k π ) 2 t , k = 1 , 2 , . . . , solve the 1D heat equation PDE with correct boundary values and for the initial conditions: u k ( x , 0 ) = sin ( k π x ) , x ∈ ( 0 , 1 ) . 1 1 1 0.8 0.5 0.5 0.6 0 0 0.4 −0.5 −0.5 0.2 0 −1 −1 1 1 1 0.5 0.5 0.5 0.1 0.1 0.1 0.08 0.08 0.08 0.06 0.06 0.06 0.04 0.04 0.04 x 0.02 x 0.02 x 0.02 0 0 0 0 0 0 t t t • To ensure the initial conditions, we use the Fourier sine series: ∞ � g ( x ) = c k sin ( k π x ) k = 1 Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 8

  9. Lehrstuhl Informatik V Limits of Analytical Solutions • For the 1D heat equation, we obtain an analytical solution ∞ c k e − ( k π ) 2 t sin ( k π x ) . � u ( x , t ) = k = 1 • Solutions for the heat equation with heat sources? Separation of variables yields X ( x ) · T t ( t ) = T ( t ) · X xx ( x ) + f ( x ) . ⇒ No analytical solution computable for arbitrary f . • need for numerical methods!!! Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 9

  10. Lehrstuhl Informatik V Part II: Numerical Solution Numerical Solution 1 – An Explicit Scheme Discretisation Accuracy of the Explicit Scheme Stability of the Explicit Scheme Numerical Solution 2 – An Implicit Scheme Implicit Time-Stepping Stability of the Implicit Scheme Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 10

  11. Lehrstuhl Informatik V Numerical Solution 1 – An Explicit Scheme Discretisation similar to ODEs: • compute approximations v ( m ) ≈ u ( x j , t m ) j at grid points x j and time points t m : x N x ? x 2 x 1 x j := j · h t m := m · τ, t 2 t M t 1 t v ( m ) j − 1 − 2 v ( m ) + v ( m ) j j + 1 approximate u xx at ( x j , t m ) by , h 2 v ( m + 1 ) − v ( m ) j j u t at ( x j , t m ) by . τ Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 11

  12. Lehrstuhl Informatik V An Explicit Scheme • This approximation results in the equations v ( m + 1 ) − v ( m ) v ( m ) j − 1 − 2 v ( m ) + v ( m ) j j j j + 1 = . (3) τ h 2 for j = 1 , . . . , N − 1 and m ≥ 0. • Add initial and boundary conditions: v ( m ) v ( m ) = = 0 , for all m ≥ 0 , n 0 v ( 0 ) = g ( x j ) , for j = 1 , . . . , n − 1 j • and obtain an explicit scheme: + τ � � v ( m + 1 ) = v ( m ) v ( m ) j − 1 − 2 v ( m ) + v ( m ) . j j j j + 1 h 2 • We can, step by step, compute v ( m ) starting with v ( 0 ) = g ( x j ) . j j Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 12

  13. Lehrstuhl Informatik V Accuracy of the Explicit Scheme Observations: • first order accurate in τ : u t ( x , t ) = u ( x , t + τ ) − u ( x , t ) + O ( τ ) τ • second order accurate in h : u xx ( x , t ) = u ( x − h , t ) − 2 u ( x , t ) + u ( x + h , t ) + O ( h 2 ) h 2 Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 13

  14. Lehrstuhl Informatik V Stability of the Explicit Scheme • exact solution: u k ( x , t ) = e − ( k π ) 2 t sin ( k π x ) • explicit scheme: v ( m + 1 ) − v ( m ) v ( m ) j − 1 − 2 v ( m ) + v ( m ) j j j j + 1 = τ h 2 Is there a stability condition for τ as observed in the ODE case? • assumption: there are solutions of the form = ( a k ) m sin ( π kx j ) , v ( m ) where x j := jh . (4) j • compare with the exact solution ⇒ ( a k ) m should decrease ⇒ | a k | < 1 necessary Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 14

  15. Lehrstuhl Informatik V Stability of the Explicit Scheme (2) • Inserting v ( m ) = ( a k ) m sin ( π kx j ) into j v ( m + 1 ) − v ( m ) v ( m ) j − 1 − 2 v ( m ) + v ( m ) j j j j + 1 = h 2 τ leads to � π kh � a k = 1 + τ h 2 ( 2 cos ( π kh ) − 2 ) = 1 − 4 τ h 2 sin 2 . 2 • for stability: | a k | < 1 for all k if τ < h 2 4 τ h 2 < 2 or 2 . Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 15

  16. Lehrstuhl Informatik V Stability of the Explicit Scheme (3) • solution for general initial conditions: n − 1 v ( 0 ) � = f j = c k sin ( π kx j ) , j k = 1 n − 1 c k ( a k ) m sin ( π kx j ) v ( m ) � := j k = 1 • stability for τ < h 2 2 because | a k | < 1 for all k Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 16

  17. Lehrstuhl Informatik V Numerical Solution 2 – An Implicit Scheme • apply implicit Euler: v ( m + 1 ) − v ( m ) v ( m + 1 ) − 2 v ( m + 1 ) + v ( m + 1 ) j j j − 1 j j + 1 = h 2 τ for j = 1 , . . . , n − 1, and m ≥ 0. • boundary conditions: v ( m ) = v ( m ) = 0 , m ≥ 0 , for all n 0 • initial conditions: v ( 0 ) = g ( x j ) , for j = 1 , . . . , n − 1 . j Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 17

  18. Lehrstuhl Informatik V An Implicit Scheme • implicit equations for v ( m + 1 ) : j − τ � � v ( m + 1 ) v ( m + 1 ) − 2 v ( m + 1 ) + v ( m + 1 ) = v ( m ) . j j − 1 j j + 1 j h 2 • With the ratio r := τ/ h 2 , we can write it as − rv ( m + 1 ) + ( 1 + 2 r ) v ( m + 1 ) − rv ( m + 1 ) = v ( m ) j − 1 j j + 1 j for j = 1 , . . . , n − 1, and m ≥ 0. • solution: v ( m ) = ( I + rA ) − 1 v ( m − 1 ) , where A = tridiag ( − 1 , 2 , − 1 ) . • Use a solver for systems of linear equations to obtain v ( m + 1 ) in j every step. Miriam Mehl based on Slides by Michael Bader (Winter 09/10): Scientific Computing I Module 7: Solving the 1D Heat Equation, Winter 2011/2012 18

Recommend


More recommend