ODE - Computational Methods - Lecture Slides, Slides of Computational Methods

Download ODE - Computational Methods - Lecture Slides and more Slides Computational Methods in PDF only on Docsity!

Chp9: ODE Solns

By MATLAB

Learning Goals

• Use MATLAB’s ODE Solvers to find Solutions to Ordinary Differential Eqns
• When Possible use Analytical ODE Solutions to Check the ACCURACY of the MATLAB Numerical Soln
• Make Approximations to perform a REALITY CHECK on MATLAB Solutions to NONLinear ODEs

Solver Summary

Solver Problem Type Order of Accuracy When to Use ode45 Nonstiff Medium Most of the time. This should be the first solver you try. ode23 Nonstiff Low For problems with crude error tolerances or for solving moderately stiff problems. ode113 Nonstiff Low to high For problems with stringent error tolerances or for solving computationally intensive problems.

ode15s Stiff Low to medium If ode45 is slow because the problem is stiff.

ode23s Stiff Low If using crude error tolerances to solve stiff systems and the mass matrix is constant. ode23t Moderately Stiff Low For moderately stiff problems if you need a solution without numerical damping. ode23tb Stiff Low If using crude error tolerances to solve stiff systems.

MATLAB ODE Format - 1

• MATLAB ODE Solver Form for Multiple Dependent Variables (MultiVar Probs)

1 2 0

2 1 2 2 0 2

2

1 1 2 1 0 1

1

m m m m

m

m

m

f t y y y y t b dt

dy

f t y y y y t b dt

dy

f t y y y y t b dt

dy

m-Eqns (1st^ order ODEs) in m-Unknowns

MATLAB ODE Format - 2

• In Vector Form (saves Writing time; does NOT make Solution Easier)

y = f y y = b

• ( t , ), ( t 0 )

 If we have written ROW vectors for the y & b quantities can Transpose to Column-Vectors

yT = [ y 1 (^) ; y 2 ;; ym ]

b T =^ [ b^1 (^) ; b 2 ;; bm ]

MATLAB Solution

• To solve this system using MATLAB the functions f 1 , f 2 , …, fm must be provided to the computer, along with the initial values of the variables; i.e., t 0 and b 1 ,b 2 , …, bm
• The functions f 1 , f 2 , …, fm are input using a function which we have to name, say fcn_vec. Then stored in an m-file called in this case fcn_vec.m
• The initial values of the y’s are stored in a one dimensional array, say y

Plotting the Solution

• MATLAB can also be used to plot the ODE Solution results.
• For example: plot(t,y) gives a plot of ALL components of the solution y 1 , y 2 , …, y (^) m , as a function of t
• Alternatively plot(t,y(:,1))gives a plot of y 1 as a function of t

Comments on ode

• Note that in its simplest form ode45 chooses its OWN time step and VARIES the time step according to how fast the solution is changing.
• Thus ode45 generates solution values at a sequence of times t 1 , t 2 , t 3 , … given by t (^) k+1= t (^) k +Δt (^) k , with Δtk selected by ode45. - Thus EACH component of the solution y 1 , y 2 , …, y (^) n is ITSELF a vector containing values such as y 1 (t1), y 1 (t2), y 1 (t3), …, then y 2 (t1), y 2 (t 2 ), y 2 (t3), …, then y 3 (t1), y 3 (t 2 ), y 3 (t3), etc.

Example  ode45 (2)

• With the Xform (^)  Subbing in the x 1 & x (^2)
  • •• • x 1 = y x 2 = y y = x 2  Find that

2

1 (^1) dt y x

dy dt

dx x = = = =

• •

 ReArranging the ODE to isolate Highest order term

x 2 (^) = sin t − 5 x 1 − 2 x 2

• • • y = sin t − 5 y − 2 y

 Thus the 1 st^ Order Eqn System in 2 Vars

1 2

2 2

2

1 1

sin t 5 x 2 x dt

dx x

x dt

dx x

= = − −

= =

0 0 73

2

1 = =

= =

y x

y x

Example  ode45 (3)

• Note that applying this Xform

 And also dx 2 /dt

x 1 = y x 2 = y

 Converted the SINGLE 2 nd^ Order ODE to a LINEAR System of TWO 1 st Order ODEs

 Sub into ODE

• •• •  = = 

  

=  2

(^2) y y x dt

d dt

dx

2

(^1) y x dt

dy dt

dx = = =

y + 2 y + 5 y = sin t

• • •

x x t dt

dx (^2) + 2 2 + 5 1 = sin

Example  ode45 (5)

• Thus the Transformation to State-Var Form of Two 1st Order ODEs y + 2 y + 5 y = sin t
• • • &

0

=

y

dt

dy

y

t

( )

sin 5 2 ( 0 ) 2

0 73

1 2 2

2

2 1

1

= − − = =

= = =

t x x x t dt

dx

x x t dt

dx

 The final xForm

Example  ode45 (6)

• Compare xForm to Slide-

( )

sin 5 2 (^0 )^2

1 2 2

2

2 1

1

t x x x t dt

dx

x x t dt

dx

( , , ), ( ) ,

( , , ), ( ) ,

2 1 2 2 0 2 2

1 1 2 1 0 1 1

f t y y y t b dt

dy

f t y y y t b dt

dy

= =

= =

Example  ode45 (8)

% Bruce Mayer, PE * 05Nov11 % ENGR25 * Lec24 on MATLAB ODE solvers % file = Demo_ODE_Lec24.m % Revised to include a set of non-zero ICs % %This script file calls FUNCTION xdot_lec % clear % clear memory % % CASE-I => set the IC's y(0) & dy(0)/dt as COL Vector % y0=[0; 0]; % comment-out if Not Used % CASE-II => set the IC's y(0) & dy(0)/dt as COL Vector y0=[-0.19; -0.73]; % Comment-Out if Not Used % % Default Time Interval of 20 Time-Units; user can change this tmax = input('input tmax = ') trng = [0, tmax]; % % %Call the ode45 routine with the above data inputs [t,y]=ode45('xdot_lec24', trng, y0); % %Plot the first column of the solution “matrix” %giving x1 or y. plot(t,y, 'LineWidth', 2), xlabel('t'), ylabel('ODE Solution y(t) & dy/dt'),... title('ODE Example - Lecture24'), grid, legend('y(t)','dy/dt')

ODE Example Result (0 for ICs)

0 2 4 6 8 10 12 14 16 18 20

-0.

-0.

-0.

-0.

-0.

0

Time, t

ODE Solution, y(t)