Using ODE45 to solve a system of ODEs

I have system of ODEs that look like:

dcndt=a(t)*cn(t)+b(t)*cm(t)
dcndt=c(t)*cn(t)+b(t)*cm(t)

I know the values of a(t), b(t), and c(t) for a given time t=0:1:1000. I have created a function called:

function dy = myODE(t,y,a,b,c)
cn = y(1);
cm = y(2);
dy(1,1) = a(t)*cn+b(t)*cm;
dy(2,1) = c(t)*cm+b(t)*cn;

When I ran the following command:

[T,Y] = ode45(@(t,y) myODE(t,y,a,b,c),t,[1 0]);

I got the error message:

??? Attempted to access a(1); index must be a positive integer or logical.
Error in ==> myODE at 5
dy(1,1) = a(t)*cn+b(t)*cm;
Error in ==> @(t,y)myODE(t,y,a,b,c)
Error in ==> ode45 at 324
    f(:,2) = feval(odeFcn,t+hA(1),y+f*hB(:,1),odeArgs{:});

I guess that it means that need to access the a1 in any time value, but I only know a, b, and c for the specific times 0:1:1000. Does anybody knows how to solve this problem?

[Walter Roberson]

Are you attempting to find an ODE that is only defined over a discrete grid, or are you attempting to find an ODE that is continuously defined but which you have only sampled at specific points and the values at other points are not (feasibly) available?

Answer 1

[Kelly Kearney]

You could try interpolating within your a, b, and c datasets.

 function dy = myODE(t,y,a,b,c)

 at = interp1(0:1000, a, t);
 bt = interp1(0:1000, b, t);
 ct = interp1(0:1000, c, t);

 cn = y(1);
 cm = y(2);
 dy(1,1) = at*cn+bt*cm;
 dy(2,1) = ct*cm+bt*cn;

Here I'm demonstrating a linear interpolation, but you'll have to decide what is most appropriate to your problem.

Answer 2

[Daniel]

Kelly, I cannot interpolate because a, b, and c are stochastic variables.

[Kelly Kearney]

ODE45 expects that your function is defined for all t. So you have to somehow figure out how you want to define a, b, and c over the entire domain. I don't see how stochasticity excludes interpolation, but then I'm not sure what your parameters represent. Perhaps a nearest-neighbor interpolation would be more appropriate than a linear one?

The only way to restrict the calculation to your set of integer t values would be to use a first-order solver, I think.