大家帮我看看这个ODE45解微分方程程序哪里出错
M文件如下
function dx=f2(t,x);
m1=1.8491 m2=1.3950 me=0.7951
a1=0.03 wg=0.46we=0.02 Fm=0.04 Fin=0.2
if x(1)>1
f_x(1)=x(1)-1;
else if x(1)<-1
f_x(1)=x(1)+1;
elsef_x(1)=0;
end
dx(1)=x(2);
dx(2)=-2*m1/m2*a1*x(2)-me/m2*f_x(1)+Fm*(1-0.1*sin(we*t))+Fin*wg^2*sin(wg*t);
%dx(2)=-2*me/m2*a1*x(2)-me/m2*(1+a2*sin(wg*t))*f_x(1)+Ft
%Ft=Fm*(1-0.1*sin(we*t))+Fin*wg^2*sin(wg*t)
dx=;
运行如下
clear all;
clc;
H=;x0=;
=ode45('f2',H,x0);
plot(x(:,1),x(:,2),'m-');
错误如下
??? Error: File: f2.m Line: 2 Column: 11
Unexpected MATLAB expression.
Error in ==> funfun\private\odearguments at 80
if (nargin(ode) == 2)
Error in ==> ode45 at 173
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in ==> Untitled6 at 4
=ode45('f2',H,x0);
着急 自己弄两天了 还是弄不出来 请高手帮忙 正好最近也在解微分方程，试着弄了一下。其实没多大问题，主要是matlab语句的写法上。
1、function 定义时，不需要分号；
2、参数赋值时，要么之间加分号，要么另起一行；
3、条件语句的写法问题，好像是：if...elseif...else...end，elseif连写
最后，加上运行的结果：
1stOpt求解：
Constant m1=1.8491,m2=1.3950,me=0.7951,
a1=0.03,wg=0.46,we=0.02,Fm=0.04,Fin=0.2;
ConstStr f_x=if(x1>1,x1-1,if(x1<-1, x1+1,0));
Variable t=,x1=0,x2=0;
Plot x1,x2;
ODEFunction x1'=x2;
x2'=-2*m1/m2*a1*x2-me/m2*f_x+Fm*(1-0.1*sin(we*t))+Fin*wg^2*sin(wg*t);
%%%下面是微分方程文件
function dx=f2(t,x)
m1=1.8491;
m2=1.3950;
me=0.7951;
a1=0.03;
wg=0.46;
we=0.02;
Fm=0.04;
Fin=0.2;
if x(1)>1
f_x(1)=x(1)-1;
elseif x(1)<-1
f_x(1)=x(1)+1;
else
f_x(1)=0;
end
dx(1)=x(2);
dx(2)=-2*m1/m2*a1*x(2)-me/m2*f_x(1)+Fm*(1-0.1*sin(we*t))+Fin*wg^2*sin(wg*t);
%dx(2)=-2*me/m2*a1*x(2)-me/m2*(1+a2*sin(wg*t))*f_x(1)+Ft
%Ft=Fm*(1-0.1*sin(we*t))+Fin*wg^2*sin(wg*t)
dx=;
%%下面是主函数
clear all;
clc;
H=0.1:0.1:100;
x0=;
=ode45('f2',H,x0);
plot(x(:,1),x(:,2),'m-');
顺便问一个技术性问题：论坛里的程序怎么贴上去的？可以显示行数，又可以直接点复制。
jixiantianzi 发表于 2012-10-12 22:30 static/image/common/back.gif
%%%下面是微分方程文件
function dx=f2(t,x)
m1=1.8491;
我也不知道 同求 的确我是因为 else 和 if 没有连在一起写 而程序出错 jixiantianzi 发表于 2012-10-12 22:30 static/image/common/back.gif
%%%下面是微分方程文件
function dx=f2(t,x)
m1=1.8491;
编辑区选择插入代码文字就行，如图红圈标识位置
happy 发表于 2012-10-18 11:38 static/image/common/back.gif
编辑区选择插入代码文字就行，如图红圈标识位置
学到了，哈哈！
试验一把
%%%下面是微分方程文件
function dx=f2(t,x)
m1=1.8491;
m2=1.3950;
me=0.7951;
a1=0.03;
wg=0.46;
we=0.02;
Fm=0.04;
Fin=0.2;
if x(1)>1
f_x(1)=x(1)-1;
elseif x(1)<-1
f_x(1)=x(1)+1;
else
f_x(1)=0;
end
dx(1)=x(2);
dx(2)=-2*m1/m2*a1*x(2)-me/m2*f_x(1)+Fm*(1-0.1*sin(we*t))+Fin*wg^2*sin(wg*t);
%dx(2)=-2*me/m2*a1*x(2)-me/m2*(1+a2*sin(wg*t))*f_x(1)+Ft
%Ft=Fm*(1-0.1*sin(we*t))+Fin*wg^2*sin(wg*t)
dx=;
%%下面是主函数
clear all;
clc;
H=0.1:0.1:100;
x0=;
=ode45('f2',H,x0);
plot(x(:,1),x(:,2),'m-');
哈哈 我也来试验一把 成功了 a=1;
b=2;
页:
[1]