24小时热门版块排行榜    

查看: 560  |  回复: 1

xuanshan

新虫 (初入文坛)

[求助] 程序问题 已有1人参与

function dy=rigid4(t,y)
w=78;%轮胎的转速
R=0.317;%轮胎半径
u=10;%外界风速
ls=1.225;%空气密度
lg=2.5*10^3;%磨粒密度
d=10^(-5);%磨粒直径
g=9.806;%重力加速度
m=1.79*10^(-5);%空气动力粘度
n=m/ls;%运动粘度
k=3*ls/(4*lg*d);
g1=(1-ls/lg)*g;
dy=zeros(4,1);%一个列向量
V=sqrt(y(2)^2+y(4)^2);
u0=w*R^2/(R+y(1));
U1=(0.32*y(1)*u^1.5*n*y(3)-0.005*y(1)^3*u^2.5)/(n^1.5*y(3)^1.5)+u0;
U2=(0.08*n*y(1)^2*u^1.5*y(3)-0.002*y(1)^4*u^2.5)/(n^1.5*y(3)^2.5);
U=sqrt(U1^2+U2^2);
Re=ls*d*abs(U-V)/m;
Cd=19.5/Re^0.633;
dy(1)=y(2);
dy(3)=y(4);
dy(2)=k*Cd*abs(U2-y(2))*(U2-y(2));
dy(4)=-g1+k*Cd*abs(U1-y(4))*(U1-y(4));
dy=[dy(1) dy(2) dy(3) dy(4)]';
end
上面是主程序,下面是命令窗口中的程序
tspan=linspace(0,1,1000);
[t,y]=ode45(@rigid4,tspan,[0 2.3 0 2.1]);
plot(y(:,1),y(:,3))
axis([0 0.003 0 1])
得出的结果只显示坐标轴不显示图像,麻烦哪位大神给指导指导,在这里谢过!
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

文刀星

木虫 (小有名气)

【答案】应助回帖

感谢参与,应助指数 +1
function dy=rigid4(t,y)
w=78;%轮胎的转速
R=0.317;%轮胎半径
u=10;%外界风速
ls=1.225;%空气密度
lg=2.5*10^3;%磨粒密度
d=10^(-5);%磨粒直径
g=9.806;%重力加速度
m=1.79*10^(-5);%空气动力粘度
n=m/ls;%运动粘度
k=3*ls/(4*lg*d);
g1=(1-ls/lg)*g;
dy=zeros(4,1);%一个列向量
V=sqrt(y(2)^2+y(4)^2);
u0=w*R^2/(R+y(1));
U1=(0.32*y(1)*u^1.5*n*y(3)-0.005*y(1)^3*u^2.5)/(n^1.5*y(3)^1.5)+u0;%%%% U1=NAN 有问题,问题在于y(3)=0,做除数了
U2=(0.08*n*y(1)^2*u^1.5*y(3)-0.002*y(1)^4*u^2.5)/(n^1.5*y(3)^2.5);%%% 同样,U2=NAN,也由于y(3)=0,在分母上
U=sqrt(U1^2+U2^2); % NAN
Re=ls*d*abs(U-V)/m;% NAN
Cd=19.5/Re^0.633;% NAN
dy(1)=y(2);
dy(3)=y(4);
dy(2)=k*Cd*abs(U2-y(2))*(U2-y(2)); %NAN
dy(4)=-g1+k*Cd*abs(U1-y(4))*(U1-y(4));% NAN

end


你求出的数值解y都是NAN,当然画不出图来啦!出现NANd的的情况我已经在上面的代码中注明了。还是好好检查检查的的微分方程组是否有错。把初值改为非零,但是也还积不出来。感觉可能是你的微分方程有问题!
2楼2015-12-16 14:09:03
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 xuanshan 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 欢迎发来filecode的Mz6后的代码验证其规律 +24 医学老男孩 2026-08-13 50/2500 2026-08-15 06:38 by 学员8dgXkO
[公派出国] 售SCI一区T0P文章,我:8.O55.1.O.54,科目全,可十急 +3 7lpolszZVXgi 2026-08-14 5/250 2026-08-15 04:41 by 4wMiSEwB6436
[找工作] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +4 k0dTPqJtl0jt 2026-08-14 4/200 2026-08-15 02:21 by 4wMiSEwB6436
[硕博家园] 售SCI一区文章,我:8O5.5.1.O5.4,科目全,可伽急 +3 k0dTPqJtl0jt 2026-08-14 4/200 2026-08-15 01:28 by 4wMiSEwB6436
[考研] 售SCI一区T0P文章,我:8.O.55.1.O.54,科目齐全,可+急 +4 HFw0lei2R37i 2026-08-14 8/400 2026-08-15 00:16 by 4wMiSEwB6436
[基金申请] 是这周出结果还是下周出结果? +4 yuleib84 2026-08-11 4/200 2026-08-14 23:05 by lfy8008
[硕博家园] 请教兼职经验 +3 是阿文鸭 2026-08-09 3/150 2026-08-14 12:07 by HER12025
[基金申请] filecode +15 documentary 2026-08-10 17/850 2026-08-14 10:08 by kissu88
[基金申请] 我的国基提前知道中了,可是同事的操作让我实在接受不了,怎么会有这样的人 +10 家与远方 2026-08-10 15/750 2026-08-14 02:08 by 绵羊哥哥
[基金申请] 关于Filecode分析方法 +9 majunge000 2026-08-10 12/600 2026-08-13 23:20 by iwuli
[文学芳草园] 阿姨 +4 汪汪锅 2026-08-09 4/200 2026-08-13 19:43 by arzu_hma
[硕博家园] 一作与独作在应聘高校教师时区别大吗 +3 mbygzh 2026-08-08 4/200 2026-08-13 19:31 by 龙-樱
[基金申请] 静等基金结果 +8 gjjjzhong 2026-08-10 21/1050 2026-08-13 17:56 by 且听虎啸
[基金申请] 重要来源:本周末出结果 +10 瞬息宇宙 2026-08-12 10/500 2026-08-13 15:46 by likettle
[基金申请] Filecode 又变了,巨变 +3 WH3796 2026-08-12 4/200 2026-08-13 14:13 by 小木虫6752397
[基金申请] 结合人工智能,周易传统文化,filecode打分制来了,3分以上希望很大。 +3 Tide man 2026-08-12 4/200 2026-08-13 08:35 by ZJTJZ
[基金申请] 帮忙看看fileCode +7 wwncly 2026-08-10 13/650 2026-08-11 19:36 by 冰心玉壶晴
[基金申请] 国自然结果 +4 Vierhys 2026-08-10 8/400 2026-08-10 15:06 by Vierhys
[基金申请] 据悉今年马上要出结果了 +7 瞬息宇宙 2026-08-10 8/400 2026-08-10 12:42 by Vivilian
[基金申请] 这样的filecode谁见过 +11 布布和一二 2026-08-08 22/1100 2026-08-10 11:10 by wmfsnow
信息提示
请填处理意见