%
% Example B.1: Carbon Removal Process
%
% - constant volume
% - well mixed
% - Monod kinetics
% - no biomass in effluent
% - no clarifier dynamics
% - same substrate in effluent, recycle & wastage
% - no oxygen in recycle & wastage
%
% Bob Newell, March 1996
%
%-----------------------------------------------------------
% parameters
muhat = 3.733;
Ks = 20;
Koh = 0.2;
Yh = 0.6;
bh = 0.4;
fp = 0.1;
a = 166;
b = 80;
So_sat = 10;
%
global p
p = [ muhat Ks Koh Yh bh fp a b So_sat ];
%-----------------------------------------------------------
% initial states
V = 12.5;
Ss = 4.4;
So = 3;
Xh = 850;
%
xo = [ V Ss So Xh ];
%-----------------------------------------------------------
% integrate the state equation
to = clock;
[ t, x ] = ode45( 'Carbon_Removal', [ 4/24 10/24 ], xo );
secs = etime( clock, to )
%-----------------------------------------------------------
% evaluate the outputs
y = Carbon_Removal_Outputs( t, x );
%-----------------------------------------------------------
% plot the results
figure;
[ ax, h1, h2 ] = plotyy( 24*t, y(:,6), 24*t, y(:,8) );
axesyy( ax, [ 4 10 0 10 0 1200 ], ...
   'Time in hours', ...
   'Nutrients in mg/l', ...
   'Biomass in mg/l', 2, ...
   [ 0 2 4 6 ], [ 800 1000 1200 ] );
set(h2,'linestyle','--');
hold on;
plot( 24*t, y(:,7), ':' );
hold off;
h = text( 8.7, 7.7, 'Xh' );
set( h, 'fontsize', 24 );
h = text( 8.7, 6.3, 'Ss' );
set( h, 'fontsize', 24 );
h = text( 8.7, 2.3, 'So' );
set( h, 'fontsize', 24 );
%
printsto -f1 -deps2 figb-2a.eps
%
%-----------------------------------------------------------
% integrate the state equation
to = clock;
days = 5;
[ t, x ] = ode45( 'Carbon_Removal', [ 0 days ], xo );
secs = etime( clock, to )
%-----------------------------------------------------------
% evaluate the outputs
y = Carbon_Removal_Outputs( t, x );
%-----------------------------------------------------------
% plot the results
figure;
[ ax, h1, h2 ] = plotyy( t, y(:,6), t, y(:,8) );
axesyy( ax, [ 0 days 0 10 0 1200 ], ...
   'Time in days', ...
   'Nutrients in mg/l', ...
   'Biomass in mg/l', 2, ...
   [ 0 2 4 6 ], [ 800 1000 1200 ] );
set(h2,'linestyle','--');
hold on;
plot( t, y(:,7), ':' );
hold off;
h = text( 4, 8.3, 'Xh' );
set( h, 'fontsize', 24 );
h = text( 4, 5.7, 'Ss' );
set( h, 'fontsize', 24 );
h = text( 4, 2, 'So' );
set( h, 'fontsize', 24 );
%
printsto -f2 -deps2 figb-2b.eps
%-----------------------------------------------------------
% the end!   
