G=zpk([],[1 1 1],1);
Ts=0.1;
Gd=c2d(G,Ts);
[A,B,C,D]=ssdata(Gd);
T=100;							% horizont optimalizace
Q=100*eye(T);
R=eye(T);
x=zeros(size(A,1),1);
r=ones(T,1);
for i=1:T
   AA(2*i-1,i)=1;
   AA(2*i,i)=-1;
   b(2*i-1)=100;
   b(2*i)=100;
end

for i=1:T
   S(i,i)=D;
   V(i,:)=C*A^(i-1);
   for j=1:i-1
      S(i,j)=C*A^(i-j-1)*B;
   end
end

W=S'*Q*S+R;

for i=1:100
	w=S'*Q*(V*x-r);
   u=quadprog(W,w,AA,b);
   y(i)=C*x+D*u(1);
   x=A*x+B*u(1);
   uu(i)=u(1);
end
figure(1);plot(y);
figure(2);plot(uu);