function xdot = vanvusse(t,x);
%  
%  (c) b.w. bequette
%  25 July 94
%  revised 20 July 96
%
%  dynamic model for the van de vusse reaction
%
%  to use this function file with ode45:
%
%   [t,x] = ode45('vanvusse',t0,tf,x0]
%
%  where t0 is the initial time (usually 0)
%        tf is the final time
%        x0 is the initial condition for the states
%
%  x    = states
%  x(1) = concentration of A
%  x(2) = concentration of B
%  t    = time
%
%  fov  = dilution rate
%  ca   = concentration of a
%  cb   = concentration of b
%  caf  = feed concentration of a
%
%  the reaction scheme is
%
%  A ---> B ---> C
%    (k1)   (k2)
% 2A ---> D
%    (k3)
%
   k1 = 5/6; % rate constant for A --> B (min^-1)
   k2 = 5/3; % rate constant for B --> C (min^-1)
   k3 = 1/6; % rate constant for A + A --> D (mol/(liter min)
   caf = 10; % feed concentration of A, mol/liter
   fov = 4/7;% flowrate/volume (min^-1)
%
%  for step changes in fov, use the following
%   fov  = 4/7 + delf;
% 
%  use familiar notation for the states
%
   ca   = x(1);
   cb   = x(2);
%
%  modeling equations
%
   dcadt = fov*(caf - ca) -k1*ca - k3*ca*ca;
   dcbdt = -fov*cb + k1*ca - k2*cb;
%
%  derivatives placed in vector notation
%
   xdot = [dcadt;dcbdt];