| 查看: 557 | 回复: 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
已经有7人回复
-
大家散了吧，后缀研究没有意义，别浪费时间了，过好目前的每一天，不要焦虑
已经有4人回复
-
filecode
已经有5人回复
-
求各位大神看下
已经有31人回复
-
咨询面上基金
已经有7人回复
-
影响面上的因素
已经有11人回复
-
8月时间戳变的，举个手。玩一下，释放压力
已经有11人回复
-
一个有机合成实验室都需要哪些设备？
已经有11人回复
-
面上再次挂了，太难了，躺也躺不了，倦也卷不过，小学校之殇！
已经有30人回复
-
好消息？这个有何含义？？？
已经有10人回复
文刀星
木虫 (小有名气)
- 应助: 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
100