How to prevent substrate concentration going below minimum limit in fed batch fermentation model.

clear all;
clc;
close all;

tic
global miu_max K_S Y_SX

% C.I.
S0 = 300;
X0 = 105;
V0 = 1936;

% miu_max   =  maximum specific growth, time^-1
% K_S       =  substrate affinity constant, [S]
% Y_XS      =  substrate yield in biomass, kg/kg

% Kinetic parameters values

miu_max = 0.44;
K_S = 0.180;
Y_SX = 1/0.48;

tin = 0;
tfinal = 20;

time_interval = 1;

num_iterations = tfinal / time_interval;

all_t = [];
all_y = [];


y0 = [S0; X0; V0];

final_biomass_target = 3200;

while true
    
    tspan = [tin, tin + time_interval];

    [t, y] = ode23s(@Fb_function, tspan, y0);

    % Extract substrate concentration and biomass
    S_conc = y(:, 1);

    % Set negative substrate concentrations to zero
    S_conc(S_conc < 0) = 0;

    b_mass = y(:, 2);
    b_mass(:, 1) = (y(:, 2) * V0) / 1000;

    all_t = [all_t; t];
    all_y = [all_y; y];

    % Check if substrate concentration is below 100
    if any(S_conc < 100)
        % Add substrate to bring the concentration back to 200
        y0 = [200; y(end, 2); V0];
    else
        y0 = [S_conc(end); y(end, 2); V0];
    end

    % Check if the final biomass target is reached
    if all_y(end, 2) >= final_biomass_target
        break;
    end

    tin = tin + time_interval;
end

figure
yyaxis left
plot(all_t, all_y(:, 1));
ylabel('Substrate Concentration, g/l');

yyaxis right
plot(all_t, all_y(:, 2));
ylabel('Biomass, kg');

xlabel('Time, h');
legend('Substrate Concentration', 'Biomass');

function dy = Fb_function(t, Z)

global miu_max K_S Y_SX miu

dy = zeros(length(Z), 1);

% State variables
S = Z(1);
X = Z(2);

% Calculate growth rate and rates of production
miu = miu_max * (S / (K_S + S));
rX = miu * X;
rS = Y_SX * miu * X;

% ODEs
dy(1) = -rS;
dy(2) = rX;
end

Answer 1

[Torsten]

You can use the "Events" option of the ode integrators to interrupt and restart integration with different settings if a certain condition is met:

https://uk.mathworks.com/help/matlab/math/ode-event-location.html

I don't know if this is what you are looking for.