| 查看: 254 | 回复: 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 ] |
回复此楼» 猜你喜欢
静等基金结果
已经有8人回复
-
关于Filecode分析方法
已经有3人回复
-
国自然结果
已经有8人回复
-
FileCode能看出啥？
已经有8人回复
-
应该是93bebmhtak前后十一个字符比较关键
已经有22人回复
-
奇怪，两个人的filecode固定段从头到尾一模一样
已经有4人回复
-
8月时间戳变的，举个手。玩一下，释放压力
已经有13人回复
-
我的国基提前知道中了，可是同事的操作让我实在接受不了，怎么会有这样的人
已经有9人回复
-
好奇怪的filecode
已经有5人回复
-
综述论文作为代表作会不会影响评审专家的印象分？
已经有10人回复
zajage
金虫 (著名写手)
- 应助: 7 (幼儿园)
- 金币: 4648.5
- 散金: 1
- 红花: 1
- 帖子: 2078
- 在线: 160.1小时
- 虫号: 244937
- 注册: 2006-04-22
- 性别: GG
- 专业: 药剂学
2楼2009-02-17 19:21:46
100