getting NaN in ODE45

Question

Hossein 2014-10-18

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/159161-getting-nan-in-ode45

编辑： Star Strider 2014-10-18

i have a set of nonlinear odes to solve but im getting NaN in output matrix r

i dont know what is wrong with my code, it should work,

ill be more than happy if you share your knowledge on such issue if you've faced before

cheers

here is my function and solver

function dr=Thermal_Bubble2(t,r)

dr=zeros(3,1);

global pa f r0 T0; a=0.150e-6;

l=0.0184;

kb=1.38e-23;

ps=1e5;

c=1500;

d=998;

s=0.0725;

cv=1460;

u=0.000001;

pv=2.33e3;

w=2*pi*f;

b=5.1e-29;

%Equations n=4*pi*(r0^3)*(ps+2*(s/r0)+pv)/(3*kb*T0+(ps+2*(s/r0)+pv)*b) ;

lth=min(r(1)/pi,sqrt(a*r(1)/abs(r(2))));

A=4*pi*(r(1)^3)-3*n*b;

pg=(3*n*kb*r(3))/A;

dr(1)=r(2);

dr(2)=(((1+r(1)/c)/d)*(pg-(2*s/r(1))-(4*u*r(2)/r(1))-ps+pv-pa*sin(w*t))+(r(1)/(d*c))*(((3*n*kb*((4*pi*((r(1))^2)/cv)*((l*(T0-r(3))/lth)-pg*r(2))))/A)-(36*pi*(r(1)^2)*r(2)*n*kb*(r(3))/(A^2))-(2*s*r(2)/(r(1)^2))-(4*u*(((r(2))/(r(2)))^2))-w*pa*cos(w*t))-(1.5*((r(2))^2)*(1-(r(2))/(3*c))))/((r(1))*(1-(r(2))/c));

dr(3)=(4*pi*((r(1))^2)/cv)*((l*(T0-r(3))/lth)-pg*r(2));

solver

clear;

global f pa r0 T0

f=25e3;

pa=50e3;

r0=10e-6;

T0=300;

[t,r]=ode45(@Thermal_Bubble2,(0:0.01/f:20/f),[r0 0 T0]);

plotmatrix(t,r);

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Star Strider 2014-10-18

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/159161-getting-nan-in-ode45#answer_155693

编辑：Star Strider 2014-10-18

在 MATLAB Online 中打开

I believe the problem is here:

dr(2)=(((1+r(1)/c)/d)*(pg-(2*s/r(1))-(4*u*r(2)/r(1))-ps+pv-pa*sin(w*t))+(r(1)/(d*c))*(((3*n*kb*((4*pi*((r(1))^2)/cv)*((l*(T0-r(3))/lth)-pg*r(2))))/A)-(36*pi*(r(1)^2)*r(2)*n*kb*(r(3))/(A^2))-(2*s*r(2)/(r(1)^2))-(4*u*(((r(2))/(r(2)))^2))-w*pa*cos(w*t))-(1.5*((r(2))^2)*(1-(r(2))/(3*c))))/((r(1))*(1-(r(2))/c));
NOTE: r(2)/r(2) in the expression ‘(4*u*(((r(2))/(r(2)))^2))-w*pa*cos(w*t))’

Since ‘r(2)’ is initially zero, this will cause the entire expression to be NaN, which will then propagate throughout your calculations.