
function DYDT=f(t,y)
% This is file cstr3.m
% INPUT: Time t and state vector y
% OUTPUT: derivatives DYDT
% Usage with ode15s:   [T,Y]=  ode15s('cstr3',[0 5],[1.0 0.9 8.1])

% I. Data (parameters and independent variables)
cAf=10; cBf=0; cCf=0;  % inlet concentrations
V = 0.9;               % reactor volume
q = 0.1*1.2;           % 20% increase in q
k1 = 1; k2 = 1;        % rate constants

% II. Extract present value of states
cA=y(1);
cB=y(2);
cC=y(3);

% III. Evaluate derivatives of states
f1= (q*cAf - q*cA - k1*cA*V)           /V;
f2= (q*cBf - q*cB + k1*cA*V - k2*cB*V) /V;
f3= (q*cCf - q*cC           + k2*cB*V) /V;
DYDT=[f1; f2; f3;];


