% programa 17.3, luyben, pag. 607
% "TLC for Ogunnaike and Ray column"
% modif. E.Quintero.M 23nov'98
%---------------------------------------------------------
clear
format short e
% GM process transfer function data
kp=[0.66, -0.61, -0.0049; 1.11, -2.36, -0.012; ...
    -34.68, 46.2, 0.87]; 
tau1=[0, 0, 0; 0, 0, 0; 0, 0, 11.61];
tau2=[6.7, 8.64, 9.06; 3.25, 5.0, 7.09; 8.15, 10.9, 3.89];
tau3=[0, 0, 0; 0, 0, 0; 0, 0, 18.8];
D=[2.6, 3.5, 1; 6.5, 3, 1.2; 9.2, 9.4, 1];
% GL process transfer function data
kl=[0.14, 0.53, -11.54];
taul=[6.2, 6.9, 7.01];
Dl=[12, 10.5, 0.6];
% controller data - empirical settings
%kc=[1.2, -0.15, 0.6];
%reset=[5, 10, 4];
% controller data - BLT
kc=[1.51, -0.29, 2.63];
reset=[16.38, 17.95, 6.61];
n=3; xx=[];
% subroutine TLC
w=0.01;
dw=10^(1/10);
while(w<=10.01)
   s=i*w;
   for j=1:n
      for k=1:n
         gm(j,k)=kp(j,k)*exp(-D(j,k)*s)*(tau1(j,k)*s+1)/ ...
         ((tau2(j,k)*s+1)*(tau3(j,k)*s+1));
      end
   end
   b=zeros(n);
   for j=1:n
      b(j,j)=kc(j)*(reset(j)*s+1)/(reset(j)*s); 
   end
   q=gm*b;
   h=diag(ones(n,1),0)+q;
   hin=inv(h);
   for j=1:n    
      gl(j)=kl(j)*exp(-Dl(j)*s)/(taul(j)*s+1);  
   end
   xtlc=hin*gl';
   db=20*log10(abs(xtlc));
   xx=[xx;w db'];
   w=w*dw;
end
subplot(2,2,1),semilogx(xx(:,1),xx(:,2)),grid
subplot(2,2,2),semilogx(xx(:,1),xx(:,3)),grid
subplot(2,2,3),semilogx(xx(:,1),xx(:,4)),grid
