MATLAB: How to solve an ODE with time-dependent parameters in MATLAB

dependentMATLABodetime

Consider the following ODE with time-dependent parameters:
 y'(t) + f(t)y(t) = g(t)
and the given initial condition:
 y(0) = 1
This is an example of an ODE with time-dependent terms. Suppose the time-dependent terms are defined only through the set of data points given in two vectors. Which of the MATLAB ODE solvers should I use, and how do I set up this problem?

Best Answer

To solve this ODE, you must pass the data sets to the derivative function as additional parameters. When the ODE solver calls the derivative function, it will pass a specified time as the first input argument. You must then interpolate the datasets to obtain the value of the time-dependent terms at the specified time. This is performed within the following function, called myODE.
 
function dydt = myODE(t, y, ft, f, gt, g)
f = interp1(ft, f, t); % Interpolate the data set (ft, f) at times t
g = interp1(gt, g, t); % Interpolate the data set (gt, g) at times t
dydt = -f.*y + g; % Evalute ODE at times t
The function f is defined through the n-by-1 vectors tf and f, and the function g is defined through the m-by-1 vectors tg and g.
Now, you can refer to myODE within a call to a MATLAB ODE solver. Assume that the time-dependent parameters ft and gt are defined within the data sets generated by the following code.
ft = linspace(0, 5, 25); % Generate t for f 
f = ft.^2 - ft - 3; % Generate f(t)
gt = linspace(1, 6, 25); % Generate t for g
g = 3*sin(gt - 0.25); % Generate g(t)
The following code uses the ODE45 function to solve this time-dependent ODE.
TSPAN = [1 5]; % Solve from t=1 to t=5

IC = 1; % y(t=0) = 1

[T Y] = ode45(@(t,y) myODE(t, y, ft, f, gt, g), TSPAN, IC); % Solve ODE

Note that if you are using a version of MATLAB prior to MATLAB 7.0 (R14), you will need to pass the four additional parameters ft, f, gt, and g, into the ODE solver as follows.
TSPAN = [1 5]; % Solve from t=1 to t=5
IC = 1; % y(t=0) = 1
[T Y] = ode45(@myODE, TSPAN, IC, [], ft, f, gt, g); % Solve ODE
Now you can plot the solution y(t) as a function of time.
plot(T, Y);
title('Plot of y as a function of time');
xlabel('Time'); ylabel('Y(t)');
Related Question