How do I fix 'Error using odearguments @.... returns a vector of length 0, but the length of initial conditions vector is 24...' error?

20 次查看（过去 30 天）
KIPROTICH KOSGEY
KIPROTICH KOSGEY 2023-7-23
移动：Torsten 2023-8-10
Ran in:
Hi there,
I keep getting this error and I can't find where the problem is. Any hint would be appreciated.
Actual error: error using odearguments @(T,Y)MBBRPDE(T,Y) returns a vector of length 0, but the length of initial conditions vector is 24. The vector returned by @(T,Y)MBBRPDE(T,Y) and the initial conditions vector must have the same number of elements.
%SCRIPT:
%initial conditions
y0=[1;1;1;1;1;1;1;1;1;1;1;1;1;1;1;1;0.4;0.2;0.1;0.3;0.1;1;1;1]';
tSpan=[0 180*24];
options = odeset('RelTol',1e-20,'AbsTol',1e-10,'InitialStep',1e-10, 'MStateDependence','strong');
[tSol, ySol]=ode15s(@(t,y) MBBRPDE(t,y), tSpan, y0, options);
Error using odearguments
@(T,Y)MBBRPDE(T,Y) returns a vector of length 0, but the length of initial conditions vector is 24. The vector returned by @(T,Y)MBBRPDE(T,Y) and the initial conditions vector must have the same number of elements.

