Azzera filtri
Azzera filtri

manually run ode45 step by step - and impose filter to the half-way solution

1 visualizzazione (ultimi 30 giorni)
Hi everybody~
I'm solving a ode:
%%code 1
t = [tt1 tt2 tt3 ... tt100];
[t, y] = ode45(@equation, t, y1, options); % y(t=1)=y1
Now I wanna to run it step by step. I wanna verify the above code is equivalent to
%%code 2
t = [tt1 tt2 tt3 ... tt100];
for i = 1:99
[t(i+1) y(i+1)] = ode45(@equation, [t(i) t(i+1)], y(i), options);
end
The reason I want step-by-step is, I want to impose filters to y in the middle. Say, at t=t(5), I do y(5) = y(5)*filter, and then continue to solve y(6).
Let me know whether code 1 is equivalent to code 2.
Any thoughts of other ways to do this work is appreciated. I think this is called "numerical filtering" - but didn't find anything to read about it. Any recommended reading is appreciated. if true % code endThanks~

Risposta accettata

Jan
Jan il 20 Giu 2013
The calculations are equivalent, but not exactly identical: While in the first version ODE45 uses the stepsize from the former step after reaching one of the intermediate points, the 2nd method restarts the integrations with a new estimated stepsize. Due to rounding and local discretization errors, the results wil be slightly different - except if the solution is instable and the small deviations are amplified, which can cause large differences. But then the "result" of the integration is doubtful at all.
  2 Commenti
Yuji Zhang
Yuji Zhang il 21 Giu 2013
Hi Jan~
Thank u for your help.
I learn from "help" that, t= t(i) is just the "record points" I want y(i) to be recorded. code45 automatically determines internal step-length of delta_t. You were talking about how the step-length is determined. Is there something I can read about it?
Let me know. Thank u!
Jan
Jan il 24 Giu 2013
This is a perfect question for an internet serach engine: Simply ask your favorite engine for "Runge Kutta stepsize control"

Accedi per commentare.

Più risposte (0)

Categorie

Scopri di più su Mathematics in Help Center e File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by