% LEP 9-5 : Bacterial growth in a Batch Reactor
%

tspan = [0 12]; % Range for the independent variable time
y0 = [1; 250; 0]; % Initial values for the dependent variables 
kd=0.01;
Ysc=12.5;
Ypc=5.6;
Ks=33.5;
m=0.03;
umax=0.46;
warning('off','all')
Cso=250;
[t,y]=ode45(@(t,y)ODEfun(t,y,kd,Ysc,Ypc,Ks,m,umax),tspan,y0);
tiledlayout(3,2)
set (gcf,'Position',[0,0.1,600,400])
nexttile
plot(t,y(:,1),'LineWidth',2)
leg=legend('C_c');
leg.ItemTokenSize=[10,10];
xlabel('t(hr)')
ylabel('C_c (g/dm^3)')
axis([0 12 0 20])
grid on
nexttile
plot(t,y(:,2),t,y(:,3),'LineWidth',2)
leg=legend('C_S','C_P');
leg.ItemTokenSize=[10,10];
xlabel('t(hr)')
ylabel('C (g/dm^3)')
axis([0 12 0 250])
grid on
z=size(y);
% Calculating values of selectivity
for i=1:z(1,1)
rsm(i)=m*y(i,1);
kobs(i)=(umax*(1-y(i,3)/93)^0.52);
rg(i)=kobs(i)*y(i,1)* y(i,2)./(Ks+y(i,2)); 
rd(i)=kd*y(i,1);
Xs(i)=(Cso-y(i,2))/Cso;
end
nexttile
plot(t,kobs(1:z(1,1)),'LineWidth',2)
axis([0 12 0 0.6])
xlabel('t(hr)')
ylabel('k_{obs} (h^-1)')
leg=legend('k_{obs}');
leg.ItemTokenSize=[10,10];
grid on
nexttile
plot(t,rg(1:z(1,1)),t,rd(1:z(1,1)),t,rsm(1:z(1,1)),'LineWidth',2)
axis([0 12 0 3])
xlabel('t(hr)')
ylabel('Rates (g/dm^3.hr)')
leg=legend('r_g','r_d','r_{sm}');
leg.ItemTokenSize=[10,10];
grid on
nexttile
plot(t,Xs(1:z(1,1)),'LineWidth',2)
axis([0 12 0 1])
xlabel('t(hr)')
ylabel('Conversion')
leg=legend('X_{S}');
leg.ItemTokenSize=[10,10];
grid on

%

function dYfuncvecdt = ODEfun(t,Yfuncvec,kd,Ysc,Ypc,Ks,m,umax); 
Cc = Yfuncvec(1); 
Cs = Yfuncvec(2); 
Cp = Yfuncvec(3); 

% Explicit equations
rd=Cc*kd;
rsm=m*Cc;
kobs=(umax*(1-Cp/93)^0.52);
rg=kobs*Cc* Cs/(Ks+Cs);

% Differential equations
dCcdt = rg-rd; 
dCsdt = Ysc *(-rg)-rsm; 
dCpdt = rg * Ypc; 
dYfuncvecdt = [dCcdt; dCsdt; dCpdt]; 
end

%