ode45 输出均为 NaN
ode45 outputs are all NaN
我正在尝试求解这个简单的 ODE 系统:
dydpdt = 1*(-f{2}-f{1}*dydp);
存储在一个名为 funsensitivity
的函数中(f
是一个元胞数组,具有 765x765 稀疏矩阵和 765x1 向量,可下载为 .MAT 文件 here )。我称它为:
dydp0 = zeros(size(f{2}));
[t2,JJ]=ode45(@(t,y)funsensitivity(t,y,f),0:4000:100000,dydp0);
JJ 大小合适,我没有发现任何错误,但 JJ 中的所有值都是 NaN
。我不知道为什么会这样。我做错了什么?
您正在使用 {}
进行索引,这是针对元胞数组的。你应该使用 ()
:
dydpdt = 1*(-f(2)-f(1)*dydp);
还有
dydp0 = zeros(size(f(2)));
[t2,JJ]=ode45(@(t,y)funsensitivity(t,y,f),0:tstep:tfinal,dydp0);
根据评论进行编辑:
我已经在 Octave 中尝试了 运行 你的代码(没有 MATLAB),通过查看解决方案看起来你的问题不稳定:
JJ =
Columns 1 through 10:
0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000
5.1645e+001 1.0181e+004 1.2727e+003 -1.1492e+004 -1.2900e+001 -7.2862e+001 7.2502e+001 7.7228e-003 1.1269e-002 1.4324e-003
1.5631e+031 1.3173e+033 1.6466e+032 -1.4936e+033 -3.7860e+030 -2.2204e+031 2.2012e+031 7.4674e+026 1.0897e+027 1.3851e+026
7.0857e+060 3.7159e+062 4.6449e+061 -4.2334e+062 -1.6685e+060 -1.0125e+061 1.0005e+061 1.9564e+056 2.8548e+056 3.6287e+055
3.4854e+090 1.3800e+092 1.7250e+091 -1.5787e+092 -7.9162e+089 -5.0174e+090 4.9378e+090 7.0373e+085 1.0269e+086 1.3053e+085
-5.3460e+120 5.8857e+121 7.3572e+120 -6.2253e+121 1.3755e+120 7.4917e+120 -7.4825e+120 3.1352e+115 4.5748e+115 5.8151e+114
-2.3670e+152 2.7260e+151 3.4075e+150 1.4565e+152 5.7833e+151 3.3559e+152 -3.3303e+152 9.3359e+145 1.3623e+146 1.7316e+145
-7.6120e+183 1.3332e+181 1.6682e+180 5.6622e+183 1.8358e+183 1.0823e+184 -1.0724e+184 3.3770e+177 4.9277e+177 6.2636e+176
-2.4360e+215 6.7970e+210 9.2304e+209 1.8189e+215 5.8022e+214 3.4726e+215 -3.4359e+215 1.4223e+209 2.0755e+209 2.6381e+208
-7.7953e+246 4.4692e+240 3.6568e+240 5.8278e+246 1.8337e+246 1.1142e+247 -1.1008e+247 6.0021e+240 8.7582e+240 1.1133e+240
-2.4947e+278 4.0992e+271 1.3587e+272 1.8673e+278 5.7955e+277 3.5749e+278 -3.5270e+278 2.5328e+272 3.6959e+272 4.6979e+271
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
您可能想尝试使用 ode15s
之类的刚性求解器,看看它是否可以改善问题,或者在您的时间向量中使用更小的时间步长,但看起来问题从根本上是错误的(数值为最少)。
我正在尝试求解这个简单的 ODE 系统:
dydpdt = 1*(-f{2}-f{1}*dydp);
存储在一个名为 funsensitivity
的函数中(f
是一个元胞数组,具有 765x765 稀疏矩阵和 765x1 向量,可下载为 .MAT 文件 here )。我称它为:
dydp0 = zeros(size(f{2}));
[t2,JJ]=ode45(@(t,y)funsensitivity(t,y,f),0:4000:100000,dydp0);
JJ 大小合适,我没有发现任何错误,但 JJ 中的所有值都是 NaN
。我不知道为什么会这样。我做错了什么?
您正在使用 {}
进行索引,这是针对元胞数组的。你应该使用 ()
:
dydpdt = 1*(-f(2)-f(1)*dydp);
还有
dydp0 = zeros(size(f(2)));
[t2,JJ]=ode45(@(t,y)funsensitivity(t,y,f),0:tstep:tfinal,dydp0);
根据评论进行编辑:
我已经在 Octave 中尝试了 运行 你的代码(没有 MATLAB),通过查看解决方案看起来你的问题不稳定:
JJ =
Columns 1 through 10:
0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000 0.0000e+000
5.1645e+001 1.0181e+004 1.2727e+003 -1.1492e+004 -1.2900e+001 -7.2862e+001 7.2502e+001 7.7228e-003 1.1269e-002 1.4324e-003
1.5631e+031 1.3173e+033 1.6466e+032 -1.4936e+033 -3.7860e+030 -2.2204e+031 2.2012e+031 7.4674e+026 1.0897e+027 1.3851e+026
7.0857e+060 3.7159e+062 4.6449e+061 -4.2334e+062 -1.6685e+060 -1.0125e+061 1.0005e+061 1.9564e+056 2.8548e+056 3.6287e+055
3.4854e+090 1.3800e+092 1.7250e+091 -1.5787e+092 -7.9162e+089 -5.0174e+090 4.9378e+090 7.0373e+085 1.0269e+086 1.3053e+085
-5.3460e+120 5.8857e+121 7.3572e+120 -6.2253e+121 1.3755e+120 7.4917e+120 -7.4825e+120 3.1352e+115 4.5748e+115 5.8151e+114
-2.3670e+152 2.7260e+151 3.4075e+150 1.4565e+152 5.7833e+151 3.3559e+152 -3.3303e+152 9.3359e+145 1.3623e+146 1.7316e+145
-7.6120e+183 1.3332e+181 1.6682e+180 5.6622e+183 1.8358e+183 1.0823e+184 -1.0724e+184 3.3770e+177 4.9277e+177 6.2636e+176
-2.4360e+215 6.7970e+210 9.2304e+209 1.8189e+215 5.8022e+214 3.4726e+215 -3.4359e+215 1.4223e+209 2.0755e+209 2.6381e+208
-7.7953e+246 4.4692e+240 3.6568e+240 5.8278e+246 1.8337e+246 1.1142e+247 -1.1008e+247 6.0021e+240 8.7582e+240 1.1133e+240
-2.4947e+278 4.0992e+271 1.3587e+272 1.8673e+278 5.7955e+277 3.5749e+278 -3.5270e+278 2.5328e+272 3.6959e+272 4.6979e+271
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
您可能想尝试使用 ode15s
之类的刚性求解器,看看它是否可以改善问题,或者在您的时间向量中使用更小的时间步长,但看起来问题从根本上是错误的(数值为最少)。