% main_reak_cl.m
% Lukket sløyfe simulering av reaktor med PI-regulatror.
% Skisse til loesning av oving 7, punkt f).

clear all

% Prosessparametre.
x1s=2.5;x2s=1; ths=10; us=25;
k1=50; k2=100; k3=10;
d=[0,1];

% Tidshorisont og samplingsintervall
dt=0.001;                 % Samplingsintervall
t=0:dt:1; t=t';
N=length(t);

% Prosess og maalestoey.
w=randn(N,1)*0.0001;
v=randn(N,2)*0.1;

% Referansesignal
dR=prbs1(N,30,150);
R=ones(N,1)+0.1*dR;

% regulatorparametre
Kp=50; Ti=1/75; 

z=25;        % Initialverdi paa regulatortilstanden
x=[x1s;x2s]; % Initialverdi for tilstandsvektoren

for k=1:N
   % "Henter" maalingen fra prosessen
   y=d*x+w(k);
   
   % PI-regulator.
   e=R(k)-y;
   u=Kp*e+z;
   z=z+dt*Kp*e/Ti;
   
   % Lagrer variable
   X(k,:)=x';
   Y(k,:)=y';
   U(k,:)=u';
  
   % Integrerer prosessmodellen
   fx1=-k1*x(1)-k3*x(1)^2+(ths-x(1))*u+v(k,1);
   fx2=k1*x(1)-k2*x(2)-x(2)*u +v(k,2);
   x=x+dt*[fx1;fx2];   
   
end

% Identifikasjon vha DSR
L=3; J=4;
[a,b,d,e,c,f,x0]=dsr(Y,U,L,0,J,1,3);
ym=dlsim(a,b,d,e,U,x0);

figure(1)
subplot(411), plot(t,R),
title('Referansesignal: r_k'), grid
subplot(412), plot(t,U)
title('Paadragssignal: u_k'), grid
subplot(413), plot(t,[Y ym])
title('Utgang og estimert utgang'), grid

