| 查看: 559 | 回复: 1 | |||
[求助]
程序问题 已有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]) 得出的结果只显示坐标轴不显示图像，麻烦哪位大神给指导指导，在这里谢过！ |
回复此楼» 猜你喜欢
奇怪，两个人的filecode固定段从头到尾一模一样
已经有9人回复
-
是这周出结果还是下周出结果？
已经有3人回复
-
应该是93bebmhtak前后十一个字符比较关键
已经有27人回复
-
小木虫上这么多卖论文的，真有人买论文么？感觉没必要啊
已经有9人回复
-
有时候，自然基金真的不能太认真 （我的申报经验）
已经有4人回复
-
帮忙看看fileCode
已经有13人回复
-
为什么网上很多人说本周 12号出结果
已经有7人回复
-
读博的好处
已经有3人回复
-
什么时候出结果，有咨询渠道？？？
已经有3人回复
-
分享一下我之前已中青C的计划书的filecode
已经有4人回复
文刀星
木虫 (小有名气)
- 应助: 13 (小学生)
- 金币: 3651.2
- 散金: 500
- 红花: 5
- 帖子: 60
- 在线: 100.9小时
- 虫号: 3693767
- 注册: 2015-02-22
- 性别: GG
- 专业: 导航、制导与传感技术
【答案】应助回帖
感谢参与，应助指数 +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都是ＮＡＮ，当然画不出图来啦！出现NANd的的情况我已经在上面的代码中注明了。还是好好检查检查的的微分方程组是否有错。把初值改为非零，但是也还积不出来。感觉可能是你的微分方程有问题！ |
2楼2015-12-16 14:09:03
