请各位大虾看看我的程序，错在什么地方！

lc622503

function asm21
clear all,clc
x0=[20,30,16,0,3.6,25,125,30,50,30,1,1,20,30,16,0,3.6,25,125,30,50,30,1,1]
[t,x]=ode45(@f,[0:1:15],x0)
%-----------------------
function dxdt=f(t,x)
% parameters
par(1)=3;par(2)=0.06;par(3)=0.1;par(4)=0.2;par(5)=0.5;
par(6)=0.1;par(7)=6;par(8)=3;par(9)=0.8;par(10)=0.4;
par(11)=0.2;par(12)=4;par(13)=20;par(14)=4;par(15)=0.5;
par(16)=0.05;par(17)=0.01;par(18)=0.1;par(19)=3;par(20)=1.5;
par(21)=1;par(22)=0.2;par(23)=0.2;par(24)=0.2;par(25)=0.2;
par(26)=4;par(26)=0.05;par(28)=0.2;par(29)=0.01;par(30)=0.1;
par(31)=0.01;par(32)=0.34;par(33)=0.02;par(34)=0.01;par(35)=1;
par(36)=0.15;par(37)=0.5;par(38)=1;par(39)=0.5;par(40)=0.01;
%influent concentrations
i1=20; i2=30; i3=16; i4=0; i5=3.6
i6=25; i7=125; i8=30;i9=50; i10=30
i11=1; i12=1
%tank1 concentrations
x(1); x(2); x(3); x(4); x(5)
x(6); x(7); x(8); x(9); x(10)
x(11); x(12);oxygen1=0
%tank2 concentrations
x(13); x(14); x(15);x(16);x(17)
x(18); x(19); x(20); x(21); x(22)
x(23); x(24);oxygen2=3
%recycle concentrations
r1=x(13);r2=x(14);r3=x(15);r4=x(16);r5=x(17)
r6=3*x(18);r7=3*x(19);r8=3*x(20);r9=3*x(21)
r10=3*x(22);r11=3*x(23);r12=3*x(24)
%waste sludge concentrations
w1=x(13);w2=x(14);w3=x(15);w4=x(16);w5=x(17)
w6=3*x(18);w7=3*x(19);w8=3*x(20);w9=3*x(21)
w10=3*x(22);w11=3*x(23);w12=3*x(24)
%efluent concentrations
e1=x(13);e2=x(14);e3=x(15);e4=x(16);e5=x(17)
%flow rates and volumes
qi=10;q1=20;v1=3;q2=10;v2=3
qr=10;qw=2;q0=8
% processes
proc1=par(1)*(oxygen1/(par(4)+oxygen1))*((x(7)/x(8))/(par(6)+x(7)/x(8)))*x(8)
proc2=par(1)*par(2)*(par(4)/(par(4)+oxygen1))*(x(4)/(par(5)+x(4)))*((x(7)/x(8))/(par(6)+x(7)/x(8)))*x(8)
proc3=par(1)*par(3)*(par(4)/(par(4)+oxygen1))*(par(5)/(par(5)+x(4)))*((x(7)/x(8))/(par(6)+x(7)/x(8)))*x(8)
proc4=par(7)*(oxygen1/(par(11)+oxygen1))*(x(2)/(par(12)+x(2)))*(x(2)/(x(1)+x(2)))*(x(3)/(par(16)+x(3)))...
       *(x(5)/(par(17)+x(5)))*x(8)
proc5=par(7)*(oxygen1/(par(11)+oxygen1))*(x(1)/(par(14)+x(1)))*(x(1)/(x(1)+x(2)))*(x(3)/(par(16)+x(3)))...
       *(x(5)/(par(17)+x(5)))*x(8)
proc6=par(7)*par(9)*(par(11)/(par(11)+oxygen1))*(x(2)/(par(12)+x(2)))*(x(2)/(x(1)+x(2)))*(x(3)/(par(16)+x(3)))...
       *(x(4)/(par(15+x(4))))*(x(5)/(par(17)+x(5)))*x(8)
proc7=par(7)*par(9)*(par(11)/(par(11)+oxygen1))*(x(1)/(par(14)+x(1)))*(x(1)/(x(1)+x(2)))*(x(3)/(par(16)+x(3)))...
       *(x(4)/(par(15+x(4))))*(x(5)/(par(17)+x(5)))*x(8)
proc8=par(8)*(par(11)/(par(11)+oxygen1))*(par(15)/(par(15)+x(4)))*(x(2)/(par(13)+x(2)))*x(8)
proc9=par(10)*x(8)
proc10=par(19)*(x(1)/(par(26)+x(1)))*((x(11)/x(10))/(par(31)+x(11)/x(10)))*x(10)
proc11=par(20)*(oxygen1/(par(25)+oxygen1))*(x(5)/(par(28)+x(5)))*((x(12)/x(10))/(par(34)+x(12)/x(10)))...
    *((par(32)-(x(11)/x(10)))/(par(33)+par(32)-(x(11)/x(10))))*x(10)
