首页 > 解决方案 > 带有差分方程系统的 ODE45 时间因变量

问题描述

我有一个需要插值的时间相关变量,并发送到我的 ODE。我的 EOM 存储在一个函数 .m 文件中,我需要将插值发送给它,因此每个值都可以在集成中使用。我该怎么做呢?

我需要插值的变量是 k。其中 k 遵循这个公式

k = (h^2*(a+b)-b*h*(a+2*b)*cos(a.*time./b)+b^3)./...
((a+b)*(-2*b*h*cos(a.*time./b)+b^2+h^2).^2);

function [ sv ] = precessing_orbit_sim( t,R ) 
%%% skipped sv 1 and 2
mu = earth.mukm; 
rmag = norm([R(1),R(2)]); 
vmag = norm([R(3),R(4)]);
sv(3)=-R(4)/vmag^2*(vmag^3*k-mu/rmag^3*(-R(4)*R(1)+R(3)*R(2)));
sv(4) =  R(3)/vmag^2*(vmag^3*k-mu/rmag^3*(-R(4)*R(1)+R(3)*R(2)));
end

%%% 在积分器文件中

[t,out] = ode45(@precessing_orbit_sim,tspan,IC,options);

标签: matlab

解决方案


包含k在您的 ODE 函数中?就像是:

function [ sv ] = precessing_orbit_sim( t,R ) 

  %%% skipped sv 1 and 2
  mu = earth.mukm; 
  rmag = norm([R(1),R(2)]); 
  vmag = norm([R(3),R(4)]);

  % Add definition of a, b and h variables
  k = (h^2*(a+b)-b*h*(a+2*b)*cos(a.*t./b)+b^3)./...
      ((a+b)*(-2*b*h*cos(a.*t./b)+b^2+h^2).^2);

  sv(3)=-R(4)/vmag^2*(vmag^3*k-mu/rmag^3*(-R(4)*R(1)+R(3)*R(2)));
  sv(4) =  R(3)/vmag^2*(vmag^3*k-mu/rmag^3*(-R(4)*R(1)+R(3)*R(2)));

end

然后像往常一样调用 ODE 求解器:

[t,out] = ode45(@precessing_orbit_sim,tspan,IC,options);

推荐阅读