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

조회 수: 1 (최근 30일)
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~

채택된 답변

Jan
Jan 2013년 6월 20일
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
Yuji Zhang
Yuji Zhang 2013년 6월 21일
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 2013년 6월 24일
This is a perfect question for an internet serach engine: Simply ask your favorite engine for "Runge Kutta stepsize control"

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Mathematics에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by