% programa 5.4, W.L.Luyben, pag. 123
% "three isothermal cstr (closedloop)" 
% modif. E.Quintero.M 18nov'98
%------------------------------------------------------
% (with proportional feedback controller)
% disturbance (cad) is a step change at time equal zero
% from 0 to 0.2
% time is in minutes
clear
format short e
% initial conditions
time=0; ca1=0.4; ca2=0.2; 
ca3=0.1; ca3set=0.1; erint=0;
% disturbance
cad=0.2;
% parameter values
tau=2; k=0.5; kc=30; taui=5; 
delta=0.01; tprint=0; xx=[];
while(time<=10)
   % feedback controller
   e=ca3set-ca3;
   cam=0.8+kc*(e+erint/taui);
   cao=cam+cad;
   % evaluate derivatives
   ca1dot=(cao-ca1)/tau-k*ca1;
   ca2dot=(ca1-ca2)/tau-k*ca2;
   ca3dot=(ca2-ca3)/tau-k*ca3;
   if(time>=tprint)
      xx=[xx;time,ca1,ca2,ca3,cam];
      tprint=tprint+0.1;
   end
   ca1=ca1+ca1dot*delta;
   ca2=ca2+ca2dot*delta;
   ca3=ca3+ca3dot*delta;
   erint=erint+e*delta;
   time=time+delta;
end
xx
subplot(2,1,1),plot(xx(:,1),xx(:,5)),ylabel('Cam'),grid
subplot(2,1,2),plot(xx(:,1),xx(:,4)),ylabel('Ca3'),grid
xlabel('Time, min')