声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 9244|回复: 10

[编程技巧] 大神求教这个警告是什么原因造成的,如何修改程序?

[复制链接]
发表于 2014-12-27 13:24 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
大神求教这个警告是什么原因造成的,如何修改程序?
用ode45求解系统方程
for j=1:L
j
[t ,Y]=ode45(@weifenfangcheng,[0:(2*pi/200):(288*2*pi),[0.01 0.001 0.01 0.001 0.01 0.001 0.01 0.001]);
.....
当j=2时出现警告:“在t=1.007596e-2处失败。在时间t处,若不将步长降至允许的最小值(2.775558e-17)一下,积分公差要求无法满足。”
in ode45 at 308
下标索引必须为正整数类型或逻辑类型。

另外,我是在求解方程里面加了一项函数,在不添加前方程可以运行得到结果,由于多考虑其他因素添加了一项后就出现了错误。请问出现这种错误的原因是什么啊,怎么解决啊?

回复
分享到:

使用道具 举报

 楼主| 发表于 2014-12-27 13:30 来自手机 | 显示全部楼层
有没有知道,指点一下
发表于 2014-12-29 10:07 | 显示全部楼层
好奇问下, 报错讯息是中文?

另提醒下
求助完整格式:出错代码和出错提示
 楼主| 发表于 2014-12-29 10:39 | 显示全部楼层
ChaChing 发表于 2014-12-29 10:07
好奇问下, 报错讯息是中文?

另提醒下

用的是最新的中文版matlab2014b,报错提示是中文的。
 楼主| 发表于 2014-12-29 11:07 | 显示全部楼层
本帖最后由 牛小贱 于 2015-1-4 09:34 编辑
ChaChing 发表于 2014-12-29 10:07
好奇问下, 报错讯息是中文?

另提醒下

程序如下,主要问题就是我在程序里面未添加Fa这一项之前是可以进行计算得到图,添加Fa之后就出现警告,终止运行。
  1. function fangchengfenchatu
  2. tic
  3. N=256;
  4. M=200;
  5. RelTol = 1e-6;                %相对误差Relative tolerance
  6. AbsTol = 1e-6;                %绝对误差Absolute tolerance
  7. options = odeset('RelTol',RelTol,'AbsTol',ones(1,8)*AbsTol);  %需要改变的参数‘ones(1,4)’
  8. ommiga=200:10:2500;
  9. L=length(ommiga);
  10. for j=1:L
  11. j
  12.   [t,Y] = ode45( @fangcheng001,[0:(2*pi/M):(288*2*pi)],[0.01 0.001 0.01 0.001 0.01 0.001 0.01 0.001],options,ommiga(j));
  13.   y1=Y(:,1);
  14. y11=y1((end-N*M):end);

  15. for i=1:N
  16.   
  17.    Point1(i,j)=y11(1+(i-1)*M);
  18. end
  19.   
  20. end

  21. ommiga = ommiga';

  22. plot(ommiga,Point1,'.k','MarkerSize',1)

  23. toc
复制代码
  1. function dy=fangcheng001(t,y,ommiga)  
  2. %ommiga=450;%转速
  3. m1=4;
  4. m2=32.1;
  5. c1=1050;
  6. c2=2100;
  7. c=0.00011;
  8. R=0.025;
  9. L=0.012;
  10. e=0.00005;
  11. miu=0.018;
  12. k=2.5*10^7;
  13. s=miu*ommiga*R*L*(R/c)^2*(L/2/R)^2;
  14. M=s/(m1*ommiga^2*c);
  15. %求解Fa;
  16. xi=1.2e-3;
  17. Rt=0.5;
  18. Rb=0.37;
  19. B1=35/360*2*pi;
  20. B2=40/360*2*pi;
  21. xii=0.83;
  22. V=30;
  23. p0=11.8;
  24. C=V^2*sin(B1)*p0*(cos(B1)+xii*B2);
  25. A1=(Rt^2-Rb^2)^2*pi*C*Rt/(Rt^2-Rb^2+2*Rt*xi)^2;
  26. A3=3*(Rt^2-Rb^2)^2*pi*C*Rt^3/(Rt^2-Rb^2+2*Rt*xi)^4;
  27. E=sqrt(y(5)^2+y(7)^2)/xi;
  28. Fa=A1*E+A3^3*E^3;
  29. %无量纲方程

  30. dy=[ y(2);
  31.     -c1/(m1*ommiga)*y(2)-k/(m1*ommiga^2)*(y(1)-y(5))+fx(y(1),y(3),y(2),y(4))*M;
  32.     y(4);
  33.     -c1/(m1*ommiga)*y(4)-k/(m1*ommiga^2)*(y(3)-y(7))-9.8/(c*ommiga^2)+fy(y(1),y(3),y(2),y(4))*M;
  34.     y(6);
  35.     -c2/(m2*ommiga)*y(6)-2*k/(m2*ommiga^2)*(y(5)-y(1))+e*cos(t)/c+Fa*cos(t)/(c*m2*ommiga^2);
  36.     y(8);
  37.     -c2/(m2*ommiga)*y(8)-2*k/(m2*ommiga^2)*(y(7)-y(3))+e*sin(t)/c-9.8/(c*ommiga^2)+Fa*sin(t)/(c*m2*ommiga^2);];
  38. end
复制代码








发表于 2014-12-30 09:01 | 显示全部楼层
Rotor2014 发表于 2014-12-29 11:07
程序如下,主要问题就是我在程序里面未添加Fa这一项之前是可以进行计算得到图,添加Fa之后就出现警告,终 ...

fx,fy表示什么
发表于 2014-12-30 11:14 | 显示全部楼层
新人报到,向前辈致敬!
 楼主| 发表于 2014-12-30 11:32 | 显示全部楼层
本帖最后由 牛小贱 于 2015-1-4 09:34 编辑

忘记写了  不好意思
  1. %油膜力程序
  2. function out=fx(y1,y3,y2,y4)
  3. alpha=atan((y3+2*y2)/(y1-2*y4))-pi/2*sign((y3+2*y2)/(y1-2*y4))-pi/2*sign(y3+2*y2);
  4. G=2*(pi/2+atan((y3*cos(alpha)-y1*sin(alpha))/(sqrt(1-y1^2-y3^2))))/(sqrt(1-y1^2-y3^2));
  5. V=(2+(y3*cos(alpha)-y1*sin(alpha))*G)/(1-y1^2-y3^2);
  6. S=(y1*cos(alpha)+y3*sin(alpha))/(1-(y1*cos(alpha)+y3*sin(alpha))^2);
  7. K=sqrt((y1-2*y4)^2+(y3+2*y2)^2)/(1-y1^2-y3^2);
  8. out=K*(3*y1*V-sin(alpha)*G-2*cos(alpha)*S);
  9. end
  10. function out=fy(y1,y3,y2,y4)
  11. alpha=atan((y3+2*y2)/(y1-2*y4))-pi/2*sign((y3+2*y2)/(y1-2*y4))-pi/2*sign(y3+2*y2);
  12. G=2*(pi/2+atan((y3*cos(alpha)-y1*sin(alpha))/(sqrt(1-y1^2-y3^2))))/(sqrt(1-y1^2-y3^2));
  13. V=(2+(y3*cos(alpha)-y1*sin(alpha))*G)/(1-y1^2-y3^2);
  14. S=(y1*cos(alpha)+y3*sin(alpha))/(1-(y1*cos(alpha)+y3*sin(alpha))^2);
  15. K=sqrt((y1-2*y4)^2+(y3+2*y2)^2)/(1-y1^2-y3^2);
  16. out=K*(3*y3*V+cos(alpha)*G-2*sin(alpha)*S);
  17. end
复制代码


发表于 2014-12-30 12:45 | 显示全部楼层
y1=Y(:,1);  
y11=y1((end-N*M):end);
y1只有一个值,所以y11=y1((end-N*M):end)相当于y11=y1((1-256*200):1)=y1(-51199:1),
y1下标索引只能是正值或者逻辑数值

评分

1

查看全部评分

 楼主| 发表于 2014-12-30 14:45 | 显示全部楼层
chybeyond 发表于 2014-12-30 12:45
y1=Y(:,1);  
y11=y1((end-N*M):end);
y1只有一个值,所以y11=y1((end-N*M):end)相当于y11=y1((1-256*200 ...

非常感谢您的指导,我没有添加Fa之前是可以得到结果的。还有就是提示“当j=2时出现警告:“在t=1.007596e-2处失败。在时间t处,若不将步长降至允许的最小值(2.775558e-17)一下,积分公差要求无法满足。”应该怎么处理啊?
发表于 2017-11-20 12:00 | 显示全部楼层
请问楼主解决了没有?
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-11-28 11:18 , Processed in 0.080330 second(s), 19 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表