马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
本帖最后由 komazsl 于 2011-5-22 22:01 编辑
第一个M文件
function dy=rigid(t,y,flag,Ce,Ta,Tb,Tc,Td,Te,Tf,Tg,Th,Ti,Tj,Tk,Tl,Tm,Tn,To,Tp,Tq,Tr);
dy=zeros(9,1);
dy(1)=-Tb*y(1)*y(4)-Tc*y(3)*y(1)+Td*y(3)*y(4)+Tg*y(3)*y(5)-Th*y(1)*y(4)-Ti*y(1)*y(6)-Tj*y(2)*y(2)*y(5)+2*Tk*y(3)*y(7)+Tl*y(7)*y(4)-Tp*y(8)*y(1);
dy(2)=Tc*y(3)*y(1);
dy(3)=Tb*y(1)*y(4)-Tc*y(3)*y(1)-Td*y(3)*y(4)-Tg*y(3)*y(5)-Tk*y(3)*y(7);
dy(4)=2*Ta*y(5)*Ce-Tb*y(1)*y(4)+Tc*y(3)*y(1)-Td*y(3)*y(4)-Te*y(4)*y(5)-Tf*y(4)*y(6)+Tg*y(3)*y(5)-Th*y(1)*y(4)-Tl*y(7)*y(4)-Tm*y(7)*y(4)-To*y(8)*y(4);
dy(5)=-Ta*y(5)*Ce+Tb*y(1)*y(4)-Te*y(4)*y(5)+2*Tf*y(4)*y(6)-Tg*y(3)*y(5)+Ti*y(1)*y(6)-Tj*y(2)*y(2)*y(5)+Tl*y(7)*y(4)+Tn*y(7)*y(6)+To*y(8)*y(4);
dy(6)=Te*y(4)*y(5)-Tf*y(4)*y(6)-Ti*y(1)*y(6)-Tn*y(7)*y(6);
dy(7)=Th*y(1)*y(4)+Ti*y(1)*y(6)-Tk*y(3)*y(7)-Tl*y(7)*y(4)-Tm*y(7)*y(4)-Tn*y(7)*y(6)+To*y(8)*y(4)+2*Tp*y(8)*y(1)-Tq*y(8)*y(7)+Tr*y(9);
dy(8)=Tm*y(7)*y(4)+Tn*y(7)*y(6)-To*y(8)*y(4)-Tp*y(8)*y(1)-Tq*y(8)*y(7)+Tr*y(9);
dy(9)=Tq*y(8)*y(7)-Tr*y(9);
end
第二个
tspan=[0,6];
y0=[0.001,0.999,0,0,0,0,0,0,0];
Ce=1.00E+13;
Ta=1.29E+15;
Tb=2.77E+9;
Tc=2.04E+13;
Td=1.45E+11;
Te=1.69E+12;
Tf=4.81E+9;
Tg=5.36E+7;
Th=8.43E+11;
Ti=1.09E+10;
Tj=7.09E+9;
Tk=1.81E+12;
Tl=5.84E+12;
Tm=4.38E+7;
Tn=1.97E+9;
To=1.02E+13;
Tp=1.20E+13;
Tq=4.57E+12;
Tr=1.32E+11;
options=odeset('RelTol',1e-6);
t0=cputime;
[t,y]=ode23s('rigid',tspan,y0,options,Ce,Ta,Tb,Tc,Td,Te,Tf,Tg,Th,Ti,Tj,Tk,Tl,Tm,Tn,To,Tp,Tq,Tr);
time=cputime-t0;
figure;
plot(t,y(:,1),'*');
hold on;
plot(t,y(:,2),'o');
set(gca,'Fontsize',12);
ylim([0,0.006]);
xlabel('\itt','Fontsize',16);
set(gcf,'color','w');
算出来出现
Warning: Matrix is close to singular or badly scaled.
Results may be inaccurate. RCOND = 7.319604e-029.
> In ode23s at 456
In Untitled11 at 24,
然后能出现图但是y值还是初值,根本没有变化,谢谢各位帮我看看,我毕设遇到难题了
|