How to make ode45 store/return only results of every n-th time step to reduce memory use?
Afficher commentaires plus anciens
Hi,
I am trying to run a relatively long simulation with a small time step with ode45, meaning the resulting y matrix is very big. I want to keep the time step small for accuracy, but I do not need the result at each time step. When running a longer simulation the memory becomes a problem so makin the output sparser would be very beneficial.
Thank you!
1 commentaire
I assumed that with "keeping the time step small" the OP meant that "RelTol" and "AbsTol" in his simulation were chosen small, not that the large output was a consequence of setting the "tspan" vector accordingly. The latter will not enhance precision as you correctly pointed out.
Réponse acceptée
You can specify the times when the solution should be saved in the input variable "tspan" to "ode45".
And you should test whether "ode15s" can make larger time steps - maybe your system is stiff so that "ode45" is not suited as solver.
Plus de réponses (2)
Walter Roberson
le 21 Fév 2025
1 vote
Set your tspan to be a vector of three elements (more correctly, use more than 2 elements). When you call ode45(), set it up to ignore the t and y outputs.
Use odeset to create ode*() options that include OutputFcn . Have the output function use a persistent variable to count how many times it has been invoked. "most" of the time throw away the data. Every N'th time it is called, store the data in an accessible place.
OutputFcn will be called for each successful time step; you ignore most of the time steps and record every N'th.
Meanwhile the internal recording of the outputs will only happen at the times in your tspan, and you will have stripped down your tspan to the minimum needed so that it does not choose to record everything (tspan with 2 entries.)
11 commentaires
Seems OutputFcn is only called at the three times specified in "tspan". Or is there something wrong with the example code ?
f = @(t,y) exp(-0.05*t)*sin(t);
tspan = [0 0.5 1];
options = odeset('OutputFcn',@myOutputFcn);
y0 = 0;
[T,Y] = ode45(f,tspan,y0,options);
function status = myOutputFcn(t,y,flag)
status = 0;
t
end
Mike Hosea
le 21 Fév 2025
Modifié(e) : Mike Hosea
le 21 Fév 2025
Clearly the interface is designed to require the user to supply the time points in situations like this, but there's no adaptivity of the output points in that case.
I don't know if the need for this is common enough to merit supporting it more direclty, but I could imagine an enhancement to support it directly, albeit in a slightly kludgy way, if refine < 0 is the number of consecutive steps to throw away.
In this case, I don't understand what "OutputFcn" is good for if the function is only called at the times specified in "tspan" (assuming "tspan" is a vector with more than two elements). Postprocessing the results of the integration could be done easily after the simulation has finished.
Mike Hosea
le 21 Fév 2025
Torsten is right. I was thinking that outputFcn was called after every successful step. In fact, it is called after every successful step in which there are any output points. If there are none, it will not get called.
If the problem isn't stiff, I guess you could minimize the number of time steps by using the highest order RK method available, i.e. ode89 and setting Refine=1.
Other than that, I think you have to specify output points up front to minimize storage.
What time points do you imagine that OUTPUTFCN shoudl be called at: all function evaluations including the invalid solutions and backtracking?
"OutputFcn" should be called after each successful integration step - independent of what is set in "tspan".
Mike Hosea
le 21 Fév 2025
@Stephen23, I was assuming that the user was using length(tspan) == 2 and then complaining that the number of output points was too large. I thought their reference to the time step was made with respect to tolerance setting. As I read your comment, I see the other interpretation that you suspected. It is plausible, but I don't know which interpretation is correct at this point.
Torsten
le 22 Fév 2025
Isn't your suggested behavior already provided by specifying TSPAN with only the end points?
No, it doesn't solve the OP's problem. Specifying TSPAN with only the end points means that all the solutions at successful integration steps will be written in the Y-matrix returned by the ode-solver - and this matrix can become very large.
In my opinion, plotting can be done after the solution is returned by the ODE-solver.
Walter Roberson
le 22 Fév 2025
"OutputFcn" should be called after each successful integration step
That's how it is documented.
- independent of what is set in "tspan".
Unfortunately the documentation does not specify the relationship between recording output and OutputFcn calls. I think it is a valid interpretation to assume that OutputFcn is called after every successful integration step regardless of whether the step is one to be recorded in the output arrays. That might turn out to be inaccurate in practice but I think it is a valid reading.
For the given example, OutputFcn is called for the three times specified in "tspan". You can see the difference if you replace "tspan" by the vector that only contains the endpoints of the integration interval.
f = @(t,y) exp(-0.05*t)*sin(t);
options = odeset('OutputFcn',@myOutputFcn);
y0 = 0;
tspan = [0 0.5 1];
[T,Y] = ode45(f,tspan,y0,options);
tspan = [0 1];
[T,Y] = ode45(f,tspan,y0,options);
T
function status = myOutputFcn(t,y,flag)
status = 0;
t
end
Mike Hosea
le 22 Fév 2025
Modifié(e) : Mike Hosea
le 22 Fév 2025
Yes, it turns out that the outputFcn is only called if there is a non-empty set of outputs to transmit to it. I guess that is in keeping with it being literally an "output function" and not, more generically, a step callback function. Probably most use cases involve length(tspan)==2.
I have reported the documentation issue.
Mike Croucher
le 21 Fév 2025
Modifié(e) : Mike Croucher
le 21 Fév 2025
I would suggest trying the new ode interface The new solution framework for Ordinary Differential Equations (ODEs) in MATLAB R2023b » The MATLAB Blog - MATLAB & Simulink
One approach that might work for you is the new solutionFcn. This is covered in the blog post above but, for completeness, here's the example I used
myODE = ode(ODEFcn=@(t,y) exp(-0.05*t)*sin(t),InitialTime=0,InitialValue=0); % Define the ODE problem
solFcn = solutionFcn(myODE,0,20) % Use it to create a solutionFcn
I can now evaluate this wherever I like within the interval defined in the solutionFcn. E.g.
solFcn(15)
This way you can keep only the exact time points you want.
2 commentaires
The OP's problem is that the ode integrator internally saves the solution for too many time steps such that available memory becomes low during a computation. I think this problem is not solved by your suggestion. Somewhere, the solution of the ODE must be saved at discrete times in order to make interpolation using the function handle "solFcn" after the computation has finished - and I don't see that this solution is somehow reduced in size by the code.
Mike Croucher
le 22 Fév 2025
Yes, you are right! Sorry for the noise
Catégories
En savoir plus sur Ordinary Differential Equations dans Centre d'aide et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
