Code covered by the BSD License

# 5 ways to Solve IVP ODEs in MATLAB

### Sulaymon Eshkabilov (view profile)

The scripts show 5 ways of solving IVP ODEs in MATLAB environment

IVP_ODEsols_5_ways.m
```%% Solving IVP ODEs using Symbolic MATH "dsolve", Laplace Transforms and
% MuPAD for analytic solutions and ODE45, ODE113 and Simulink for numeric
% solutions of the given 2nd order ODE is shown in this script.
%
% Given a 2nd order ODE: 3*y''+3*y'+9*y=15*t with ICs: y(0)=3, y'(0)=-2;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                by Sulaymon L. ESHKABILOV, Ph.D
%                         December, 2011
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

close all
clear
%% 1st Solution method: DSOLVE
Solution_dsolve=dsolve('3*D2y+3*Dy+9*y-15*t', 'y(0)=3', 'Dy(0)+2=0');
display('Symbolic MATH dsolve function output is: ')
disp(Solution_dsolve)
% display plot of the solution of dsolve and compare it with LT
figure;
subplot(211)
ezplot(Solution_dsolve)
title('Solution found with DSOLVE '); grid
%% 2nd Solution Method: Laplace Transforms
% Define symbolic variables' names
syms t s Y
ODE2nd='3*D(D(y))(t)+3*D(y)(t)+9*y(t)-15*t';

% Laplace Transforms
lt_A=laplace(ODE2nd, t, s);

% Substitute Initial Conditions and initiate an unknown Y
lt_A=subs(lt_A,{'laplace(y(t),t,s)','y(0)','D(y)(0)'},{Y,3,-2});

% Solve for Y unknown
Y=solve(lt_A, Y);

display('Laplace Transforms of the given 2nd Order ODE with ICs')
disp(Y)

Solution_Laplace=ilaplace(Y);
display('Solution found using Laplace Transforms ')
disp(Solution_Laplace)
subplot(212)
ezplot(Solution_Laplace); grid
title('Solution found with Laplace Transforms')
hold off

% Plot two solutions in one plot area

t=0:pi/100:6*pi;
SDsol=eval(vectorize(Solution_dsolve));
LTsol=eval(vectorize(Solution_Laplace));
figure;
plot(t, SDsol, 'bx-'); grid
hold on; xlabel('t'); ylabel('solution values')
plot(t, LTsol, 'ro-');

% % % ooO:=ode({3*y''(t)+3*y'(t)+9*y(t)=15*t, y(0)=3, y'(0)=-2}, y(t))
% % % solutions:=ode::solve(ooO)
% % % plot(solutions, #G, #L)

%% 4th Solution Method: ODE45, ODE113 solvers to get numeric solutions
% Help: to get this cell executed the following function file needs to be
% created in separate (e.g. compareODEsols.m):
% % % % % function dydt=compareODEsols(t, y)
% % % % % % This is a 2nd order ODE IVP example
% % % % % %    3*y''+3*y'+9*y=15*t with ICs: y(0)=3, y'(0)=-2;
% % % % %
% % % % % dydt=[y(2);-(y(2)+3*y(1))+5*t ];
% % % % % return

ICs=[3, -2];
timeSPAN=0:pi/100:6*pi;
options_ODE=odeset('RelTol', 1e-6, 'AbsTol', 1e-8);
[t, y_ODE45]=ode45(@compareODEsols, timeSPAN, ICs, options_ODE);
plot(t, y_ODE45(:,1), 'm+-');
[t, y_ODE113]=ode113(@compareODEsols, timeSPAN, ICs, options_ODE);
plot(t, y_ODE113(:,1), 'k<-');
title('Comparison analysis of ODE sols with ODE45, ODE23, dsolve, Laplace Transforms & Simulink')

% NB: SIMULINK models uses ODE solver ODE23s in simulation of its *.mdl
% model. A solver type can be altered if needed for a better/faster solution search.
% This cell first recalls and executes Simulink model ODEsols_COMPARE.mdl
% that is located and it to be put in currect directory of MATLAB from the archived folder.
sim('ODEsols_COMPARE.mdl');
SIM_sols=ODEsols_SIM.signals.values; time=tout;
plot(time, SIM_sols, 'cd-');