% ECE504 optimal control example
% DRB 08-Dec-2009

A = [0 1; 0 0];
B = [0;1];
C = [1 0];
D = 0;
rho = 1;
Q = C'*C;
x0 = [1;1];
Tfinal = 15;

% solve algebraic Riccati equation
[P,L,G] = care(A,B,Q,rho); 

% make feedback rule
F = (1/rho)*B'*P';

%simulate
sys = ss(A-B*F,B,C,D);
[y,t,x] = initial(sys,x0,Tfinal);
u = -F*x';
subplot(2,1,1);
plot(t,u,t,y);
xlabel('time');
legend('input','output');
grid on
subplot(2,1,2);
plot(t,cumsum(u.^2),t,cumsum(y.^2));
xlabel('time');
legend('cumulative input energy','cumulative output energy','Location','SouthEast');
grid on
