Hello,
I am trying to use ode45 to implement 4th order Runge-Kutta on a system of first order differential equations. Within the body of my code, I call ode45 as follows:
[times, T] = ode45(@(time, T) interpolate(time, A, ... int_data, tmax, temp_dist, a, b), ... subt, temp_dist);
where @(time, T) interpolate(…) is my function handle, subt is my interval of integration, and temp_dist is my initial condition. interpolate.m is defined as follows:
function dTdt = interpolate(t, A, int_data, tmax, temp_dist, a, b) m = interp1([1:tmax], int_data, t, 'linear', 'extrap'); % interpolate int_data at time t
m(m < 0) = 0; % none of the values in m can be negative
% define matrix A at interpolated times (the first row of A depends on
% t)
A(1,1) = -1*b-a*m(3)-a*m(4); A(1,2) = b; A(1,end) = a*(m(1) + m(2) + m(3)*m(5) + ... m(4)*m(5)); % define output (evaluate system of ODE's at time t)
dTdt = A*temp_dist;
where t is the time input by ode45, A is a matrix of coefficients whose first row changes with t, int_data is a matrix of information used to determine the first row of A, temp_dist is a vector of initial conditions, and tmax, a and b are used to determine the first row of A.
As you probably know, ode45 calls interpolate.m each time it needs a new value of dTdt for a specific time t. My problem is that temp_dist needs to be updated in between calls. After dTdt is calculated, temp_dist should be updated to:
temp_dist = temp_dist + dTdt;
This new value of temp_dist should be the one used by ode45 the next time interpolate.m is called.
I'm stuck trying to implement this change. The obvious thing might be to redefine temp_dist within ode45.m, but the source code for this function is pretty incomprehensible to me. I've tried redefining temp_dist as shown above within interpolate.m, but when ode45 calls interpolate.m the next time, it's back to its original value.
Does anyone know how to solve this problem? My code is written in R2014a.
Thanks!!
Best Answer