I am struggling to set a maximum value for a output for a function. I have nested my function in another but it is giving an error.

2 次查看（过去 30 天）
KIPROTICH KOSGEY
KIPROTICH KOSGEY 2023-7-13
编辑： Torsten 2023-7-13
I have the below script and main function. After nesting so as to introduce a limit for maxSnh4, it is giving an error. Any hint on how to restructure it would be highly appreciated.
%Script:
tSpan=[0 300];
y0=[1;1;1;1;1;1;1;1;1;1;1];
M1=eye(55);
options = odeset('RelTol',1e-20,'AbsTol',1e-10,'Mass',M1,'MassSingular', 'maybe','MStateDependence','strong','InitialStep',1e-10, 'MStateDependence','strong');
[tSol, ySol]=ode15s(@(t,y) SBR(t,y), tSpan, y0);
%plotting
figure (1)
plot(tSol,ySol(:,7),'-r'); hold on %Sse
plot(tSol,ySol(:,8),'-g'); hold on %Sno3
plot(tSol,ySol(:,9),'--r'); hold on %Sno2
plot(tSol,ySol(:,10),'-k'); hold on %nh4
plot(tSol,ySol(:,11),'-b'); hold on %So2
figure (2)
plot(tSol,ySol(:,1),'-r'); hold on %Xaob
plot(tSol,ySol(:,2),'-g'); hold on %Xnb
plot(tSol,ySol(:,3),'--r'); hold on %Xnsp
plot(tSol,ySol(:,4),'-k'); hold on %amx
plot(tSol,ySol(:,5),'-b'); hold on %Xhan
plot(tSol,ySol(:,6),'--g'); hold off %Xs
xlabel('Time (days)'); ylabel('Concentration (g/m^3)');
%Main function
% SBR simulation in MATLAB
function dy=SBR(t,y)
maxSnh4=200;
function fval=SBR1(t,y)
% Initialize arrays
Xaob=y(1);
Xnb = y(2);
Xnsp= y(3);
Xamx= y(4);
Xhan= y(5);
Xs= y(6);
Sse= y(7);
Sno3= y(8);
Sno2= y(9);
Snh4= y(10);
So2= y(11);
T=303;
%aob
O1=68/(0.00831*293*(273+T)); u1=0.054*exp(O1*(T-293)); Ko2aob=0.6; Knh4aob=2.4; Yaob=0.15; Kno3aob=0.5;
%nb
O2=44/(0.00831*293*(273+T)); u2=0.047*exp(O2*(T-293)); Ko2nb=2.2; Kno2nb=0.5*1.375; Ynb=0.041;Kno3nb=0.5;
%nsp
O3=44/(0.00831*293*(273+T)); u3=0.69*exp(O3*(T-293)); Ko2nsp=0.33; Kno2nsp=0.52; Ynsp=0.14;Kno3nsp=0.5;
%amx
O5=71/(0.00831*293*(273+T)); u5=0.0022*exp(O5*(T-293)); Kno2amx=0.05; Knh4amx=0.07; Yamx=0.159; Ko2amx=0.01; Kno3amx=0.5;
%han
O6=0.069;u6=0.5*7.2*exp((T-293)*0.069); KSsehan=2; Kno2han=0.5;Kno3han=0.5; YNo2han=0.3; Ko2han=0.2;
YNo3han=0.43; Yhaer=0.54;
%other constants
baob=0.0054*exp(O6*(T-293)); bnb=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; Qo=116.67; Xso=250;
Snh4o=650; Sno2o=0.2; So2o=0;
%oxygen transfer parameters
KaL=1920;
V=1400; SRT=100;
Xaobo=2; Xnbo=2; Xnspo=2; Xhano=2; Xamxo=2;
Vo=V-0.3*V; % volume in the reactor at the start of refilling
% SBR operation cycle
refilling_duration = 1.59; % Duration of refilling phase (min)
aeration_duration = 15; % Duration of aeration phase (min)
mixing_duration = 45; % Duration of mixing phase (min)
settling_duration = 30; % Duration of settling phase (min)
withdrawal_duration = 0.58; % Duration of withdrawal phase (min)
cycle_duration=refilling_duration+aeration_duration+mixing_duration+settling_duration+withdrawal_duration;
t0=0+cycle_duration;
t1=t0+15;
t2=t1+15;
t3=t2+45;
t4=t3+60;
t5=t4+15;
%refilling
while t>t0 && t<=t1
% input parameters
Sseo=300; Sno3o=0.2; Qo=220; Xso=250;
Snh4o=650; Sno2o=0.2; So2o=0; KaL=0;
Xaobo=2; Xnbo=2; Xnspo=2; Xhano=2; Xamxo=2;
break
end
%aeration
while t>t1 && t<=t2
% input parameters
Sseo=0; Sno3o=0; Qo=0; Xso=0; Snh4o=0; Sno2o=0; So2o=0; Xaobo=0; Xnbo=0; Xnspo=0; Xhano=0; Xamxo=0;
if Snh4<10 || Sno2>10 || So2>0.8
KaL=0;
elseif Snh4>10 && So2<=0.5
KaL=1920;
elseif Sno2<10 && So2<0.5
KaL=1920;
end
break
end
%mixing & settling
while t>t2 && t<=t4
% input parameters
Sseo=0; Sno3o=0; Qo=0; Xso=0; Snh4o=0; Sno2o=0.2; So2o=0; KaL=0; Xaobo=0; Xnbo=0; Xnspo=0; Xhano=0; Xamxo=0;
break
end
%withdrawal
while t>t4 && t<=t5
% input parameters
Sseo=0; Sno3o=0; Qo=600; Xso=0; Snh4o=0; Sno2o=0.2; So2o=0; KaL=0; Xaobo=0; Xnbo=0; Xnspo=0; Xhano=0; Xamxo=0;
break
end
fval(1)=((Vo*Xaobo-Vo*Xaob/SRT)/V + u1*(So2/(Ko2aob+So2))*(Snh4/(Knh4aob+Snh4))*Xaob - baob*(So2/(Ko2aob+So2))*Xaob - baob*naob*((Ko2aob/(Ko2aob+So2))*(Sno2+Sno3)/(Kno3aob+Sno2+Sno3))*Xaob);
fval(2)=((Vo*Xnbo-Vo*Xnb/SRT)/V + u2*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2))*Xnb - bnb*(So2/(Ko2nb+So2))*Xnb - bnb*nnb*((Ko2nb/(Ko2nb+So2))*(Sno2+Sno3)/(Kno3nb+Sno2+Sno3))*Xnb);
fval(3)=((Vo*Xnspo-Vo*Xnsp/SRT)/V + u3*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2))*Xnsp - bnsp*(So2/(Ko2nsp+So2))*Xnsp - bnsp*nnsp*((Ko2nsp/(Ko2nsp+So2))*(Sno2+Sno3)/(Kno3nsp+Sno2+Sno3))*Xnsp);
fval(4)=((Vo*Xamxo-Vo*Xamx/SRT)/V + u5*((Snh4/(Knh4amx+Snh4))*(Sno2/(Kno2amx+Sno2))*(Ko2amx/(Ko2amx+So2)))*Xamx-bamx*(So2/(Ko2amx+So2))*Xamx - bamx*namx*((Ko2amx/(Ko2amx+So2))*(Sno2+Sno3)/(Kno3amx+Sno2+Sno3))*Xamx);
fval(5)=((Vo*Xhano-Vo*Xhan/SRT)/V + (u6*(Sse/(KSsehan+Sse))*(Sno2/(Kno2han+Sno2))*(Ko2han/(Ko2han+So2)) + u6*(Sse/(KSsehan+Sse))*(Sno3/(Kno3han+Sno3))*(Ko2han/(Ko2han+So2))+u6*((Sse/(KSsehan+Sse))*(So2/(Ko2han+So2)))*Xhan)- bhan*(So2/(Ko2han+So2))*Xhan - bhan*nhan*((Ko2han/(Ko2han+So2))*(Sno2+Sno3)/(Kno3han+Sno2+Sno3))*Xhan);
fval(6)=((Qo*Xso-Qo*Xs)/V - KH*((Xs/Xhan)/(KX+(Xs/Xhan)))*Xhan);
fval(7)=((Qo*Sseo-Qo*Sse)/V - (u6*nhan*(1/YNo2han)*(Sse/(KSsehan+Sse))*(Sno2/(Kno2han+Sno2))*(Ko2han/(Ko2han+So2)))*Xhan -(u6*nhan*(1/YNo3han)*(Sse/(KSsehan+Sse))*(Sno3/(Kno3han+Sno3))*(Ko2han/(Ko2han+So2)))*Xhan - ((u6/Yhaer)*(Sse/(KSsehan+Sse))*(So2/(Ko2han+So2))*Xhan) + KH*((Xs/Xhan)/(KX+(Xs/Xhan)))*Xhan);
fval(8)=((Qo*Sno3o-Qo*Sno3)/V - (u6*nhan*((1-YNo3han)/(2.86*YNo3han))*(Sse/(KSsehan+Sse))*(Sno3/(Kno3han+Sno3))*(Ko2han/(Ko2han+So2)))*Xhan + (u5/1.14)*(Snh4/(Knh4amx+Snh4))*(Sno2/(Kno2amx+Sno2))*(Ko2amx/(Ko2amx+So2))*Xamx + (u2/Ynb*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2))*Xnb) + ((u3/Ynsp)*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2))*Xnsp) - ((1-fI)/2.86)*baob*naob*(Ko2aob/(Ko2aob+So2))*(Sno2+Sno3)/(Kno3aob+Sno2+Sno3)*Xaob - ((1-fI)/2.86)*bnb*nnb*(Ko2nb/(Ko2nb+So2))*(Sno2+Sno3)/(Kno3nb+Sno2+Sno3)*Xnb - ((1-fI)/2.86)*bnsp*nnsp*(Ko2nsp/(Ko2nsp+So2))*(Sno2+Sno3)/(Kno3nsp+Sno2+Sno3)*Xnsp - ((1-fI)/2.86)*bamx*namx*(Ko2amx/(Ko2amx+So2))*(Sno2+Sno3)/(Kno3amx+Sno2+Sno3)*Xamx-((1-fI)/2.86)*bhan*nhan*(Ko2han/(Ko2han+So2))*(Sno2+Sno3)/(Kno3han+Sno2+Sno3)*Xhan);
fval(9)=((Qo*Sno2o-Qo*Sno2)/V - (((1-YNo2han)/(1.71*YNo2han))*u6*nhan*(Sse/(KSsehan+Sse))*(Sno2/(Kno2han+Sno2))*(Ko2han/(Ko2han+So2))*Xhan) - (((u5/Yamx)+(u5/1.14))*(Snh4/(Knh4amx+Snh4))*(Sno2/(Kno2amx+Sno2))*(Ko2amx/(Ko2amx+So2)))*Xamx + ((u1/Yaob)*(So2/(Ko2aob+So2))*(Snh4/(Knh4aob+Snh4)))*Xaob - ((u2/Ynb)*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2))*Xnb) - ((u3/Ynsp)*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2)))*Xnsp);
fval(10)=((Qo*Snh4o-Qo*Snh4)/V - (u1*(INBM+1/Yaob)*(So2/(Ko2aob+So2))*(Snh4/(Knh4aob+Snh4)))*Xaob -(u5*(INBM+1/Yamx)*(Snh4/(Knh4amx+Snh4))*(Sno2/(Kno2amx+Sno2))*(Ko2amx/(Ko2amx+So2)))*Xamx - (u2*INBM*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2)))*Xnb -(u3*INBM*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2)))*Xnsp - INBM*u6*nhan*((Sse/(KSsehan+Sse))*(Sno2/(Kno2han+Sno2))*(Ko2han/(Ko2han+So2)))*Xhan - u6*INBM*nhan*((Sse/(KSsehan+Sse))*(Sno3/(Kno3han+Sno3))*(Ko2han/(Ko2han+So2)))*Xhan - u6*INBM*((Sse/(KSsehan+Sse))*(So2/(Ko2han+So2)))*Xhan + (INBM-(fI*INXI))*baob*naob*((Ko2aob/(Ko2aob+So2))*(Sno2+Sno3)/(Kno3aob+Sno2+Sno3))*Xaob + (INBM-(fI*INXI))*bnb*nnb*((Ko2nb/(Ko2nb+So2))*(Sno2+Sno3)/(Kno3nsp+Sno2+Sno3))*Xnb +(INBM-(fI*INXI))*bnsp*nnsp*((Ko2nsp/(Ko2nsp+So2))*(Sno2+Sno3)/(Kno3nsp+Sno2+Sno3))*Xnsp + (INBM-(fI*INXI))*bamx*namx*((Ko2amx/(Ko2amx+So2))*(Sno2+Sno3)/(Kno3amx+Sno2+Sno3))*Xamx + (INBM-(fI*INXI))*bhan*nhan*((Ko2han/(Ko2han+So2))*(Sno2+Sno3)/(Kno3han+Sno2+Sno3))*Xhan +(INBM-(fI*INXI))*baob*(So2/(Ko2aob+So2))*Xaob + (INBM-(fI*INXI))*bnb*(So2/(Ko2nb+So2))*Xnb +(INBM-(fI*INXI))*bnsp*(So2/(Ko2nsp+So2))*Xnsp + (INBM-(fI*INXI))*bamx*(So2/(Ko2amx+So2))*Xamx + (INBM-(fI*INXI))*bhan*(So2/(Ko2han+So2))*Xhan);
fval(11)=((Qo*So2o-Qo*So2)/V + KaL*(8-So2) -(((1-Yhaer)/Yhaer)*u6*(Sse/(KSsehan+Sse))*(So2/(Ko2han+So2)))*Xhan -(((3.43-Yaob)/Yaob)*u1*(So2/(Ko2aob+So2))*(Snh4/(Knh4aob+Snh4)))*Xaob - (((1.14-Ynb)/Ynb)*u2*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2)))*Xnb - (((1.14-Ynsp)/Ynsp)*u3*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2)))*Xnsp - (1-fI)*baob*(So2/(Ko2aob+So2))*Xaob - (1-fI)*bnb*(So2/(Ko2nb+So2))*Xnb -(1-fI)*bnsp*(So2/(Ko2nsp+So2))*Xnsp - (1-fI)*bamx*(So2/(Ko2amx+So2))*Xamx - (1-fI)*bhan*(So2/(Ko2han+So2))*Xhan);
if fval(10)>maxSnh4
dy=maxSnh4;
else
dy=fval(10);
end
end
end
  3 个评论
