| 查看: 259 | 回复: 1 | |||
| 当前主题已经存档。 | |||
chshy1234银虫 (小有名气)
|
[交流]
【求助】，MATLAB程序错在哪了？？？
|
||
|
Error in ==> ode45 at 173 [neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ... Error in ==> chemical_kinetic at 17 [t,x]=ode45(@myfun,[t0:0.1:t1],x0,optimset('display','iter'),A); 运行下面程序出现上述错误：各位高手帮忙分析一下错在哪了？？？ function chemical_kinetic %计算反应物随浓度的变化关系 clear all clc %输入初值 A0=1.50; B0=0;C0=0;D0=0; k1=1.00;k2=0.50;k3=0.30;k4=0.10; x0=[B0 C0 D0]; %赋浓度初值 t0=0;t1=20; %赋时间起始值 for f=0.5:0.1:1.0 switch f case 0.5 %计算反应速率方程中的系数矩阵 A=[-k1*f*A0 k2 -k3*f*A0 0;k1*f*A0 -k2 -k3*f*A0 0 ;0 0 k3*f*A0 -k4;0 0 0 k4]; %调用ode45函数计算常微分方程组 [t,x]=ode45(@myfun,[t0:0.1:t1],x0,optimset('display','iter'),A); %作反应物浓度随时间的变化图 figure(1) plot(t,x(:,1),'k-') xlabel('t') ylabel('B浓度') title('B的浓度变化曲线(f=0.5)') figure(2) plot(t,x(:,2),'k-') xlabel('t') ylabel('C浓度') title('C的浓度变化曲线(f=0.5)') figure(3) plot(t,x(:,3),'k-') xlabel('t') ylabel('D浓度') title('D的浓度变化曲线(f=0.5)') fprintf('\n f=0.5\n') fprintf('\n B=%.4g\n',x(length(x),1)) fprintf('\n C=%.4g\n',x(length(x),2)) fprintf('\n D=%.4g\n',x(length(x),3)) case 0.75 A=[-k1*f*A0 k2 -k3*f*A0 0;k1*f*A0 -k2 -k3*f*A0 0 ;0 0 k3*f*A0 -k4;0 0 0 k4]; [t,x]=ode45(@myfun,[t0:0.1:t1],x0,optimset('display','iter'),A); figure(4) plot(t,x(:,1),'k-') xlabel('t') ylabel('B浓度') title('B的浓度变化曲线(f=1.6)') figure(5) plot(t,x(:,2),'k-') xlabel('t') ylabel('C浓度') title('C的浓度变化曲线(f=1.6)') figure(6) plot(t,x(:,3),'k-') xlabel('t') ylabel('D浓度') title('D的浓度变化曲线(f=1.6)') fprintf('\n f=1.6\n') fprintf('\n B=%.4g\n',x(length(x),1)) fprintf('\n C=%.4g\n',x(length(x),2)) fprintf('\n D=%.4g\n',x(length(x),3)) case 1.0 A=[-k1*f*A0 k2 -k3*f*A0 0;k1*f*A0 -k2 -k3*f*A0 0 ;0 0 k3*f*A0 -k4;0 0 0 k4]; [t,x]=ode45(@myfun,[t0:0.1:t1],x0,optimset('display','iter'),A); figure(7) plot(t,x(:,1),'k-') xlabel('t') ylabel('B浓度') title('B的浓度变化曲线(f=2.3)') figure(8) plot(t,x(:,2),'k-') xlabel('t') ylabel('C浓度') title('C的浓度变化曲线(f=2.3)') figure(9) plot(t,x(:,3),'k-') xlabel('t') ylabel('D浓度') title('D的浓度变化曲线(f=2.3)') fprintf('\n f=2.3\n') fprintf('\n B=%.4g\n',x(length(x),1)) fprintf('\n C=%.4g\n',x(length(x),2)) fprintf('\n D=%.4g\n',x(length(x),3)) end end %------------------ function y=myfun(t,x,A) y=A*[1;x(1);x(1);x(2)]; [ Last edited by sunxiao on 2009-3-9 at 08:30 ] |
回复此楼» 猜你喜欢
梦想
已经有5人回复
-
中青基了要发朋友圈吗？
已经有7人回复
-
有多少人是今天查系统知道结果的？
已经有19人回复
-
29号明天会评吗
已经有4人回复
-
国社科又开始会评了，不知道这次命运如何
已经有10人回复
-
基金不中，共勉
已经有12人回复
-
能否申诉？
已经有7人回复
-
面上意见出来了
已经有19人回复
-
为什么到现在没收到通知？
已经有5人回复
-
有没有仍没收到信息的
已经有8人回复
zajage
金虫 (著名写手)
- 应助: 7 (幼儿园)
- 金币: 4648.5
- 散金: 1
- 红花: 1
- 帖子: 2078
- 在线: 160.1小时
- 虫号: 244937
- 注册: 2006-04-22
- 性别: GG
- 专业: 药剂学
2楼2009-02-17 19:21:46
投票:
