Figure 4.9 (page 122):

Batch reactor with nth-order kinetics, r=kc_A^n, k_0=kc_{A0}^{n-1}=1, n<1.

Code for Figure 4.9

Text of the GNU GPL.

main.m

%% Copyright (C) 2001, James B. Rawlings and John G. Ekerdt
%%
%% This program is free software; you can redistribute it and/or
%% modify it under the terms of the GNU General Public License as
%% published by the Free Software Foundation; either version 2, or (at
%% your option) any later version.
%%
%% This program is distributed in the hope that it will be useful, but
%% WITHOUT ANY WARRANTY; without even the implied warranty of
%% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
%% General Public License for more details.
%%
%% You should have received a copy of the GNU General Public License
%% along with this program; see the file COPYING.  If not, write to
%% the Free Software Foundation, 59 Temple Place - Suite 330, Boston,
%% MA 02111-1307, USA.

npts = 100;
tfin = 10;
K = 1;
orders = [-2;-1;-0.5;0;0.5;1;2];
norders = length(orders);
store = [];
for i = 1: norders
  t = linspace(0,tfin,npts)';
  n=orders(i);
  if (n ==1 )
    c = exp(-K*t);
  else
    %% test here for nonnegativity of c
    if (n < 1)
      tstop =  1/((-n+1)*K);
      t = linspace(0,tstop,npts)';
    end
    arg = (1 - (-n+1)*K*t);
    c = arg.^(1/(-n+1));
  end
  store = [store t c];
end
plot (store(:,1), store(:,2), store(:,3), store(:,4), ...
      store(:,5), store(:,6), store(:,7), store(:,8), ...
      store(:,9), store(:,10), store(:,11), store(:,12), ...
      store(:,13), store(:,14));
title ('Figure 4.9')