Error in ode15s (line 153)
odearguments(odeIsFuncHandle, odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);
%Function:
function fval=MBBRPDE(t,y)
%Define the three variables
SseBL=y(1);
Sno3BL=y(2);
Sno2BL=y(3);
Snh4BL=y(4);
So2BL=y(5);
Xs=y(6);
muaob=y(7);
munb=y(8);
munsp=y(9);
muamx=y(10);
muhan=y(11);
SseBF=y(12);
Sno3BF=y(13);
Sno2BF=y(14);
Snh4BF=y(15);
So2BF=y(16);
faob=y(17);
fnb=y(18);
fnsp=y(19);
famx=y(20);
fhan=y(21);
muaver=y(22);
UL=y(23);
VL=y(24);
T=301;
%aob
O1=68/(0.00831*293*(273+T)); u1=1.443*exp(O1*8); Ko2aob=0.594; Knh4aob=1.5*2.4; Yaob=0.15; Kno3aob=0.5;
%nb
O2=44/(0.00831*293*(273+T)); u2=0.01*0.8*0.48*exp(O2*8); Ko2nb=0.98*0.435; Kno2nb=0.5*1.375; Ynb=1.0*0.13;Kno3nb=0.5;
%nsp
O3=44/(0.00831*293*(273+T)); u3=0.01*0.69*exp(O3*8); Ko2nsp=0.435; Kno2nsp=1.5*0.76; Ynsp=1.0*0.14;Kno3nsp=0.5;
%amx
O5=71/(0.00831*293*(273+T)); u5=5*0.0664*exp(O5*8); Kno2amx=0.2*0.175; Kno3amx=0.5; Knh4amx=0.2*0.185; Yamx=0.159; Ko2amx=0.1;
%han
O6=0.069;u6=7.2*exp(8*0.069); KSsehan=0.5*2; Kno2han=0*0.5;Kno3han=0.5*0.32; YNo2han=2*0.49; Ko2han=0.2;
YNo3han=0.5*0.79; Yhaer=0.37;
%other constants
baob=0.0054*exp(O6*(T-293)); bnb=3*0.0025*exp(O6*(T-293));bnsp=bnb; bamx=0.00013*exp(O6*(T-293)); bhan=0.008*exp(O6*(T-293));
INBM=0.07; fI=0.1; namx=0.5; nhan=0.6;naob=0.5; nnb=0.5; nnsp=0.5;
KX=1; KH=0.125*exp(0.069*(T-293)); INXI=0.02;
% input parameters
Sseo=300; Sno3o=0.2; Xso=355; Sno2o=0.2; So2o=0;
VR=350;
z1=10^-6;
N=length(UL);
hz=N/z1;
SseBL=zeros(N,1); Sno3BL=zeros(N,1); Sno2BL=zeros(N,1); Snh4BL=zeros(N,1); So2BL=zeros(N,1); Xs=zeros(N,1);
muaob=zeros(N,1); munb=zeros(N,1); munsp=zeros(N,1); muamx=zeros(N,1); muhan=zeros(N,1); SseBF=zeros(N,1);
Sno3BF=zeros(N,1); Sno2BF=zeros(N,1); Snh4BF=zeros(N,1); So2BF=zeros(N,1); faob=zeros(N,1); fnb=zeros(N,1);
fnsp=zeros(N,1); famx=zeros(N,1); fhan=zeros(N,1); muaver=zeros(N,1); UL=zeros(N,1); VL=zeros(N,1);
sigma=0; L=0.0001;
DLSse=4.8*10^-5; DLSno3=1.8*10^-4; DLSno2=1.8*10^-4; DLSnh4=1.8*10^-4; DLSo2=1.8*10^-4;
if t>=0 && t<5
Snh4o=916; Qo=0.0452*VL*10^3/(916*24);
elseif t>=5 && t<9
Snh4o=1078; Qo=0.073*VL*10^3/(916*24);
elseif t>=9 && t<12
Snh4o=909;Qo=0.132*VL*10^3/(909*24);
elseif t>=12 && t<=15
Snh4o=1078;Qo=0.191*VL*10^3/(1078*24);
elseif t>15 && t<=21
Snh4o=937;Qo=0.185*VL*10^3/(937*24);
elseif t>21 && t<=33
Snh4o=937;Qo=0.2286*VL*10^3/(937*24);
elseif t>33 && t<=35
Snh4o=937;Qo=0.3095*VL*10^3/(937*24);
elseif t>35 && t<=37
Snh4o=937;Qo=0.3158*VL*10^3/(937*24);
elseif t>37 && t<=49
Snh4o=937;Qo=0.348*VL*10^3/(937*24);
elseif t>49 && t<=77
Snh4o=937;Qo=0.448*VL*10^3/(937*24);
elseif t>77 && t<91
Snh4o=1078; Qo=0.3928*VL*10^3/(1078*24);
elseif t>=91 && t<=119
Snh4o=959;Qo=0.4582*VL*10^3/(959*24);
elseif t>119 && t<=127
Snh4o=959;Qo=0.1738*VL*10^3/(959*24);
elseif t>127 && t<=134
Snh4o=959;Qo=0.6246*VL*10^3/(959*24);
elseif t>134 && t<=140
Snh4o=1190;Qo=0.625*VL*10^3/(1190*24);
elseif t>140 && t<=148
Snh4o=1293;Qo=0.3929*VL*10^3/(1293*24);
elseif t>148 && t<=155
Snh4o=1087;Qo=0.8252*VL*10^3/(1087*24);
elseif t>155 && t<=161
Snh4o=1003;Qo=0.47*VL*10^3/(1003*24);
elseif t>161 && t<=169
Snh4o=1078;Qo=0.564*VL*10^3/(1078*24);
else
Snh4o=928; Qo=0.47*VL*10^3/(928*24);
end
KaL=80;
AL=120400; LL=1;
DSse=5.8*10^-5; DSno3=1.4*10^-5; DSno2=1.4*10^-4; DSnh4=1.5*10^-4; DSo2=2.2*10^-4;
%% PDEs
for i=2:hz:N-1
%BC-1
SseBF(1)=0; Sno3BF(1)=0; Sno2BF(1)=0; Snh4BF(1)=0; So2BF(1)=0; faob(1)=0.4;fnb(1)=0.1;fnsp(1)=0.1;famx(1)=0.3;fhan(1)=0.1;
%BC-2
SseBF(end)=L*DLSse.*(SseBL(end)-SseBF(end))./LL; Sno3BF(end)=L*DLSno3.*(Sno3BL(end)-Sno3BF(end))/LL; Sno2BF(end)=L*DLSno2.*(Sno2BL(end)-Sno2BF(end))/LL; Snh4BF(end)=L*DLSnh4.*(Snh4BL(end)-Snh4BF(end))/LL;...
So2BF(end)=L*DLSo2.*(So2BL(end)-So2BF(end))./LL; SseBF(end)=SseBL(end); Sno3BF(end)=Sno3BL(end); Sno2BF(end)=Sno2BL(end); Snh4BF(end)=Snh4BL(end); So2BF(end)=So2BL(end);
SseBL(i)=Qo.*(Sseo-SseBL(i))./VL(i)-AL.*((DSse/LL)-UL(i)).*(SseBL(i)-SseBF(i));
Sno3BL(i)=Qo.*(Sno3o-Sno3BL(i))/VL(i)-AL.*((DSno3/LL)-UL(i)).*(Sno3BL(i)-Sno3BF(i));
Sno2BL(i)=Qo.*(Sno2o-Sno2BL(i))/VL(i)-AL.*((DSno2/LL)-UL(i)).*(Sno2BL(i)-Sno2BF(i));
Snh4BL(i)=Qo.*(Snh4o-Snh4BL(i))/VL(i)-AL.*((DSnh4/LL)-UL(i)).*(Snh4BL(i)-Snh4BF(i));
So2BL(i)=Qo.*(So2o-So2BL(i))/VL(i) + KaL.*(8-So2BL(i)) - AL.*((DSo2)-UL(i)).*(So2BL(i)-So2BF(i));
Xs(i)=Qo.*(Xso-Xs(i))/VL(i) - KH.*((Xs(i)./fhan(i))./(KX+(Xs(i)/fhan(i)))).*fhan(i);
muaob(i)=u1.*(So2BF(i)./(Ko2aob+So2BF(i))).*(Snh4BF(i)./(Knh4aob+Snh4BF(i))).*faob(i) - baob.*(So2BF(i)./(Ko2aob+So2BF(i))).*faob(i) - baob.*naob.*((Ko2aob./(Ko2aob+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3aob+Sno2BF(i)+Sno3BF(i))).*faob(i);
munb(i)=u2.*(Sno2BF(i)./(Kno2nb+Sno2BF(i))).*(So2BF(i)./(Ko2nb+So2BF(i))).*fnb(i) - bnb.*(So2BF(i)./(Ko2nb+So2BF(i))).*fnb(i) - bnb*nnb.*((Ko2nb./(Ko2nb+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3nb+Sno2BF(i)+Sno3BF(i))).*fnb(i);
munsp(i)=u3.*(Sno2BF(i)./(Kno2nsp+Sno2BF(i))).*(So2BF(i)./(Ko2nsp+So2BF(i))).*fnsp(i) - bnsp.*(So2BF(i)./(Ko2nsp+So2BF(i))).*fnsp(i) - bnsp*nnsp.*((Ko2nsp./(Ko2nsp+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3nsp+Sno2BF(i)+Sno3BF(i))).*fnsp(i);
muamx(i)=u5.*((Snh4BF(i)./(Knh4amx+Snh4BF(i))).*(Sno2BF(i)./(Kno2amx+Sno2BF(i))).*(Ko2amx/(Ko2amx+So2BF(i)))).*famx(i) - bamx.*(So2BF(i)./(Ko2amx+So2BF(i))).*famx(i) - bamx.*namx.*((Ko2amx./(Ko2amx+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3amx+Sno2BF(i)+Sno3BF(i))).*famx(i);
muhan(i)=(u6.*(SseBF(i)./(KSsehan+SseBF(i))).*(Sno2BF(i)./(Kno2han+Sno2BF(i))).*(Ko2han/(Ko2han+So2BF(i))) + u6.*(SseBF(i)./(KSsehan+SseBF(i))).*(Sno3BF(i)./(Kno3han+Sno3BF(i))).*(Ko2han./(Ko2han+So2BF(i)))+u6.*((SseBF(i)./(KSsehan+SseBF(i))).*(So2BF(i)./(Ko2han+So2BF(i)))).*fhan(i))- bhan.*(So2BF(i)./(Ko2han+So2BF(i))).*fhan(i) - bhan.*nhan.*((Ko2han/(Ko2han+So2BF(i))).*(Sno2BF(i+1)+Sno3BF(i))./(Kno3han+Sno2BF(i)+Sno3BF(i))).*fhan(i);
SseBF(i)=DSse.*(SseBF(i+1)-2.*SseBF(i)+(SseBF(i-1)))./(hz.*i)^2 +i.*UL(i).*(SseBF(i+1)-SseBF(i))./L - (u6*nhan.*(1/YNo2han).*(SseBF(i)./(KSsehan+SseBF(i))).*(Sno2BF(i)./(Kno2han+Sno2BF(i))).*(Ko2han./(Ko2han+So2BF(i)))).*fhan(i) -(u6*nhan.*(1/YNo3han).*(SseBF(i)./(KSsehan+SseBF(i))).*(Sno3BF(i)./(Kno3han+Sno3BF(i))).*(Ko2han./(Ko2han+So2BF(i)))).*fhan(i) - (u6/Yhaer).*(SseBF(i)./(KSsehan+SseBF(i))).*(So2BF(i)./(Ko2han+So2BF(i))).*fhan(i);
Sno3BF(i)=DSno3.*(Sno3BF(i+1)-2.*Sno3BF(i)+(Sno3BF(i-1)))./(hz.*i)^2 +i.*UL(i).*(Sno3BF(i+1)-Sno3BF(i))./L - (u6*nhan*((1-YNo3han)/(2.86*YNo3han)).*(SseBF(i)./(KSsehan+SseBF(i))).*(Sno3BF(i)./(Kno3han+Sno3BF(i))).*(Ko2han./(Ko2han+So2BF(i)))).*fhan(i) + (u5/1.14).*(Snh4BF(i)./(Knh4amx+Snh4BF(i))).*(Sno2BF(i)./(Kno2amx+Sno2BF(i))).*(Ko2amx./(Ko2amx+So2BF(i))).*famx(i) + (u2/Ynb*(Sno2BF(i)./(Kno2nb+Sno2BF(i))).*(So2BF(i)./(Ko2nb+So2BF(i))).*fnb(i)) + ((u3/Ynsp).*(Sno2BF(i)./(Kno2nsp+Sno2BF(i))).*(So2BF(i)./(Ko2nsp+So2BF(i))).*fnsp(i)) - ((1-fI)/2.86)*baob*naob.*(Ko2aob./(Ko2aob+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3aob+Sno2BF(i)+Sno3BF(i)).*faob(i) - ((1-fI)/2.86)*bnb*nnb.*(Ko2nb/(Ko2nb+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3nb+Sno2BF(i)+Sno3BF(i)).*fnb(i) - ((1-fI)/2.86)*bnsp*nnsp.*(Ko2nsp/(Ko2nsp+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3nsp+Sno2BF(i)+Sno3BF(i)).*fnsp(i) - ((1-fI)/2.86)*bamx*namx.*(Ko2amx/(Ko2amx+So2BF(i))).*(Sno2BF(i)+Sno3(i))./(Kno3amx+Sno2BF(i)+Sno3BF(i)).*famx(i)-((1-fI)/2.86)*bhan*nhan.*(Ko2han/(Ko2han+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3han+Sno2BF(i)+Sno3BF(i)).*fhan(i);
Sno2BF(i)=DSno2.*(Sno2BF(i+1)-2.*Sno2BF(i)+(Sno2BF(i-1)))./(hz.*i)^2 +i.*UL(i).*(Sno2BF(i+1)-Sno2BF(i))./L - (((1-YNo2han)/(1.71*YNo2han))*u6*nhan.*(SseBF(i)./(KSsehan+SseBF(i))).*(Sno2BF(i)./(Kno2han+Sno2BF(i))).*(Ko2han/(Ko2han+So2BF(i))).*fhan(i)) - (((u5/Yamx)+(u5/1.14)).*(Snh4BF(i)./(Knh4amx+Snh4BF(i))).*(Sno2BF(i)/(Kno2amx+Sno2BF(i))).*(Ko2amx/(Ko2amx+So2BF(i)))).*famx(i) + ((u1/Yaob).*(So2BF(i)/(Ko2aob+So2(i))).*(Snh4BF(i)./(Knh4aob+Snh4BF(i)))).*faob(i) - ((u2/Ynb).*(Sno2BF(i)./(Kno2nb+Sno2BF(i))).*(So2BF(i)./(Ko2nb+So2BF(i))).*fnb(i)) - ((u3/Ynsp).*(Sno2BF(i)./(Kno2nsp+Sno2BF(i))).*(So2BF(i)./(Ko2nsp+So2BF(i)))).*fnsp(i);
Snh4BF(i)=DSnh4.*(Snh4BF(i+1)-2.*Snh4BF(i)+(Snh4BF(i-1)))./(hz.*i)^2 +i.*UL(i).*(Snh4BF(i+1)-Snh4BF(i))./L - (u1*(INBM+1/Yaob).*(So2BF(i)./(Ko2aob+So2BF(i))).*(Snh4BF(i)./(Knh4aob+Snh4BF(i)))).*faob(i) -(u5*(INBM+1/Yamx).*(Snh4BF(i)./(Knh4amx+Snh4BF(i))).*(Sno2BF(i)./(Kno2amx+Sno2BF(i))).*(Ko2amx./(Ko2amx+So2BF(i)))).*famx(i) - (u2*INBM.*(Sno2BF(i)./(Kno2nb+Sno2BF(i))).*(So2BF(i)./(Ko2nb+So2BF(i)))).*(i) -(u3*INBM.*(Sno2BF(i)./(Kno2nsp+Sno2BF(i))).*(So2BF(i)./(Ko2nsp+So2BF(i)))).*fnsp(i) - INBM*u6*nhan.*((SseBF(i)./(KSsehan+SseBF(i))).*(Sno2BF(i)./(Kno2han+Sno2BF(i))).*(Ko2han/(Ko2han+So2BF(i)))).*fhan(i) - u6*INBM*nhan.*((SseBF(i)./(KSsehan+SseBF(i))).*(Sno3BF(i)./(Kno3han+Sno3BF(i))).*(Ko2han./(Ko2han+So2BF(i)))).*fhan(i) - u6*INBM.*((SseBF(i)./(KSsehan+SseBF(i))).*(So2BF(i)./(Ko2han+So2BF(i)))).*fhan(i) + (INBM-(fI*INXI))*baob*naob.*((Ko2aob/(Ko2aob+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3aob+Sno2BF(i)+Sno3BF(i))).*faob(i) + (INBM-(fI*INXI))*bnb*nnb.*((Ko2nb/(Ko2nb+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3nsp+Sno2BF(i)+Sno3BF(i))).*fnb(i) +(INBM-(fI*INXI))*bnsp*nnsp.*((Ko2nsp./(Ko2nsp+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3nsp+Sno2BF(i)+Sno3BF(i))).*fnsp(i) + (INBM-(fI*INXI))*bamx*namx.*((Ko2amx./(Ko2amx+So2BF(i))).*(Sno2BF(i+1)+Sno3BF(i))./(Kno3amx+Sno2BF(i)+Sno3BF(i))).*famx(i) + (INBM-(fI*INXI))*bhan*nhan.*((Ko2han/(Ko2han+So2BF(i))).*(Sno2BF(i)+Sno3BF(i))./(Kno3han+Sno2BF(i)+Sno3BF(i))).*fhan(i) +(INBM-(fI*INXI))*baob.*(So2BF(i)./(Ko2aob+So2BF(i))).*faob(i) + (INBM-(fI*INXI))*bnb.*(So2BF(i)./(Ko2nb+So2BF(i))).*fnb(i) +(INBM-(fI*INXI))*bnsp.*(So2BF(i)./(Ko2nsp+So2BF(i))).*fnsp(i) + (INBM-(fI*INXI))*bamx.*(So2BF(i)./(Ko2amx+So2BF(i))).*famx(i) + (INBM-(fI*INXI))*bhan.*(So2BF(i)./(Ko2han+So2BF(i))).*fhan(i);
So2BF(i)=DSo2.*(So2BF(i+1)-2.*So2BF(i)+(So2BF(i-1)))./(hz.*i)^2 +i.*UL(i).*(So2BF(i+1)-So2BF(i))./L - -(((1-Yhaer)/Yhaer)*u6.*(SseBF(i)./(KSsehan+SseBF(i))).*(So2BF(i)./(Ko2han+So2BF(i)))).*fhan(i) -(((3.43-Yaob)/Yaob)*u1.*(So2BF(i)./(Ko2aob+So2BF(i))).*(Snh4BF(i)./(Knh4aob+Snh4BF(i)))).*faob(i) - (((1.14-Ynb)/Ynb)*u2.*(Sno2BF(i)/(Kno2nb+Sno2BF(i))).*(So2BF(i)./(Ko2nb+So2BF(i)))).*fnb(i) - (((1.14-Ynsp)/Ynsp)*u3.*(Sno2BF(i)/(Kno2nsp+Sno2BF(i))).*(So2BF(i)./(Ko2nsp+So2BF(i)))).*fnsp(i) - (1-fI)*baob.*(So2BF(i)./(Ko2aob+So2BF(i))).*faob(i) - (1-fI)*bnb.*(So2BF(i)./(Ko2nb+So2BF(i))).*fnb -(1-fI)*bnsp.*(So2BF(i)./(Ko2nsp+So2BF(i))).*fnsp(i) - (1-fI)*bamx.*(So2BF(i)./(Ko2amx+So2BF(i))).*famx(i) - (1-fI)*bhan.*(So2BF(i)./(Ko2han+So2BF(i))).*fhan(i);
faob(i)=(muaob(i)-muaver(i)).*faob(i) - (U-UL(i).*UL(i)).*(faob(i+1)-faob(i))./L;
fnb(i)=(munb(i)-muaver(i)).*fnb(i) - (U-UL(i).*UL(i)).*(fnb(i+1)-fnb(i))./L;
fnsp(i)=(munsp(i)-muaver(i)).*fnsp(i) - (U-UL(i).*UL(i)).*(fnsp(i+1)-fnsp(i))./L;
famx(i)=(muamx(i)-muaver(i)).*famx(i) - (U-UL(i).*UL(i)).*(famx(i+1)-famx(i))./L;
fhan(i)=(muhan(i)-muaver(i)).*fhan(i) - (U-UL(i).*UL(i)).*(fhan(i+1)-fahan(i))./L;
muaver(i)=faobi(i).*muaob(i) + fnb(i).*munb(i) + fnsp(i).*munsp(i) + famx(i).*muamx(i) + fhan(i).*muhan(i);
UL(i)=L.*muaver(i).*i+sigma;
VL(i)=VR-AL.*(i + LL);
end
fval=[SseBL(i);Sno3BL(i);Sno2BL(i);Snh4BL(i);So2BL(i);Xs(i);muaob(i);munb(i);munsp(i);muamx(i);muhan(i);SseBF(i);Sno3BF(i);Sno2BF(i);Snh4BF(i);So2BF(i);faob(i);fnb(i);fnsp(i);famx(i);fhan(i);muaver(i);UL(i);VL(i)];
end

请先登录，再回答此问题。

采纳的回答

Cris LaPierre
Cris LaPierre 2023-7-23
Ran in:
Your for loop is not executing because, after plugging in variables, your loop counter is
for i=2:1000000:0
end
i
i = []
so the value of I is [], resulting in fval being empty. Check the values you are using for your loop counter and correct as necessary.
  19 个评论
显示 17更早的评论
Torsten
Torsten 2023-7-27
Ran in:
tstart = 0;
tend = 3;
nt = 10;
xstart = 0;
xend = 1;
nx = 50;
t = linspace(tstart,tend,nt);
x = linspace(xstart,xend,nx);
% Method-of-lines solution
x = x.';
xhalf = (x(1:nx-1)+x(2:nx))/2;
% Enlarge x and xhalf by a ghostcell
x = [x;x(end)+(x(end)-x(end-1))];
xhalf = [xhalf;(x(end)+x(end-1))/2];
a=arrayfun(@(i)oscic(x(i)),1:nx,'UniformOutput',false);
b=cell2mat(a);
y01 = b(1,:).';
y02 = b(2,:).';
y01(1) = 0;
y02(1) = 0;
y1=[y01; y02];
[tSol,ySol]=ode15s(@(t,y) fun(t,y,x,xhalf,nx),t,y1);
figure(1)
plot(x(1:end-1),ySol(:,1:nx))
figure(2)
plot(x(1:end-1),ySol(:,nx+1:end))
function dydt=fun(t,y,x,xhalf,nx)
dydt = zeros(length(y),1);
u1 = y(1:nx,:);
u2 = y(nx+1:2*nx,:);
u1nxplus1 = u1(nx-1) + (u1(nx)-300)*(x(nx+1)-x(nx-1));
u2nxplus1 = u2(nx-1);
u1 = [u1;u1nxplus1];
u2 = [u2;u2nxplus1];
du1dt = zeros(nx,1);
du2dt = zeros(nx,1);
for ix = 2:nx
du1dt(ix) = (u1(ix+1)-2*u1(ix)+u1(ix-1))/(x(ix+1)-x(ix))^2;
du2dt(ix) = (u2(ix+1)-2*u2(ix)+u2(ix-1))/(x(ix+1)-x(ix))^2 + u2(ix)*(1-u1(ix)-u2(ix));
end
dydt = [du1dt;du2dt];
end
function u0 = oscic(x)
u0 = [600;x*(x-2)];
end
KIPROTICH KOSGEY
KIPROTICH KOSGEY 2023-8-10
移动：Torsten 2023-8-10
Hi Torsten
Thank you so much for the assistance you have accorded me all along. I added an ode in the above simple system of equations and everything worked fine. However, when I applied in my system of equations that I captured in the attachment previously given, it is giving an error relating to the boundary conditions: Index exceeds the number of array elements. Index must not exceed 1.
Error in MBBR>fun (line 90)
SseBFmax=SseBF(nx-1)+L.*DLSse.*(SseBL(nx)-SseBF(nx))./LL; Sno3BFmax=Sno3BF(nx-1)+L.*DLSno3.*(Sno3BL(nx)-Sno3BF(nx))./LL; Sno2BFmax=Sno2BF(nx-1)+L.*DLSno2.*(Sno2BL(nx)-Sno2BF(nx))/LL; Snh4BFmax=Snh4BF(nx-1)+L.*DLSnh4.*(Snh4BL(nx)-Snh4BF(nx))/LL;
the left boundary conditions are ode (dfi/dt=(mui-muaver)*fi) i=aob,nb,nsp,am,han
and also pde(dsi/dx=0) i=SseBF,Sno3BF,Sno2BF,So2BF,Snh4BF
right boundary conditions are pde(L*DLi*(SLi-Si))/LL*Di, Li=SseBL,Sno3BL,Sno2BL,So2BL,Snh4BL and i=SseBF,Sno3BF,Sno2BF,So2BF,Snh4BF
I have attached the code .
%BC-1
SseBFmin=SseBF(2); Sno3BFmin=Sno3BF(2); Sno2BFmin=Sno2BF(2); Snh4BFmin=Snh4BF(2); So2BFmin=So2BF(2);
faobmin=exp((muaob(1)-muaver(1))*0.0003);fnbmin=exp((munb(1)-muaver(1))*0.0003);fnspmin=exp((munsp(1)-muaver(1))*0.0003);famxmin=exp((muamx(1)-muaver(1))*0.0003);fhanmin=exp((muhan(1)-muaver(1))*0.0003);
%BC-2
SseBFmax=SseBF(nx-1)+L.*DLSse.*(SseBL(nx)-SseBF(nx))./LL; Sno3BFmax=Sno3BF(nx-1)+L.*DLSno3.*(Sno3BL(nx)-Sno3BF(nx))./LL; Sno2BFmax=Sno2BF(nx-1)+L.*DLSno2.*(Sno2BL(nx)-Sno2BF(nx))/LL; Snh4BFmax=Snh4BF(nx-1)+L.*DLSnh4.*(Snh4BL(nx)-Snh4BF(nx))/LL;
So2BFmax=So2BF(nx-1)+L.*DLSo2.*(So2BL(nx)-So2BF(nx))./LL;
SseBF=[SseBFmin;SseBFmax];Sno3BF=[Sno3BFmin;Sno3BFmax];Sno2BF=[Sno2BFmin;Sno2BFmax];Snh4BF=[Snh4BFmin;Snh4BFmax];So2BF=[So2BFmin;So2BFmax];
faob=[faobmin;faob];fnb=[fnbmin;fnb];fnsp=[fnspmin;fnsp];famx=[famxmin;famx];fhan=[fhanmin;fhan];

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 1-D Partial Differential Equations 的更多信息

标签

产品


版本

R2023a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by