求解一个简单的微分方程组

Question

如何使用八度数值求解以下简单的微分方程组？

注意：

• 我使用限定词“简单”，因为据我了解，该系统是一阶且未耦合。
• 我已经尽力在线尝试解决此问题的方法和脚本，包括here，here和here。无论如何，我要么被吊死，无响应的八度，提示“重复收敛失败”，建议我手动调整的错误初始和最大步长（我确实尝试过，但没有），或者由于没有错误而最初看起来像是解决方案的东西，但绘制解决方案会显示空白图表
• [Octave为等效的Matlab例程提供的位置，我尝试了各种例程ode45，ode23，ode113，ode15s，ode23s，ode23t，ode23tb，ode15i，当然，八度音阶的lsode命令，都给出与上述相同的错误。

Answer 1

让我们首先复制香草溶液

% z = [x,y]
f = @(t,z) [ z(1).^2+t; z(1).*z(2)-2 ];
z0 = [ 2; 1];

[ T, Z ] = ode45(f, [0, 10], z0);

plot(T,Z); legend(["x";"y"]);

积分器失败，并显示警告消息

警告：解决方案未成功。迭代积分循环在到达t = 0.494898的端点之前的时间tend = 10.000000退出。如果步长变得太小，可能会发生这种情况。尝试使用命令'InitialStep'减小'MaxStep'和/或'odeset'的值。

直到关键时刻之前重复集成

opt = odeset('MaxStep',0.01);
[ T, Z ] = ode45(f, [0, 0.49], z0, opt);
clf; plot(T,Z); legend(["x";"y"]);

图形结果

可以看到第一个方程式中的二次项导致失控增长。由于某些原因，求解器只能识别不断减小的步长，而不能识别解决方案的失控值。

实际上是第一个Riccati方程，已知它在有限的时间具有极点。使用典型的参数化x(t)=-u'(t)/u(t)给出ODE

u''(t)+t*u(t)=0

具有t>0的振荡分支的Airy方程。 u的第一个根是x的极点，无法将解扩展到此点之外。

g=@(t,u) [u(2); -t.*u(1)]
u0 = [ 1; -2];

function [val,term, dir] = event(t,u)
   val = u(1);
   term = 0;
   dir = 0;
end
opt = odeset('MaxStep',0.1, 'Events', @(t,u) event(t,u));
[T,U,Te,Ue,Ie] = ode45(g,[0,4],u0,opt);
disp(Te)
clf; plot(T,U); legend(["u";"u'"]);

[将u的零列为0.4949319379979706, 2.886092605590324，再次确认警告的原因，并给出曲线图