Simultaneous differential equations with conditional change of a coefficient

Question

I want to model a Chemostat (a certain kind of bioreactor). I setup the following system which can be solved with the ode45 method:

%chemostat model, based on:
%DCc=-v0*Cc/V + umax*Cs*Cc/(Ks+Cs)-rd -->Change of cell concentration over time
%Dcs=(v0/V)*(Cs0-Cs) - Cc*(Ys*umax*Cs/(Ks+Cs)-m) -->Change of substrate concentration over time

function dydt=sys(t,y,v0,V,umax,Ks,rd,Cs0,Ys,m)
         dydt=[-(v0/V)*y(1)+(umax*y(1)*y(2))/(Ks+y(2))-rd; 
               (v0/V)*(Cs0-y(2))-(Ys*umax*y(2)*y(1))/(Ks+y(2))];

I call the function with:

[t,y]=ode45(@systemEquations, [0 40],[1 100],[],**v0**, V,umax,Ks,rd,Cs0,Ys,m);

The values of all the additional coefficients are set before running the calculation. So far everything works. Now I would like v0 to be dependent on the system state. E.g. initially v0 = 0 and when my cell concentration surpasses a certain value I want it to change.

The problem is, I absolutely have no idea how to implement this. The Matlab help of the ode solver wasn't of any help either...

Thanks for any help or suggestions!

Cheers, dahlai

Answer 1

Use your Differential System function to take care of it for you! First formulate a v0 as a function of time and state:

function v0 = funV0( t, y )
    v0 = 1;
    if y(1) > 5
        v0 = y(2); % As an example, if [Y1] > 5, then set v0 = [Y2]
    end
end

Then utilize that function in your differential system. You can either pass @funV0 as a parameter, or let your differential system function execute it directly. I will treat it as a parameter for now.

function dydt=sys(t,y,v0,V,umax,Ks,rd,Cs0,Ys,m)
     dydt=[-(v0(t,y)/V)*y(1)+(umax*y(1)*y(2))/(Ks+y(2))-rd; 
           (v0(t,y)/V)*(Cs0-y(2))-(Ys*umax*y(2)*y(1))/(Ks+y(2))];

Called using:

[t,y]=ode45(@systemEquations, [0 40],[1 100],[],@funV0, V,umax,Ks,rd,Cs0,Ys,m);

Of course, this principle can be applied to any time or state dependent parameter. This style also allows you to try various parameter functions Let's say you are approximating the parameter as a function of state: You can trial out different approximations by passing in different V0 functions.