proc12=par(21)*(oxygen1/(par(25)+oxygen1))*(x(5)/(par(29)+x(5)))*((x(12)/x(10))/(par(34)+x(12)/x(10)))...
         *(x(3)/(par(27)+x(3)))*x(10)
proc13=par(22)*x(10)
proc14=par(23)*x(11)  
proc15=par(24)*x(12)
proc16=par(35)*(oxygen1/(par(37)+oxygen1))*(x(3)/(par(38)+x(3)))*(x(5)/(par(40)+x(5)))*x(9)
proc17=par(36)*x(9)
proc18=par(1)*(oxygen2/(par(4)+oxygen2))*((x(19)/x(20))/(par(6)+x(19)/x(20)))*x(20)
proc19=par(1)*par(2)*(par(4)/(par(4)+oxygen2))*(x(16)/(par(5)+x(16)))*((x(19)/x(20))/(par(6)+x(19)/x(20)))*x(20)
proc20=par(1)*par(3)*(par(4)/(par(4)+oxygen2))*(par(5)/(par(5)+x(16)))*((x(19)/x(20))/(par(6)+x(19)/x(20)))*x(20)
proc21=par(7)*(oxygen2/(par(11)+oxygen2))*(x(14)/(par(12)+x(14)))*(x(14)/(x(13)+x(14)))*(x(15)/(par(16)+x(15)))...
       *(x(17)/(par(17)+x(17)))*x(20)
proc22=par(7)*(oxygen2/(par(11)+oxygen2))*(x(13)/(par(14)+x(13)))*(x(13)/(x(13)+x(14)))*(x(15)/(par(16)+x(15)))...
       *(x(17)/(par(17)+x(17)))*x(20)
proc23=par(7)*par(9)*(par(21)/(par(11)+oxygen2))*(x(14)/(par(12)+x(14)))*(x(14)/(x(13)+x(14)))*(x(15)/(par(16)+x(15)))...
       *(x(16)/(par(15+x(16))))*(x(17)/(par(17)+x(17)))*x(20)
proc24=par(7)*par(9)*(par(21)/(par(11)+oxygen2))*(x(13)/(par(14)+x(13)))*(x(13)/(x(13)+x(14)))*(x(15)/(par(16)+x(15)))...
       *(x(16)/(par(15+x(16))))*(x(17)/(par(17)+x(17)))*x(20)
proc25=par(8)*(par(11)/(par(11)+oxygen2))*(par(15)/(par(15)+x(16)))*(x(14)/(par(13)+x(14)))*x(20)
proc26=par(10)*x(20)
proc27=par(19)*(x(13)/(par(26)+x(13)))*((x(23)/x(22))/(par(31)+x(23)/x(22)))*x(22)
proc28=par(20)*(oxygen2/(par(25)+oxygen2))*(x(17)/(par(28)+x(17)))*((x(24)/x(22))/(par(34)+x(24)/x(22)))...
    *((par(32)-(x(23)/x(22)))/(par(33)+par(32)-(x(23)/x(22))))*x(22)
proc29=par(21)*(oxygen2/(par(25)+oxygen2))*(x(17)/(par(29)+x(17)))*((x(24)/x(22))/(par(34)+x(24)/x(22)))...
         *(x(15)/(par(27)+x(15)))*x(22)
proc30=par(22)*x(22)
proc31=par(23)*x(23)      
proc32=par(24)*x(24)
proc33=par(35)*(oxygen2/(par(37)+oxygen2))*(x(15)/(par(38)+x(15)))*(x(17)/(par(40)+x(17)))*x(21)
proc34=par(36)*x(21)
%f
f1=(qi/v1)*i1+(qr/v1)*x(13)-(q1/v1)*x(1)-1.59*(proc5+proc7)+proc8-proc10+proc15
f2=(qi/v1)*i2+(qr/v1)*x(14)-(q1/v1)*x(2)+proc1+proc2+proc3-1.59*(proc4+proc6)-proc8
f3=(qi/v1)*i3+(qr/v1)*x(15)-(q1/v1)*x(3)+0.01*(proc1+proc2+proc3)-0.022*(proc4+proc6)-0.07*(proc5+proc7+proc12)...
    +0.03*proc8+0.031*(proc9+proc13+proc17)-4.24*proc16
