MATLAB: How to save a variable computed inside ODE

odeode45odesetoutput function

Inside an ode function I perform some calculations to compute a variable. I would like to save the values of that variable in a matrix so after the ode is done I can plot it vs time. Because it would not be useful to save the values during the integration, since the function is evaluated several times for each step and a step could be rejected, I followed michio's suggestion to use an Output function.
His solution is also described here: Get variable out of ode 45
function status = myOutputFcn(time, x, flag)
% q: value of variable at each succesfull time step
% q_save: matrix of variable values I want to save for plotting after ode is done
persistent q_save
switch flag
case 'init'
q_save = q;
case ''
q_save = [q_save; q];
case 'done' % when it's done
assignin('base','q_save',q_save); % get the data to the workspace.
end
status = 0;
end
The problem I have with this approach is that after the ode finishes, the size of the time vector is different from the size of the variable 'q_save'. The Output function is called less times than the succesfull steps. For example, time=800×1 and q_save=150×1.
I have also tried declaring tspan in two ways
tspan = linspace(t0,tf,100)
tspan = [t0 tf]
So, how can I save a variable computed inside ODE for every succesfull ode step?
Meaning that I want a value of q for every (time,x) pair that the ode gives me as a solution.

Best Answer

assignin is a bad programming technique. Producing variables remotely in another workspace is not clean. A better method is to let the integration run regularily and request the wanted values afterwards:
[t, y] = ode45(@fcn, tSpan, y0);
function dy = fcn(t, y)
p = sin(y(2)); % <- The wanted parameter

dy = [2 * t; p / y(1)];
end
If this is the function to be integrated, it must be modified such, that it accepts vectors as input and replies p as output:
function [dy, p] = fcn(t, y)
p = sin(y(:, 2)); % <- The wanted parameter
dy = [2 .* t; p ./ y(:, 1)];
end
Then:
[t, y] = ode45(@fcn, tSpan, y0);
[~, p] = fcn(t.', y.');
This is much easier than letting the output function create variables remotely.
For some cases the output function is the best choice. Then provide the value as output instead of assignin
function status = myOutputFcn(time, x, flag)
persistent q_save
status = 0;
switch flag
case 'init'
q_save = q;
case ''
q_save = [q_save; q];
case 'done' % when it's done
status = q_save;
q_save = [];
end
end
But now to the actual problem: Unfortunately the ODE integrators reply different time steps, if the output is obtained as a struct. You did not show, how you call the integrator, so the most important detail is missing. (In opposite to this, the code or the output function is not a part of the problem.)
If you want to call the integrator in the form sol = ode45(___), and use a tspan defined as vector also, you have to modify ode45 such that the Output function is called regularily as for the [t,y,te,ye,ie] = ode45(odefun,tspan,y0,options) form. I consider this a bug and sent a bug report to MathWorks some years ago.
You wrote:
"I have also tried declaring tspan in two ways"
tspan = linspace(t0,tf,100)
tspan = [t0 tf]
And what did you observe as effect? How did you call which integrator?