显示 1更早的评论
KIPROTICH KOSGEY
KIPROTICH KOSGEY 2023-7-13
maxSnh4 is defined in row 3.
It is my first time coding nested functions and I am struggling to structure it properly hence the errors.
the error is: Output argument "dy" (and possibly others) not assigned a value in the execution with "SBR" function.
I just need to set fval(10) to <=200. Any hint would be greatly appreciated.
Stephen23
Stephen23 2023-7-13
"Output argument "dy" (and possibly others) not assigned a value in the execution with "SBR" function."
Because you never define the variable DY inside of the function SBR.
You also never call the nested function SBR1, so you might as well just delete it. Your code is exactly equivalent to this:
function dy=SBR(t,y)
maxSnh4=200;
end
which clearly does not define DY anywhere.

请先登录，再进行评论。

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

采纳的回答

Torsten
Torsten 2023-7-13
移动：Torsten 2023-7-13
Ran in:
Instead of programming if-statements within the function routine and setting parameters depending on the actual time t, you should integrate up to t0, return control to the calling program and save results, reset parameters, integrate from t0 to t1, return control to the calling program and save results, reset parameters, integrate ... This is the usual procedure to avoid discontinuities in the ODE equations which most probably makes ode15s stop with an error message at time t ~ 107.
%Script:
tSpan=[0 300];
y0=[1;1;1;1;1;1;1;1;1;1;1];
M1=eye(55);
options = odeset('RelTol',1e-20,'AbsTol',1e-10,'Mass',M1,'MassSingular', 'maybe','MStateDependence','strong','InitialStep',1e-10, 'MStateDependence','strong');
[tSol, ySol]=ode15s(@(t,y) SBR(t,y), tSpan, y0);
Warning: Failure at t=1.071763e+02. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (2.273737e-13) at time t.
%plotting
figure (1)
plot(tSol,ySol(:,7),'-r'); hold on %Sse
plot(tSol,ySol(:,8),'-g'); hold on %Sno3
plot(tSol,ySol(:,9),'--r'); hold on %Sno2
plot(tSol,ySol(:,10),'-k'); hold on %nh4
plot(tSol,ySol(:,11),'-b'); hold on %So2
figure (2)
plot(tSol,ySol(:,1),'-r'); hold on %Xaob
plot(tSol,ySol(:,2),'-g'); hold on %Xnb
plot(tSol,ySol(:,3),'--r'); hold on %Xnsp
plot(tSol,ySol(:,4),'-k'); hold on %amx
plot(tSol,ySol(:,5),'-b'); hold on %Xhan
plot(tSol,ySol(:,6),'--g'); hold off %Xs
xlabel('Time (days)'); ylabel('Concentration (g/m^3)');
%Main function
% SBR simulation in MATLAB
function dy=SBR(t,y)
maxSnh4=200;
%function fval=SBR1(t,y)
% Initialize arrays
Xaob=y(1);
Xnb = y(2);
Xnsp= y(3);
Xamx= y(4);
Xhan= y(5);
Xs= y(6);
Sse= y(7);
Sno3= y(8);
Sno2= y(9);
Snh4= y(10);
So2= y(11);
T=303;
%aob
O1=68/(0.00831*293*(273+T)); u1=0.054*exp(O1*(T-293)); Ko2aob=0.6; Knh4aob=2.4; Yaob=0.15; Kno3aob=0.5;
%nb
O2=44/(0.00831*293*(273+T)); u2=0.047*exp(O2*(T-293)); Ko2nb=2.2; Kno2nb=0.5*1.375; Ynb=0.041;Kno3nb=0.5;
%nsp
O3=44/(0.00831*293*(273+T)); u3=0.69*exp(O3*(T-293)); Ko2nsp=0.33; Kno2nsp=0.52; Ynsp=0.14;Kno3nsp=0.5;
%amx
O5=71/(0.00831*293*(273+T)); u5=0.0022*exp(O5*(T-293)); Kno2amx=0.05; Knh4amx=0.07; Yamx=0.159; Ko2amx=0.01; Kno3amx=0.5;
%han
O6=0.069;u6=0.5*7.2*exp((T-293)*0.069); KSsehan=2; Kno2han=0.5;Kno3han=0.5; YNo2han=0.3; Ko2han=0.2;
YNo3han=0.43; Yhaer=0.54;
%other constants
baob=0.0054*exp(O6*(T-293)); bnb=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; Qo=116.67; Xso=250;
Snh4o=650; Sno2o=0.2; So2o=0;
%oxygen transfer parameters
KaL=1920;
V=1400; SRT=100;
Xaobo=2; Xnbo=2; Xnspo=2; Xhano=2; Xamxo=2;
Vo=V-0.3*V; % volume in the reactor at the start of refilling
% SBR operation cycle
refilling_duration = 1.59; % Duration of refilling phase (min)
aeration_duration = 15; % Duration of aeration phase (min)
mixing_duration = 45; % Duration of mixing phase (min)
settling_duration = 30; % Duration of settling phase (min)
withdrawal_duration = 0.58; % Duration of withdrawal phase (min)
cycle_duration=refilling_duration+aeration_duration+mixing_duration+settling_duration+withdrawal_duration;
t0=0+cycle_duration;
t1=t0+15;
t2=t1+15;
t3=t2+45;
t4=t3+60;
t5=t4+15;
%refilling
%while t>t0 && t<=t1
if t>t0 && t<=t1
% input parameters
Sseo=300; Sno3o=0.2; Qo=220; Xso=250;
Snh4o=650; Sno2o=0.2; So2o=0; KaL=0;
Xaobo=2; Xnbo=2; Xnspo=2; Xhano=2; Xamxo=2;
%break
end
%aeration
%while t>t1 && t<=t2
if t>t1 && t<=t2
% input parameters
Sseo=0; Sno3o=0; Qo=0; Xso=0; Snh4o=0; Sno2o=0; So2o=0; Xaobo=0; Xnbo=0; Xnspo=0; Xhano=0; Xamxo=0;
if Snh4<10 || Sno2>10 || So2>0.8
KaL=0;
elseif Snh4>10 && So2<=0.5
KaL=1920;
elseif Sno2<10 && So2<0.5
KaL=1920;
end
%break
end
%mixing & settling
%while t>t2 && t<=t4
if t>t2 && t<=t4
% input parameters
Sseo=0; Sno3o=0; Qo=0; Xso=0; Snh4o=0; Sno2o=0.2; So2o=0; KaL=0; Xaobo=0; Xnbo=0; Xnspo=0; Xhano=0; Xamxo=0;
%break
end
%withdrawal
%while t>t4 && t<=t5
if t>t4 && t<=t5
% input parameters
Sseo=0; Sno3o=0; Qo=600; Xso=0; Snh4o=0; Sno2o=0.2; So2o=0; KaL=0; Xaobo=0; Xnbo=0; Xnspo=0; Xhano=0; Xamxo=0;
%break
end
fval(1)=((Vo*Xaobo-Vo*Xaob/SRT)/V + u1*(So2/(Ko2aob+So2))*(Snh4/(Knh4aob+Snh4))*Xaob - baob*(So2/(Ko2aob+So2))*Xaob - baob*naob*((Ko2aob/(Ko2aob+So2))*(Sno2+Sno3)/(Kno3aob+Sno2+Sno3))*Xaob);
fval(2)=((Vo*Xnbo-Vo*Xnb/SRT)/V + u2*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2))*Xnb - bnb*(So2/(Ko2nb+So2))*Xnb - bnb*nnb*((Ko2nb/(Ko2nb+So2))*(Sno2+Sno3)/(Kno3nb+Sno2+Sno3))*Xnb);
fval(3)=((Vo*Xnspo-Vo*Xnsp/SRT)/V + u3*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2))*Xnsp - bnsp*(So2/(Ko2nsp+So2))*Xnsp - bnsp*nnsp*((Ko2nsp/(Ko2nsp+So2))*(Sno2+Sno3)/(Kno3nsp+Sno2+Sno3))*Xnsp);
fval(4)=((Vo*Xamxo-Vo*Xamx/SRT)/V + u5*((Snh4/(Knh4amx+Snh4))*(Sno2/(Kno2amx+Sno2))*(Ko2amx/(Ko2amx+So2)))*Xamx-bamx*(So2/(Ko2amx+So2))*Xamx - bamx*namx*((Ko2amx/(Ko2amx+So2))*(Sno2+Sno3)/(Kno3amx+Sno2+Sno3))*Xamx);
fval(5)=((Vo*Xhano-Vo*Xhan/SRT)/V + (u6*(Sse/(KSsehan+Sse))*(Sno2/(Kno2han+Sno2))*(Ko2han/(Ko2han+So2)) + u6*(Sse/(KSsehan+Sse))*(Sno3/(Kno3han+Sno3))*(Ko2han/(Ko2han+So2))+u6*((Sse/(KSsehan+Sse))*(So2/(Ko2han+So2)))*Xhan)- bhan*(So2/(Ko2han+So2))*Xhan - bhan*nhan*((Ko2han/(Ko2han+So2))*(Sno2+Sno3)/(Kno3han+Sno2+Sno3))*Xhan);
fval(6)=((Qo*Xso-Qo*Xs)/V - KH*((Xs/Xhan)/(KX+(Xs/Xhan)))*Xhan);
fval(7)=((Qo*Sseo-Qo*Sse)/V - (u6*nhan*(1/YNo2han)*(Sse/(KSsehan+Sse))*(Sno2/(Kno2han+Sno2))*(Ko2han/(Ko2han+So2)))*Xhan -(u6*nhan*(1/YNo3han)*(Sse/(KSsehan+Sse))*(Sno3/(Kno3han+Sno3))*(Ko2han/(Ko2han+So2)))*Xhan - ((u6/Yhaer)*(Sse/(KSsehan+Sse))*(So2/(Ko2han+So2))*Xhan) + KH*((Xs/Xhan)/(KX+(Xs/Xhan)))*Xhan);
fval(8)=((Qo*Sno3o-Qo*Sno3)/V - (u6*nhan*((1-YNo3han)/(2.86*YNo3han))*(Sse/(KSsehan+Sse))*(Sno3/(Kno3han+Sno3))*(Ko2han/(Ko2han+So2)))*Xhan + (u5/1.14)*(Snh4/(Knh4amx+Snh4))*(Sno2/(Kno2amx+Sno2))*(Ko2amx/(Ko2amx+So2))*Xamx + (u2/Ynb*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2))*Xnb) + ((u3/Ynsp)*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2))*Xnsp) - ((1-fI)/2.86)*baob*naob*(Ko2aob/(Ko2aob+So2))*(Sno2+Sno3)/(Kno3aob+Sno2+Sno3)*Xaob - ((1-fI)/2.86)*bnb*nnb*(Ko2nb/(Ko2nb+So2))*(Sno2+Sno3)/(Kno3nb+Sno2+Sno3)*Xnb - ((1-fI)/2.86)*bnsp*nnsp*(Ko2nsp/(Ko2nsp+So2))*(Sno2+Sno3)/(Kno3nsp+Sno2+Sno3)*Xnsp - ((1-fI)/2.86)*bamx*namx*(Ko2amx/(Ko2amx+So2))*(Sno2+Sno3)/(Kno3amx+Sno2+Sno3)*Xamx-((1-fI)/2.86)*bhan*nhan*(Ko2han/(Ko2han+So2))*(Sno2+Sno3)/(Kno3han+Sno2+Sno3)*Xhan);
fval(9)=((Qo*Sno2o-Qo*Sno2)/V - (((1-YNo2han)/(1.71*YNo2han))*u6*nhan*(Sse/(KSsehan+Sse))*(Sno2/(Kno2han+Sno2))*(Ko2han/(Ko2han+So2))*Xhan) - (((u5/Yamx)+(u5/1.14))*(Snh4/(Knh4amx+Snh4))*(Sno2/(Kno2amx+Sno2))*(Ko2amx/(Ko2amx+So2)))*Xamx + ((u1/Yaob)*(So2/(Ko2aob+So2))*(Snh4/(Knh4aob+Snh4)))*Xaob - ((u2/Ynb)*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2))*Xnb) - ((u3/Ynsp)*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2)))*Xnsp);
fval(10)=((Qo*Snh4o-Qo*Snh4)/V - (u1*(INBM+1/Yaob)*(So2/(Ko2aob+So2))*(Snh4/(Knh4aob+Snh4)))*Xaob -(u5*(INBM+1/Yamx)*(Snh4/(Knh4amx+Snh4))*(Sno2/(Kno2amx+Sno2))*(Ko2amx/(Ko2amx+So2)))*Xamx - (u2*INBM*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2)))*Xnb -(u3*INBM*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2)))*Xnsp - INBM*u6*nhan*((Sse/(KSsehan+Sse))*(Sno2/(Kno2han+Sno2))*(Ko2han/(Ko2han+So2)))*Xhan - u6*INBM*nhan*((Sse/(KSsehan+Sse))*(Sno3/(Kno3han+Sno3))*(Ko2han/(Ko2han+So2)))*Xhan - u6*INBM*((Sse/(KSsehan+Sse))*(So2/(Ko2han+So2)))*Xhan + (INBM-(fI*INXI))*baob*naob*((Ko2aob/(Ko2aob+So2))*(Sno2+Sno3)/(Kno3aob+Sno2+Sno3))*Xaob + (INBM-(fI*INXI))*bnb*nnb*((Ko2nb/(Ko2nb+So2))*(Sno2+Sno3)/(Kno3nsp+Sno2+Sno3))*Xnb +(INBM-(fI*INXI))*bnsp*nnsp*((Ko2nsp/(Ko2nsp+So2))*(Sno2+Sno3)/(Kno3nsp+Sno2+Sno3))*Xnsp + (INBM-(fI*INXI))*bamx*namx*((Ko2amx/(Ko2amx+So2))*(Sno2+Sno3)/(Kno3amx+Sno2+Sno3))*Xamx + (INBM-(fI*INXI))*bhan*nhan*((Ko2han/(Ko2han+So2))*(Sno2+Sno3)/(Kno3han+Sno2+Sno3))*Xhan +(INBM-(fI*INXI))*baob*(So2/(Ko2aob+So2))*Xaob + (INBM-(fI*INXI))*bnb*(So2/(Ko2nb+So2))*Xnb +(INBM-(fI*INXI))*bnsp*(So2/(Ko2nsp+So2))*Xnsp + (INBM-(fI*INXI))*bamx*(So2/(Ko2amx+So2))*Xamx + (INBM-(fI*INXI))*bhan*(So2/(Ko2han+So2))*Xhan);
fval(11)=((Qo*So2o-Qo*So2)/V + KaL*(8-So2) -(((1-Yhaer)/Yhaer)*u6*(Sse/(KSsehan+Sse))*(So2/(Ko2han+So2)))*Xhan -(((3.43-Yaob)/Yaob)*u1*(So2/(Ko2aob+So2))*(Snh4/(Knh4aob+Snh4)))*Xaob - (((1.14-Ynb)/Ynb)*u2*(Sno2/(Kno2nb+Sno2))*(So2/(Ko2nb+So2)))*Xnb - (((1.14-Ynsp)/Ynsp)*u3*(Sno2/(Kno2nsp+Sno2))*(So2/(Ko2nsp+So2)))*Xnsp - (1-fI)*baob*(So2/(Ko2aob+So2))*Xaob - (1-fI)*bnb*(So2/(Ko2nb+So2))*Xnb -(1-fI)*bnsp*(So2/(Ko2nsp+So2))*Xnsp - (1-fI)*bamx*(So2/(Ko2amx+So2))*Xamx - (1-fI)*bhan*(So2/(Ko2han+So2))*Xhan);
if fval(10)>maxSnh4
fval(10)=maxSnh4;
%else
% dy=fval(10);
end
dy = fval.';
end
%end
  2 个评论
KIPROTICH KOSGEY
KIPROTICH KOSGEY 2023-7-13
Thank you for the insights. I will re-structure it further. However, I have coded for fval(10)<=200 but it is still going over this value as is evident in figure(1).
Please advise on this on this.
Torsten
Torsten 2023-7-13
编辑：Torsten 2023-7-13
The figure doesn't display "fval", but "y".
And the mass matrix must be of size 11x11, not 55x55.
But since you don't use the options settings, MATLAB does not throw an error.

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Programming 的更多信息

标签

产品


版本

R2023a

Community Treasure Hunt

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

Start Hunting!

Translated by