f4=(qi/v1)*i4+(qr/v1)*x(16)-(q1/v1)*x(4)-0.21*(proc6+proc7)+4.17*proc16
f5=(qi/v1)*i5+(qr/v1)*x(17)-(q1/v1)*x(5)-0.004*(proc4+proc6)-0.02*(proc5+proc7+proc12+proc16)...
    +0.01*(proc8+proc9+proc13+proc17)+0.4*proc10-proc11+proc14
f6=(qi/v1)*i6+(qr/v1)*3*x(18)-(q1/v1)*x(6)+0.1*(proc9+proc13+proc17)
f7=(qi/v1)*i7+(qr/v1)*3*x(19)-(q1/v1)*x(7)-(proc1+proc2+proc3)+0.9*(proc9+proc13+proc17)
f8=(qi/v1)*i8+(qr/v1)*3*x(20)-(q1/v1)*x(8)+proc4+proc5+proc6+proc7-proc9
f9=(qi/v1)*i9+(qr/v1)*3*x(21)-(q1/v1)*x(9)+proc16-proc17
f10=(qi/v1)*i10+(qr/v1)*3*x(22)-(q1/v1)*x(10)+proc12-proc13
f11=(qi/v1)*i11+(qr/v1)*3*x(23)-(q1/v1)*x(11)-0.4*proc10+proc11-proc14
f12=(qi/v1)*i12+(qr/v1)*3*x(24)-(q1/v1)*x(12)+proc10-0.2*proc11-1.6*proc12-proc15
f13=(q1/v2)*x(1)-(q0+qr+qw/v2)*x(13)-1.59*(proc22+proc24)+proc25-proc27+proc32
f14=(q1/v2)*x(2)-(q0+qr+qw/v2)*x(14)+proc18+proc19+proc20-1.59*(proc21+proc23)-proc25
f15=(q1/v2)*x(3)-(q0+qr+qw/v2)*x(15)+0.01*(proc18+proc19+proc20)-0.022*(proc21+proc23)-0.07*(proc22+proc24+proc29)...
    +0.03*proc25+0.031*(proc26+proc30+proc34)-4.24*proc33
f16=(q1/v2)*x(4)-(q0+qr+qw/v2)*x(16)-0.21*(proc23+proc24)+4.17*proc33
f17=(q1/v2)*x(5)-(q0+qr+qw/v2)*x(17)-0.004*(proc21+proc23)-0.02*(proc22+proc24+proc29+proc33)...
    +0.01*(proc25+proc26+proc30+proc34)+0.4*proc27-proc28+proc31
f18=(q1/v2)*x(6)-3*((qr+qw)/v2)*x(18)+0.1*(proc26+proc30+proc34)
f19=(q1/v2)*x(7)-3*((qr+qw)/v2)*x(19)-(proc18+proc19+proc20)+0.9*(proc26+proc30+proc34)
f20=(q1/v2)*x(8)-3*((qr+qw)/v2)*x(20)+proc21+proc22+proc23+proc24-proc26
f21=(q1/v2)*x(9)-3*((qr+qw)/v2)*x(21)+proc33-proc34
f22=(q1/v2)*x(10)-3*((qr+qw)/v2)*x(22)+proc29-proc30
f23=(q1/v2)*x(11)-3*((qr+qw)/v2)*x(23)-0.4*proc27+proc28-proc31
f24=(q1/v2)*x(12)-3*((qr+qw)/v2)*x(24)+proc27-0.2*proc28-1.6*proc29-proc32
dxdt=[f1,f2,f3,f4,f5,f6,f7,f8,f9,f10,f11,f12,f13,f14,f15,f16,f17,f18,f19,f20,f21,f22,f23,f24]

谢谢看到这里

icepoint

建议将运行时的错误提示或警告等也写出来，以便于找出错误原因

sysh320

建议楼上说的，这样盲目检查程序很费时间的

lc622503

谢谢楼上的各位了
昨晚已经查出错误来了
dxdt=[f1,f2,f3,f4,f5,f6,f7,f8,f9,f10,f11,f12,f13,f14,f15,f16,f17,f18,f19,f20,f21,f22,f23,f24]
上面的逗号应该改成分号
这样返回值就会是列向量

sysh320

原来这样啊，学到了一点知识，谢谢楼主啊

lc622503

还有一个问题问大家
就是上面的程序其实前面12个方程和后面12个结构和系数是一样的
集体物理过程就是有两个反应器是串联的，每个里面有12个物质在发生反应
前面的物质反应完后再流入后面一个
在模拟的过程中，可能有必要再加几个反应器
但是到时程序编写就会更加麻烦，方程组会有n*12个方程组成
请问大家要解决这方面的问题要看matlab哪方面的内容
刚刚学习matlab现学现用，所以不是很了解
有这方面的例子更好可以发给我
lc622503@163.com
谢谢大家