# Lecture 34 - Ordinary Differential Equations - BVP by tdo11445

VIEWS: 70 PAGES: 48

• pg 1
```									Lecture 34 - Ordinary Differential
Equations - BVP
CVEN 302
November 28, 2001
Systems of Ordinary Differential
Equations - BVP
• Shooting Method for Nonlinear BVP
• Finite Difference Method
• Partial Differential Equations
Shooting Method for Nonlinear
ODE-BVPs
• Nonlinear        y   f ( x , y , y  ) , a  x  b
 y( a )  y , y( b )  y
ODE
             a                b

• Consider u  f ( x , u , u ) with guessed slope t
u(a)  y , u( a )  t
        a
• Use the difference between u(b) and yb to adjust u’(a)
• m(t) = u(b, t) - yb is a function of the guessed value t
• Use secant method or Newton method to find the correct t
value with m(t) = 0
Nonlinear Shooting Based on Secant
Method

• Nonlinear  y   f ( x , y , y  ) , a  x  b
 y( a )  y , h( y( b ), y ( b ))  0
ODE                   a

Step1 : use u(a)  y a , u (a)  t(1)  Error  m(1)
Step 2 : use u(a)  y a , u (a)  t(2)  Error  m(2)
Step 3 : use secantmethodto obtain a new estimate
t(i  1)  t(i  2)
t(i)  t(i  1)                       m(i  1)
m(i  1)  m(i  2)
Step 4 : iterate until t(i) - t(i - 1)  tol
MATLAB Example in Nonlinear Shooting
Method
• Nonlinear shooting with secant method
 y   2 yy  , 0  x  1

 y( 0 )  1 , h( y( 1 ), y ( 1 ))  y( 1 )  y ( 1 )  0 .25  0
exact solution y  1 /( x  1 )
• Convert to two first-order ODE-IVPs
let z 1  y , z 2  y 

z1  z 2 ,              z1 ( 0 )  1

 z   2 z 1 z 2 ,
2                    z2 ( 0 )  t
• Update t using the secant method
Nonlinear Shooting - Secant Method

y(x)

y’(x)
Nonlinear Shooting
Based on Newton’s Method

• Nonlinear  y   f ( x , y , y  ) , a  x  b

ODE       y( a )  y a , y( b )  y b
u  f(x,u, u)                           u(a)  y a , u(a)  t k
                                          
v   vf u (x,u, u)  v f u (x,u, u)  v(a)  0, v (a)  1
Check for convergence of m(t)
if m  tol , stop

m( t k )  u(b, t k )  y b                                  m
Otherwise , t k  1  t k  v ( b , t )
                                     k
Nonlinear Shooting
Based on Newton’s Method
• Nonlinear u  f ( x , u , u ) , u( a )  y a , u( a )  t

ODE-IVP m( t )  u( b , t )  y b
0
• Chain u f ( x , u , u ) f x f u f u
                                         
Rule t          t              x t       u t u t
x and t are independent
u              u          u
let v         , then v       , v  
t              t            t
u f u f u
                      v   f u v  f u v 
t       u t u t
u( a , t )  t  v' (a)  u/t  1
Nonlinear Shooting with Newton’s Method
• Solve ODE-IVP
          ( y  )2                ( u ) 2
 y             ,      u            ,
             y                      u
 y( 0 )  1 , y( 1 )  2  u( 0 )  1 , u( 0 )  1
                          
• Construct the auxiliary equations
                      f ( u  )2
 f u ( x , u , u )                     ( u ) 2     2 u
                      u   u  2
v            v       v
                                           u 2
u
 f ( x , u , u )  f   2 u    v ( 0 )  0 , v ( 0 )  1
 u                                
                      u     u
Nonlinear Shooting with Newton’s Method

• Calculate m(t) -- deviation from the exact BC
m( t )  u(b, t )  y b

 m t ( t )  u t ( b, t )
• Update t by Newton’s method
m(t i 1 )                    m
t i  t i 1                  ; v  m t (t)      u t ( b, t )
m t (t i 1 )                 t
Finite-Difference Methods

• Divide the interval of interest into subintervals
ba
x0  a , x n  b , h       x i 1  x i
h
• Replace the derivatives by appropriate finite-difference
approximations in Chapter 11
• Solve the system of algebraic equations by methods in
Chapters 3 and 4
• For nonlinear ODEs, methods in Chapter 5 may be used
Finite-Difference Method
• General Two-Point BVPs
 y   p( x ) y   q( x ) y  r ( x ), a  x  b

 y(a)   , y(b)  
• Replace the derivatives by appropriate finite-difference
approximations
y i  1  2 y i  y i 1                y i  1  y i 1
y ( x i )                           , y ( x i ) 
h2                                   2h
h            h         h          h
xi-1       xi        xi+1
Finite-Difference Method
• Central difference approximations
y i 1  y i 1
y ( x i ) 
2h
y i 1  2 y i  y i 1
y ( x i ) 
h2
• Tridiagonal system
y i  1  2 y i  y i 1      y i 1  y i 1
2
 pi                  q i y i  ri
h                         2h


h 

1  pi  y i 1  2  h q i y i   1  pi  y i  1  h 2 ri
2
             h 
    2                                         2 
Finite-Difference Method
• Central Difference ==> Tridiagonal system
 ( 2  h 2 q1 )         0                    0                  0            y1 
                                                                             
 1  (h/2)p 2  ( 2  h 2 q 2 ) 1  (h/2)p 2                     0            y2  
       0           1  (h/2)p 3  ( 2  h q 3 ) 
2
0

 y3   
                                                                                    
                                                                          
                                                       ( 2  h 2 q n 1 )  y n 1 
       0                 0                    0                                     
 h 2 r1  ( 1  hp1 /2 ) y 0   h 2 r1  ( 1  hp1 /2 )  
                                                              
            h 2 r2                           h 2 r2           
                                                              
                2
h r3                              2
h r3             
                                                            
                                                              
h rn 1  ( 1  hpn -1 /2 ) y n  h rn 1  ( 1  hpn -1 /2 ) 

2
 
2

Finite-Difference Method for
Nonlinear BVPs
• Nonlinear ODE-BVPs
 y   f ( x , y , y  ) , a  x  b

 y(a)   , y(b)  
• Evaluate fi by appropriate finite-difference
approximations
y i 1  2 y i  y i 1
2
 fi  0
h
h           h        h          h
xi-1       xi       xi+1
Finite-Difference Method for
Nonlinear BVPs
• SOR method
1
yi             ( y i  1  2 y i  y i  1  h 2 f i )
2( 1   )
• Iterative solution
• Convergence criterion

  h Q / 2
2        0  Q  f y ( x , y , y' )  Q 

            with 
h  2 / P

 f y ( x , y , y  )  P

Example 14.12 - MATLAB

• Note error in Text
( y )2
y   f ( x , y , y  )             , y(0)  1 , y(1)  2
y
1
yi                  ( y i  1  2 y i  y i  1  h 2 f i )
2( 1   )

f       
( y i 1  y i 1 ) / 2 h 2   
( y i 1  y i 1 )2
i
yi                          4h2 yi

fi : negative sign
Chapter 15

Partial Differential Equations
Classification of PDEs
• General form of linear second-order PDEs with
two independent variables

au xx  bu xy  cu yy  du x  eu y  fu  g  0

• linear PDEs: a, b, c,….,g = f(x,y) only

b 2  4 ac  0 , Hyperbolic (2 real roots)
 2
b  4 ac  0 , Parabolic (1 double root)
b 2  4 ac  0 , Elliptic  (2 complex roots)

Heat Equation: Parabolic PDE
• Heat transfer in a one-dimensional rod
x=0                                     x=a
g1(t)                                         g2(t)

u       2u
c       ; 0  x  a, 0  t  T
t     x  2

u(x,0)  f(x),        0 xa
u(0,t)  g 1 ( t )
                       , 0 t T
u( a , t )  g 2 ( t )
Discretize the solution domain in space
and time with h = x and k = t
t
10

9

8

7

6

5

4

3

Time      2

(j index)   1

0
0   1   2   3   4   5   6   7   8   9   10

space                  x
(i index)
Initial and Boundary Conditions
Explicit Euler method
10

9

8

7

u(0, t)
6

5
u(a, t)
= g1(t)   4                                                   = g2(t)
3

2

1

0
0    1   2    3   4   5   6   7   8   9   10

Initial conditions : u(x,0) = f(x)
Heat Equation
• Finite-difference
tj+1
t
u(x,t)                                              (i,j+1)
t                   tj
(i-1,j) (i,j)     (i+1,j)

x
xi-1     xi xi+1              x
1
Forward-difference   ut  ( ui , j  1  ui , j )
k
Central-difference            c
at time level j    cu xx  2 ( ui  1 , j  2 ui , j  ui  1 , j )
h
Explicit Method
• Explicit Euler method for heat equation
h   x  a / n , x i  ih
let 
k   t  T / m , t j  jk
1                            c
ut  cu xx       ( u i , j  1  u i , j )  2 ( ui  1 , j  2 u i , j  ui  1 , j )
k                           h

• Rearrange
ck
ui , j  1    ui , j  2 ( ui  1 , j  2 ui , j  ui  1 , j )   ck   c t
h                                      r 2 
h     x2
 ru i 1 , j  ( 1  2 r )ui , j  ru i  1 , j

Stability: 0  r  0.5
Explicit Euler Method
ui , j  1  ru i  1 , j  ( 1  2 r )ui , j  ru i  1 , j
• Stable
r  0.01  ui , j  1  0.01 ui  1 , j  0.98 ui , j  0.01 ui  1 , j
r  0.1  ui , j  1  0.1 ui  1 , j  0.8 ui , j  0.1 ui  1 , j
r  0.4  ui , j  1  0.4 ui  1 , j  0.2 ui , j  0.4 ui  1 , j
r  0.5  ui , j  1  0.5 ui  1 , j  0.5 ui  1 , j
• Unstable (negative coefficients)
r  1  ui , j  1  ui  1 , j  ui , j  ui  1 , j

r  10  ui , j  1  10 ui 1 , j  19 ui , j  10 ui  1 , j
r  100  ui , j  1  100ui 1 , j  199 ui , j  100ui  1 , j

Heat Equation: Explicit Euler Method

r = 0.5
Example: Explicit Euler Method
• Heat Equation (Parabolic PDE)
ut  cu xx ; 0  x  1
u(x,0)  20  40 x
            t                  2t
u(0,t)  20e , u( 1 , t )  60e
• c = 0.5, h = 0.25, k = 0.05

2
20e   -t                                              60e -2t
1

0
0        1        2         3          4
20 + 40 x
Example
• Explicit Euler method
ck ( 0.5 )( 0.05 )
r 2            2
 0.4
h     ( 0.25 )
ui , j  1  ru i  1 , j  ( 1  2 r )ui , j  ru i  1 , j
 0.4 ui 1 , j  0.2 ui , j  0.4 ui  1 , j
• First step: t = 0.05
 u0 ,1     20e 0.05  19.02458849

 u 1 ,1    0.4 u0 ,0  0.2 u1 ,0  0.4 u 2 ,0  0.4( 20 )  0.2( 30 )  0.4 ( 40 )  30

 u 2 ,1    0.4 u1 ,0  0.2 u 2 ,0  0.4 u 3 ,0  0.4 ( 30 )  0.2( 40 )  0.4( 50 )  40
u          0.4 u 2 ,0  0.2 u 3 ,0  0.4 u4 ,0  0.4 ( 40 )  0.2( 50 )  0.4 ( 60 )  50
 3 ,1
u          60e 0.10  54.29024508
 4 ,1
• Second step: t = 0.10
 u0 , 2    20e 0.10  18.09674836

 u1 , 2    0.4 u0 ,1  0.2 u1 ,1  0.4 u2 ,1
           0.4( 19.02458849 )  0.2( 30 )  0.4( 40 )  29.6098354


u 2 ,2     0.4 u1 ,1  0.2 u2 ,1  0.4 u3 ,1  0.4( 30 )  0.2( 40 )  0.4( 50 )  40
u          0.4 u2 ,1  0.2 u3 ,1  0.4 u4 ,1
 3 ,2
           0.4( 40 )  0.2( 50 )  0.4( 54.2924508 )  47.71609803

u4 ,2
           60e 0.20  49.12384518

29.61          40         47.72
20e -t                                                               60e -2t
30             40         50

0            1              2           3             4
20 + 40 x
Heat Equation: Time-dependent BCs

r = 0.4
Numerical Stability
• Stability for Explicit Euler Method
• It can be shown by Von Neumann
analysis that
1             1x   2
r      or  t 
2             2 c

• Switch to Implicit method to avoid
instability
Explicit Euler Method: Stability

r=1

Unstable !!
Implicit Euler method
Unconditionally Stable
10

9

8

7

u(0, t)
6

5
u(a, t)
= g1(t)   4                                                   = g2(t)
3

2

1

0
0    1     2   3   4   5   6   7   8   9   10

Initial conditions : u(x,0) = f(x)
Implicit Method
• Finite-difference
tj+1
t (i-1,j+1) (i,j+1) (i+1,j+1)
T(x,t)
t                       tj
(i,j)

x
xi-1      xi xi+1            x
Forward-difference u  1 ( u
t        i , j  1  ui , j )
k
Central-difference         c
at time level j+1 cu xx  2 ( ui  1 , j  1  2 ui , j  1  ui  1 , j  1 )
h
Implicit Euler Method
• Implicit Euler method for heat equation
1                            c
( ui , j  1  ui , j )  2 ( ui  1 , j  1  2 u i , j  1  ui  1 , j  1 )
k                           h
 ru i 1 , j  1  ( 1  2 r )ui , j  1  ru i  1 , j  1  ui , j

• Tridiagonal matrix (Thomas algorithm)
1  2 r   r                                  u1 , j  1   u1 , j  ru0 , j  1 
 r     1  2r   r                          u                      u2 , j          
                                              2 , j 1  
                                       

          r   1  2r           r            u3 , j  1           u3 , j          
                                            
                                     r                                          
                                       

                                 r 1  2 r  un 1 , j  1  un1 , j  ru n , j  1 
                                        

• Unconditionally stable
Implicit Euler Method

r=2

Unconditionally stable
Example: Implicit Euler Method
• Heat Equation (Parabolic PDE)
ut  cu xx ; 0  x  1
u(x,0)  20  40 x
            t                  2t
u(0,t)  20e , u( 1 , t )  60e
• c = 0.5, h = 0.25, k = 0.1

1
20e   -t                                              60e -2t

0
0        1        2         3          4
20 + 40 x
Example
• Implicit Euler method
ck ( 0.5 )( 0.10 )
r 2            2
 0.8
h     ( 0.25 )
(  r )ui 1 , j  1  ( 1  2 r )ui , j  1  (  r )ui  1 , j  1  ui , j
( 0.8 )ui 1 , j  1  ( 2.6 )ui , j  1  ( 0.8 )ui  1 , j  1  ui , j

1  2 r         r              0   u1, 1  u1,0  ru0 , 1 
 r           1  2r                  u   u
 r   2 , 1   2 ,0

                                                                     
 0
                 r          1  2 r  u 3 , 1  u 3 ,0  ru 4 , 1 
                            
• Solve the tridiagonal matrix
 2.6  0.8    0   u1 ,1   30  0.8( 20e 0.1 )
  0.8 2.6  0.8  u   40                     
                   2 ,1                       
 0
        0.8 2.6  u3 ,1  50  0.8( 60e 0.2 )
                             
 u1 ,1   28.95515793
                     
 u2 ,1    38.50751457
u  46.19426454
 3 ,1                

1           28.96    38.51     46.19

20e -t                                         60e -2t

0
0     1        2         3        4
20 + 40 x
Crank-Nicolson method
Implicit Euler method : first-order in time
Crank-Nicolson : second-order in time
10

9

8

u(0, t)
7
u(a, t)
6

= g1(t)   5
= g2(t)
4

3

2

1

0
0     1   2   3   4   5   6   7   8   9   10

Initial conditions : u(x,0) = f(x)
Crank-Nicolson Method
• Crank-Nicolson method for heat equation
• Average between two time levels
1                            c
( ui , j  1  ui , j )     2
( ui  1 , j  2 ui , j  ui  1 , j )
k                           2h
c
 2 ( ui  1 , j  1  2 ui , j  1  ui  1 , j  1 )
2h

• Tridiagonal matrix
r                                      r                r                             r
     ui 1 , j  1  ( 1  r )ui , j  1  ui  1 , j  1  ui 1 , j  ( 1  r )ui , j  ui  1 , j
2                                      2                2                             2

• Unconditionally stable (neutrally stable)
• Oscillation may occur
General Two-Level Method
• General two-stage method for heat equation
1                           c
( ui , j  1  ui , j )  2 ( ui  1 , j  2 ui , j  ui  1 , j )
k                           h
c( 1   )
       2
( ui  1 , j  1  2 ui , j  1  ui  1 , j  1 )
h

• Weighted-average of spatial derivatives
between two time levels n and n+1

 λ  0 : im plicit Euler schem e

 λ  1 : explicit Euler schem e
 λ  1/2 : Crank - Nicolson schem e

Example: Crank-Nicolson Method
• Heat Equation (Parabolic PDE)
ut  cu xx ; 0  x  1
u(x,0)  20  40 x
            t                  2t
u(0,t)  20e , u( 1 , t )  60e
• c = 0.5, h = 0.25, k = 0.1

1
20e   -t                                              60e -2t

0
0        1        2         3          4
20 + 40 x
Example
• Crank-Nicolson method
ck ( 0.5 )( 0.10 )
r      2
                  2
 0.8
h             ( 0.25 )
r                                   r                  r                          r
 ui  1, j  1  (1  r)ui, j  1  ui  1, j  1  ui  1, j  (1  r)ui, j  ui  1, j
2                                   2                  2                          2
 0.4ui  1, j  1  1.8ui, j  1  0.4ui  1, j  1  0.4ui  1, j  0.2ui, j  0.4ui  1, j

• Tridiagonal matrix (r = 0.8)
         r                   r                          r       r      
 1 r        0                  u0 ,0  ( 1  r )u1 ,0  u2 ,0  u0 ,1 
 r
2
 u 1 ,1   2
r
2       2

r                                        r              
     1  r   u2 ,1    u1 ,0  ( 1  r )u2 ,0  u3 ,0              
   2           2              2                         2
 0                  u 3 ,1   r
 

1  r
r                                                r       r
                                u2 ,0  ( 1  r )u3 ,0  u4 ,0  u4 ,1 
         2                   2
                           2       2      

• Solve the tridiagonal matrix
 1.8  0.4    0   u1 ,1  0.4( 20 )  0.2( 30 )  0.4( 40 )  0.4( 20e 0.1 )
  0.4 1.8  0.4  u   0.4( 30 )  0.2( 40 )  0.4( 50 )                     
                   2 ,1                                                      
 0
        0.4 1.8  u3 ,1  0.4( 40 )  0.2( 50 )  0.4( 60 )  0.4( 60e 0.2 ) 
                                                            
 37.23869934                 u1 ,1   29.42144598
                                                 
 40                         u2 ,1    39.29975855
69.64953807                 u  47.42746748
                             3 ,1                

1                29.42        39.30         47.43

20e -t                                                             60e -2t

0
0          1            2             3             4
20 + 40 x
Implicit Euler method

r=2

Unconditionally stable
Heat Equation with Insulated
Boundary
• No heat flux at x = 0 and x = a
x=0                                 x=a

ux(0,t) = 0                        ux(a,t) = 0

u      2u
c       ; 0  x  a, 0  t  T
t     x 2

u(x,0)  f(x),    0 xa
u x (0,t)  0
                  ,   0 t T
u x ( a , t )  0
Insulated Boundary
• No heat flux at x = a
un  1 , j  un  1 , j
ux ( a , t )                               0  un  1  un  1
x n 1  x n 1
ux(a,t)=0

xn-1    xn     xn+1

x=a
un , j  1  ru n1 , j  ( 1  2 r )un , j  ru n 1 , j
 2 ru n1 , j  ( 1  2 r )un , j

```
To top