想用if语句编写裂纹综合模型的程序，总出错，

竹语随风

ChaChing 发表于 2013-9-9 16:01

                               
登录/注册后可看大图

可能个人在4F没说清楚吧!?

基本就是13F所说的, 不过建议LZ需谨慎使用预设值(类似13F提供的)

老师，我把我要引用那个综合模型已经贴出来了，请老师帮我看看我的整个程序都出什么问题了，具体在不在我也弄不清楚了。

ChaChing

竹语随风 发表于 2013-9-9 16:12

                               
登录/注册后可看大图

老师，我就想引用这个函数不好意思，最后一个取值范围应该是3*pi/2-alfa到后面的值

LZ是否了解atan,atan2,angle等函数间的差异了!?
有无看了下帮助文献了!?
再次提醒注意下范围!

竹语随风

竹语随风 发表于 2013-9-9 16:12

                               
登录/注册后可看大图

老师，我就想引用这个函数不好意思，最后一个取值范围应该是3*pi/2-alfa到后面的值

老师，您好，我将我的程序重新改了一下。可是还是出错，你帮忙再看看。
调用程序
function dx=crackandliner(t,x,c,w)
e=0.1;
k=0.3;
beta=0;
phi0=0;
alfa=pi/6;
phi=atan(x(1)/x(3));
theta=w*t+phi0+beta-phi;
if -pi/2+alfa<=theta & theta<=pi/2-alfa
    f=1;
elseif pi/2-alfa<theta & theta<=pi/2+alfa
    f=1/2+1/2*cos(pi*(theta-pi/2+alfa)./(2*alfa));
elseif pi/2+alfa<theta & theta<=3*pi/2-alfa
    f=0;
elseif 3*pi/2-alfa<=theta & theta<=3*pi/2+alfa
    f=1/2+1/2*cos(pi*(theta-3*pi/2-alfa)./(2*alfa));   
end
n=sin(w*t);
m=cos(w*t);
dx=zeros(4,1);
dx(1)=x(2);
dx(2)=e*w^2.*cos(w*t+phi0)-c*x(2)-x(1)+f*k*m^2*x(1)...
    +f*k*n*m*x(3);
dx(3)=x(4);
dx(4)=e*w^2.*sin(w*t+phi0)-c*x(4)-x(3)+f*k*n^2*x(3)...
    +f*k*n*m*x(1);
主程序
    j=1;
    c=0.2;
    x0=[0.01;0.01;0.01;0.01];
    tend=linspace(0,2000,5000);
    for w=0.1:0.01:3
         options=odeset;options.RelTol=1e-4;
         [t,x]=ode45(@crackandliner,tend,x0,options,c,w,f);
          p= sqrt((x(end-400:end,1).^2)+(x(end-400:end,3).^2));
          Pmax(j,1)=max(p);
          j=j+1;
    end
w=0.1:0.01:3;
h=plot(w,Pmax);
title('p-w图');
xlabel('回转速度w');
ylabel('振幅p');

ChaChing

竹语随风 发表于 2013-9-10 09:27

                               
登录/注册后可看大图

老师，您好，我将我的程序重新改了一下。可是还是出错，你帮忙再看看。

重新改了一下? 移动执行的位置有差吗?
我想楼主真的没清楚个人所说的内容!?
个人一直在强调逻辑范围的不适当, 楼主好像还是没能朝此方向思考呀
由於可能专业不同, 怕造成个人误导, 一直不敢也不想贴
但好像楼主还没有结果似的!? 就提供参考吧

dx=crackandliner.m中
phi=atan(x(1)/x(3)); 改成 phi=atan2(x(1),x(3));

fcrack.m中一开始加上
theta=mod(theta+2*pi,2*pi);
if theta>=3*pi/2+alfa, theta=theta-2*pi; end

结果如下! 再强调下, 仅针对原程序出错讯息侦错, 其他未有能力侦错

飞驰吧少年

ChaChing 发表于 2013-9-13 15:38
重新改了一下? 移动执行的位置有差吗?
我想楼主真的没清楚个人所说的内容!?
个人一直在强调逻辑范围的 ...

请问一下当时您在这篇帖子里最后对theta取mod为什么要这样做呢？再判断theta和3/2pi+alfa的关系有什么意思呢？谢谢！！

ChaChing

飞驰吧少年 发表于 2014-5-28 10:42
请问一下当时您在这篇帖子里最后对theta取mod为什么要这样做呢？再判断theta和3/2pi+alfa的关系有什么意 ...

看下atan,atan2,angle等函数的差异
取mod主要控制角度范围落在LZ所想要的范围内(我猜测)
不然那个if判断压根没作用,设断点执行下很容易发现的

飞驰吧少年

ChaChing 发表于 2014-5-28 13:24
看下atan,atan2,angle等函数的差异
取mod主要控制角度范围落在LZ所想要的范围内(我猜测)
不然那个if判 ...

谢谢~!好像有点明白了`