Plotting Kinetic Model Monod, Haldane With Parameters

Hello, I want to ask some questions regarding my code (attacthed with name CobaParZhang.m) about plotting graph. When I run the code, it shows error "Index exceeds the number of array elements (1)" like attached bellow.

I would like to plot C(1) (biomass), C(2) (substrate), and C(3) (product) into one figure from these 3 equations below:

, where dX/dt = miu*X

And here are the kinetic parameters:

I am wondering maybe it has some problem with the equation that I am using or the inconsistent unit (biomass and substrate is in g/L, product is in mg/L). Any help or comment or advice about this problem is very welcomed. Thank you so much for your attention.

function CobaParameterZhang
Rentang = linspace(0, 20, 50);
C0 = [0.143; 0.9659; 0.094406]; 

Miumax=0.5258;
Ks=0.0211;
Kxi=56.6813;
KXI=53.26;
CXm=2.9161;
CPm=55.6192;

Miupx=6.0679;
YPX=20.7502;
KPS=0.0017;
Kpi=63.9593;
KPI=72.5525;
KP=19.1078;

Miusx=1.8216;
YXS=2.6872;
YPS=60.3616;

[t,C]=ode23(@(t,C)ParameterZhang(t,C,Miumax,Ks,Kxi,KXI,CXm,CPm,Miupx,YPX,KPS,Kpi,KPI,KP,Miusx,YXS,YPS),Rentang,C0);

figure
yyaxis left
plot(t,C(:,[1 2]))
ylabel('Concentration (g/L)')
yyaxis right
plot(t, C(:,3))
ylabel('Concentration (\mu g/L')
grid
legend('C_1','C_2','C_3', 'Location','W')
end
function dCdt=ParameterZhang(t,C,Miumax,Ks,Kxi,KXI,CXm,CPm,Miupx,YPX,KPS,Kpi,KPI,KP,Miusx,YXS,YPS)
        z = [1
             2
             3
             4
             5
             6
             7
             8
             9
             10
             11
             12
             13
             14
             15
             16
             17
             18
             19
             20];
         
        j = 1;
        for j = 1 : 20   
        if z(j) <= 4 
            I = 45; 
        elseif z(j) <= 7
            I = 90;
        elseif z(j) <= 11
            I = 120;
        elseif z(j) <= 15
            I = 150;
        else 
            I = 180;
        end
mu = ((Miumax*C(2))/(Ks+C(2)+(C(2)^2/Kxi)))*(I/(KXI+I))*(1-C(1)/CXm)*(1-C(3)/CPm);
dCdt(1,:) = mu*C(1);
dCdt(2,:) = (-(1/YXS)*dCdt(1))+(-(1/YPS)*dCdt(3))+(-Miusx*C(1));
dCdt(3,:) = (YPX*dCdt(1)+Miupx*C(1))*(C(2)/(KPS+C(2)+C(2)^2/Kpi))*(C(3)/(KP+C(3)))*(I/(KPI+I));
        end
end

Answer 1

[Star Strider]

The problem is that in this assignment —

dCdt(2,:) = (-(1/YXS)*dCdt(1))+(-(1/YPS)*dCdt(3))+(-Miusx*C(1));

‘dCdt(3)’ has not yet been created, and only ‘dCdt(1)’ exists in the code at that point.

It will llikely be necessary to re-code the differential equations.

Also, ‘I’ is not a function of anything useful in the original code. I created an anonymous function to calculate it correctly as a function of ‘t’ rather than ‘j’ that does not otherwise exist in the code.

CobaParameterZhang

Index exceeds the number of array elements. Index must not exceed 1.

Error in solution>CobaParameterZhang/ParameterZhang (line 44)
dCdt(2,:) = (-(1/YXS)*dCdt(1))+(-(1/YPS)*dCdt(3))+(-Miusx*C(1));

Error in solution>@(t,C)ParameterZhang(t,C,Miumax,Ks,Kxi,KXI,CXm,CPm,Miupx,YPX,KPS,Kpi,KPI,KP,Miusx,YXS,YPS) (line 25)
[t,C]=ode23(@(t,C)ParameterZhang(t,C,Miumax,Ks,Kxi,KXI,CXm,CPm,Miupx,YPX,KPS,Kpi,KPI,KP,Miusx,YXS,YPS),Rentang,C0);

Error in odearguments (line 90)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.

Error in ode23 (line 105)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);

Error in solution>CobaParameterZhang (line 25)
[t,C]=ode23(@(t,C)ParameterZhang(t,C,Miumax,Ks,Kxi,KXI,CXm,CPm,Miupx,YPX,KPS,Kpi,KPI,KP,Miusx,YXS,YPS),Rentang,C0);

function CobaParameterZhang
Rentang = linspace(0, 20, 50);
C0 = [0.143; 0.9659; 0.094406]; 

Miumax=0.5258;
Ks=0.0211;
Kxi=56.6813;
KXI=53.26;
CXm=2.9161;
CPm=55.6192;

Miupx=6.0679;
YPX=20.7502;
KPS=0.0017;
Kpi=63.9593;
KPI=72.5525;
KP=19.1078;

Miusx=1.8216;
YXS=2.6872;
YPS=60.3616;

[t,C]=ode23(@(t,C)ParameterZhang(t,C,Miumax,Ks,Kxi,KXI,CXm,CPm,Miupx,YPX,KPS,Kpi,KPI,KP,Miusx,YXS,YPS),Rentang,C0);

figure
yyaxis left
plot(t,C(:,[1 2]))
ylabel('Concentration (g/L)')
yyaxis right
plot(t, C(:,3))
ylabel('Concentration (\mu g/L')
grid
legend('C_1','C_2','C_3', 'Location','W')
% end

function dCdt=ParameterZhang(t,C,Miumax,Ks,Kxi,KXI,CXm,CPm,Miupx,YPX,KPS,Kpi,KPI,KP,Miusx,YXS,YPS)

I = @(t) (45).*((t>=0)&(t<4)) + (90).*((t>=4)&(t<7)) + (120).*((t>=7)&(t<11)) + (150).*((t>=11)&(t<15)) + (180).*(t>=15);

    mu = ((Miumax*C(2))/(Ks+C(2)+(C(2)^2/Kxi)))*(I(t)/(KXI+I(t)))*(1-C(1)/CXm)*(1-C(3)/CPm);
    dCdt(1,:) = mu*C(1);
    dCdt(2,:) = (-(1/YXS)*dCdt(1))+(-(1/YPS)*dCdt(3))+(-Miusx*C(1));
    dCdt(3,:) = (YPX*dCdt(1)+Miupx*C(1))*(C(2)/(KPS+C(2)+C(2)^2/Kpi))*(C(3)/(KP+C(3)))*(I(t)/(KPI+I(t)));
end

end

Also, I have seen this before, and very recently, posted by a different person!

.

[Regina Maria]

Hello @Star Strider, thank you so much for the answer and the suggestion for "I" function as a function of time, I really appreciate it. For right now, I think it is best for me to change the differential equation to make the code works.

About this post, actually I have posted this 4 days ago, but I reposted again to make the question simpler and with a new title. And also, yes there is a similar question posted by a different person with a different equations and parameters, but I cannot make mine work.. so I decided to ask again, sorry if it is bothering you.

Wish you all the best in life and stay healthy!

[Star Strider]

As always, my pleasure!

You are not bothering me! I apologise if I gave you that impression.

I am simply suggesting that the code for the differential equations needs a bit of work to get it running correctly, and I described the problem.