% set up global variables for use in the DE routine.
global q uy uz Vmax K M eps delta

t0  = 0;     %initial time
t1  = 40;    % final time
z0  = 0;     % initial ?? conc
y0  = 0;     % initial nitrogen conc
x0  = 4.9;   % initial yeast concentration

% initial state for use in DE solver.

S0  = [x0;y0;z0];

q   = 0.15; % dilution rate
eps = q/1.1;

Vmax  =  12.8909;
K    =   0.0499;
M    =   5.7694;
delta =   0.0592;
uz    = 134.6042;

uy    =  68;

[t,St] = ode45('KRA2nchemostatdynamics',[t0,t1],S0);
    
subplot(311)
plot(t,St(:,1)),title('Biomass')

subplot(312)
plot(t,St(:,2)),title('nitrogen')

subplot(313)
plot(t,St(:,3)),title('nutrient 2')