Solve nonstiff differential equations — medium order method
tspan = [t0 tf], integrates the system of
differential equations from
y0. Each row in the solution
y corresponds to a value returned in column
All MATLAB® ODE solvers can solve systems of equations of
the form ,
or problems that involve a mass matrix, .
The solvers all use similar syntaxes. The
only can solve problems with a mass matrix if the mass matrix is constant.
solve problems with a mass matrix that is singular, known as differential-algebraic
equations (DAEs). Specify the mass matrix using the
ode45 is a versatile ODE solver and is the
first solver you should try for most problems. However, if the problem
is stiff or requires high accuracy, then there are other ODE solvers
that might be better suited to the problem. See Choose an ODE Solver for
finds where functions of (t,y),
called event functions, are zero. In the output,
the time of the event,
ye is the solution at the
time of the event, and
ie is the index of the triggered
For each event function, specify whether the integration is
to terminate at a zero and whether the direction of the zero crossing
matters. Do this by setting the
to a function, such as
and creating a corresponding function: [
For more information, see ODE Event Location.
Simple ODEs that have a single solution component can be specified as an anonymous function in the call to the solver. The anonymous function must accept two inputs
(t,y) even if one of the inputs is not used.
Solve the ODE
Use a time interval of
[0,5] and the initial condition
y0 = 0.
tspan = [0 5]; y0 = 0; [t,y] = ode45(@(t,y) 2*t, tspan, y0);
Plot the solution.
The van der Pol equation is a second order ODE
where is a scalar parameter. Rewrite this equation as a system of first-order ODEs by making the substitution . The resulting system of first-order ODEs is
The function file
vdp1.m represents the van der Pol equation using
. The variables
are the entries
y(2) of a two-element vector,
function dydt = vdp1(t,y) %VDP1 Evaluate the van der Pol ODEs for mu = 1 % % See also ODE113, ODE23, ODE45. % Jacek Kierzenka and Lawrence F. Shampine % Copyright 1984-2014 The MathWorks, Inc. dydt = [y(2); (1-y(1)^2)*y(2)-y(1)];
Solve the ODE using the
ode45 function on the time interval
[0 20] with initial values
[2 0]. The resulting output is a column vector of time points
t and a solution array
y. Each row in
y corresponds to a time returned in the corresponding row of
t. The first column of
y corresponds to
, and the second column to
[t,y] = ode45(@vdp1,[0 20],[2; 0]);
Plot the solutions for
plot(t,y(:,1),'-o',t,y(:,2),'-o') title('Solution of van der Pol Equation (\mu = 1) with ODE45'); xlabel('Time t'); ylabel('Solution y'); legend('y_1','y_2')
ode45 only works with functions that use two input arguments,
y. However, you can pass in extra parameters by defining them outside the function and passing them in when you specify the function handle.
Solve the ODE
Rewriting the equation as a first-order system yields
odefcn.m represents this system of equations as a function that accepts four input arguments:
function dydt = odefcn(t,y,A,B) dydt = zeros(2,1); dydt(1) = y(2); dydt(2) = (A/B)*t.*y(1);
Solve the ODE using
ode45. Specify the function handle such that it passes in the predefined values for
A = 1; B = 2; tspan = [0 5]; y0 = [0 0.01]; [t,y] = ode45(@(t,y) odefcn(t,y,A,B), tspan, y0);
Plot the results.
Consider the following ODE with time-dependent parameters
The initial condition is
. The function
f(t) is defined by the n-by-1 vector
f evaluated at times
ft. The function
g(t) is defined by the m-by-1 vector
g evaluated at times
Create the vectors
ft = linspace(0,5,25); f = ft.^2 - ft - 3; gt = linspace(1,6,25); g = 3*sin(gt-0.25);
Write a function named
myode that interpolates
g to obtain the value of the time-dependent terms at the specified time. Save the function in your current folder to run the rest of the example.
myode function accepts extra input arguments to evaluate the ODE at each time step, but
ode45 only uses the first two input arguments
function dydt = myode(t,y,ft,f,gt,g) f = interp1(ft,f,t); % Interpolate the data set (ft,f) at time t g = interp1(gt,g,t); % Interpolate the data set (gt,g) at time t dydt = -f.*y + g; % Evaluate ODE at time t
Solve the equation over the time interval
[1 5] using
ode45. Specify the function using a function handle so that
ode45 uses only the first two input arguments of
myode. Also, loosen the error thresholds using
tspan = [1 5]; ic = 1; opts = odeset('RelTol',1e-2,'AbsTol',1e-4); [t,y] = ode45(@(t,y) myode(t,y,ft,f,gt,g), tspan, ic, opts);
Plot the solution,
y, as a function of the time points,
The van der Pol equation is a second order ODE
Solve the van der Pol equation with
ode45. The function
vdp1.m ships with MATLAB® and encodes the equations. Specify a single output to return a structure containing information about the solution, such as the solver and evaluation points.
tspan = [0 20]; y0 = [2 0]; sol = ode45(@vdp1,tspan,y0)
sol = struct with fields: solver: 'ode45' extdata: [1×1 struct] x: [1×60 double] y: [2×60 double] stats: [1×1 struct] idata: [1×1 struct]
linspace to generate 250 points in the interval
[0 20]. Evaluate the solution at these points using
x = linspace(0,20,250); y = deval(sol,x);
Plot the first component of the solution.
Extend the solution to
odextend and add the result to the original plot.
sol_new = odextend(sol,@vdp1,35); x = linspace(20,35,350); y = deval(sol_new,x); hold on plot(x,y(1,:),'r')
odefun— Functions to solve
Functions to solve, specified as a function handle which defines the functions to be integrated.
dydt = odefun(t,y), for a scalar
a column vector
y, must return a column vector
corresponds to .
accept both input arguments,
even if one of the arguments is not used in the function.
For example, to solve , use the function:
function dydt = odefun(t,y) dydt = 5*y-3;
For a system of equations, the output of
a vector. Each element in the vector is the solution to one equation.
For example, to solve
use the function:
function dydt = odefun(t,y) dydt = zeros(2,1); dydt(1) = y(1)+2*y(2); dydt(2) = 3*y(1)+2*y(2);
For information on how to provide additional parameters to the
odefun, see Parameterizing Functions.
tspan— Interval of integration
Interval of integration, specified as a vector. At minimum,
be a two element vector
[t0 tf] specifying the
initial and final times. To obtain solutions at specific times between
use a longer vector of the form
The elements in
tspan must be all increasing or
The solver imposes the initial conditions,
tspan(1), then integrates from
tspan has two elements,
tf], then the solver returns the solution evaluated at each
internal integration step within the interval.
tspan contains more than two
[t0,t1,t2,...,tf], then the solver returns
the solution evaluated at the given points. This does not affect
the internal steps that the solver uses to traverse from
Thus, the solver does not necessarily step precisely to each point
tspan. However, the solutions produced
at the specified points are of the same order of accuracy as the solutions
computed at each internal step.
Specifying several intermediate points has little effect on the efficiency of computation, but for large systems it can affect memory management.
The solution obtained by the solver might be different depending
on whether you specify
tspan as a two-element vector
or as a vector with intermediate points. If
several intermediate points, then they give an indication of the scale
for the problem, which can affect the size of the initial step taken
by the solver.
3 5 7 9 10]
y0— Initial conditions
Initial conditions, specified as a vector.
be the same length as the vector output of
y0 contains an initial condition for each
equation defined in
options— Option structure
options = odeset('RelTol',1e-5,'Stats','on','OutputFcn',@odeplot) specifies
a relative error tolerance of
1e-5, turns on the
display of solver statistics, and specifies the output function
plot the solution as it is computed.
t— Evaluation points
Evaluation points, returned as a column vector.
tspan contains two elements,
t contains the internal evaluation
points used to perform the integration.
tspan contains more than two
t is the same as
Solutions, returned as an array. Each row in
to the solution at the value returned in the corresponding row of
te— Time of events
Time of events, returned as a column vector. The event times
te correspond to the solutions returned in
ie specifies which event occurred.
ye— Solution at time of events
Solution at time of events, returned as an array. The event
te correspond to the solutions returned
ie specifies which
ie— Index of vanishing event function
Index of vanishing event function, returned as a column vector.
The event times in
te correspond to the solutions
which event occurred.
sol— Structure for evaluation
Structure for evaluation, returned as a structure array. Use
this structure with the
to evaluate the solution at any point in the interval
sol structure array always includes
Row vector of the steps chosen by the solver.
Solutions. Each column
Additionally, if you specify
Events option and events are detected, then
includes these fields:
Points when events occurred.
Solutions that correspond to events in
Indices into the vector returned by the function specified
 Dormand, J. R. and P. J. Prince, "A family of embedded Runge-Kutta formulae," J. Comp. Appl. Math., Vol. 6, 1980, pp. 19–26.
 Shampine, L. F. and M. W. Reichelt, "The MATLAB ODE Suite," SIAM Journal on Scientific Computing, Vol. 18, 1997, pp. 1–22.