Time derivative of parameters within ODE solvers
4 次查看(过去 30 天)
显示 更早的评论
I have an ODE which has a parameter whose 1st and 2nd order time derivatives are also included in the ODE:
function dy = ODE(t,y)
f1 = myfun(y(t));
dy = y + f1 + df1/dt + ddf1/dt^2;
end
Unfortunately, the function of analytical derivatives of myfun is not available. Therefore, df1 and ddf1 can be computed numerically, only. Given that the time step in Matlab ode solvers is not fixed, I wonder if there is a way to numerically compute df1 and ddf1.
0 个评论
采纳的回答
Torsten
2018-11-23
编辑:Torsten
2018-11-23
function dy = ODE(t,y)
dt = 1e-8;
fm = myfun(t-dt);
f = myfun(t);
fp = myfun(t+dt);
df = (fp - fm) / (2 * dt);
ddf = (fp - 2 * f + fm) / dt^2;
dy = y + f + df + ddf;
end
6 个评论
Torsten
2018-11-23
编辑:Torsten
2018-11-23
Let
z = f(y)
the value that "myfun" returns for argument y.
Then
dz/dt = df/dy * dy/dt
d^2z/dt^2 = d^2f/dy^2 * (dy/dt)^2 + df/dy * d^2y/dt^2
Inserting into your differential equation gives
dy/dt = y + f + df/dy * dy/dt + d^2f/dy^2 * (dy/dt)^2 + df/dy * d^2y/dt^2
or
df/dy * d^2y/dt^2 + (df/dy - 1) * dy/dt + d^2f/dy^2 * (dy/dt)^2 + y + f = 0
Now you can approximate df/dy and d^2f/dy^2 as I suggested above and solve the system (z1 = y, z2 = dy/dt)
z1' = z2
z2' = -((df/dy - 1) * z2 + d^2f/dy^2 * (z2)^2 + z1 + f)/(df/dy)
using ODE45, e.g.
更多回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!