% Simulate a digital PI controller
% Plant = Gs
% Sampling Rate = T

% Plant
Gs = zpk([],[-1,-5,-10],200);

% Sampling Rate
T = 0.5;

dt = T/100

Gs_ABCD = ss(Gs);
K = 0.7;
A = Gs_ABCD.A;
B = Gs_ABCD.B;
C = Gs_ABCD.C;
D = Gs_ABCD.D;

A = expm(A*dt);
B = B*dt;
N = length(B);


Data = [];
X = zeros(N,1);

Y = 0;
R = 1;
t = 0;
E0 = 0;
U0 = 0;

while(t < 9.9)
   if(mod(t,T) < dt)
       E1 = E0;
       E0 = R - Y;
       U1 = U0;   
       U0 = U1 + 0.7*(E0 - 0.6*E1);
       end

    X = A*X + B*U0;
    Y = C*X + D*U0;
   
    t = t + dt;  
    
    Data = [Data ; [Y, U0, R]];
end

t = [0:length(Data)-1]' * dt;
plot(t, Data)
xlabel('Time (seconds)')