Solutions by MATLAB - Computational Methods - Lecture Slides, Slides of Calculus for Engineers

Download Solutions by MATLAB - Computational Methods - Lecture Slides and more Slides Calculus for Engineers in PDF only on Docsity!

Engr/Math/Physics 25

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 orfor 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 stiffsystems and the mass matrix is constant.

ode23t Moderately Stiff Low For moderately stiff problems if you need asolution 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 , …, ym , 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 tk+1= tk+Δtk, with Δtk selected by ode45. - Thus EACH component of the solution y 1 , y 2 , …, yn is ITSELF a vector containing values such as y 1 (t 1 ), y 1 (t 2 ), y 1 (t 3 ), …, then y 2 (t 1 ), y 2 (t 2 ), y 2 (t 3 ), …, then y 3 (t 1 ), y 3 (t 2 ), y 3 (t 3 ), 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

x dx

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

( , , ), ( ) ,

( , , ), ( ) ,

(^2212202)

(^1112101)

f t y y y t b dt

dy

f t y y y t b dt

dy

 

 

Example  ode45 (8) % Bruce Mayer, PE% ENGR25 * Lec24 on MATLAB ODE solvers* 05Nov % file =% Revised Demo_ODE_Lec24.m 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;% CASE-II 0]; % comment-out => set the IC's y(0) & dy(0)/dt as COL Vectorif Not Used y0=[-0.19; -0.73]; % Comment-Out if Not Used% % Defaulttmax = input('input tmax = Time Interval of 20 Time-Units; user can change this ') trng% = [0, tmax]; % %Call the [t,y]=ode45( ode45'xdot_lec24' routine (^) ,with the above data inputs trng, y0); % %Plot the first column of the solution “matrix” %giving x1 or plot(t,y, 'LineWidth y. ', 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)

ODE Example - Lecture