discrete-time state feedback LMI SDP formulation

[Moon]

Hi,

I have a discrete-time linear system and I would like to design a state-feedback controller u= - K x such that the following finite horizon cost is minimized: J_0^*(x_0)= \min_{u(0),\dots,u(N-1) } \ x^\top(N)Q_fx(N) + \sum_{k=0}^{N-1} x^\top (k)Q x(k) +u^\top (k)Ru (k).

I have two questions:

1) I tried to implement the infinite horizon using both Yalmip and cvx in a SDP formulation (see reference https://yalmip.github.io/example/lpvstatefeedback/) and compare it with its dual.

This is the code:

clear all

clc

n=2;

sys = drss(n,1,1); 

A=sys.A;

B=sys.B;

Qx = eye(n);

R = 1;

Qf = eye(n);

[Kss,Pss] = dlqr(A,B,Qx,R)

%% Yalmip LMI inf horizon 

P = sdpvar(2,2,'symmetric');

F = [P>=0, [A'*P*A+Qx-P, A'*P*B; B'*P*A, R+B'*P*B ] >=0 , R+B'*P*B>0];

optimize(F,-trace(P));

Pf = value(P)

Kf = inv(R+B'*Pf_m*B)*B'*Pf_m*A

%% SDP formulation. LMI from (A+BK)P(A+BK)'- P with X:= P^-1, L:=KP^-1

X = sdpvar(n,n);

L = sdpvar(1,n,'full' );

MAT = [(-A*X-B*L)'-A*X-B*L, X, L'; X, inv(Qx), zeros(n,1); L, zeros(1,n), inv(R)]

c = [X>=0, MAT>=0] 

optimize(c,trace(X));

Kf_sdp = value(L)*inv(value(X))

%% using cvx 

cvx_begin

    variable X(n,n)

    variable L(1,n)

    minimize (trace(X))

    subject to 

        [(-A*X-B*L)'-A*X-B*L, X, L'; X, inv(Qx), zeros(n,1); L, zeros(1,n), inv(R)]>=0

        X>=0

cvx_end

Kf_cvx = L*inv(X)

I obtained Kf = Kss but Kf_sdp != Kf_cvx != Kf. (sorry if cvx is out of the scope of this group)

There is something wrong in the formulation but I cannot find my mistake. 

2) how is it possible to solve the finite horizon problem using Yalmip? How to compute P{k} at each iteration? are iterations supported in yalmip? do you have any hints?

Thank you very much if you can help me.

[Johan Löfberg]

First, your code doesn't run as Pf_m isn't defined

Second, strict inequalities are not supported

Third, R+B'*P*B >(=)0 is redundant, as you already have that in a diagonal block in the other LMI

Fourth, your LMI for a model in X and L is not correct, looks like some variant of continuous-time model. . The correct model is

MAT = [X [A*X - B*L;L;X]';[A*X - B*L;L;X] blkdiag(X,inv(R),inv(Qx))]

optimize(MAT >= 0, -trace(X))

[Johan Löfberg]

finite horizon is a trivial special case of MPC (i.e., you shouldn't solve as full-fledged MPC, but simply derive the quadratic program in U and x(0), and then compute the explicit solution U = H*x(0) by simply setting the gradient to 0)

https://yalmip.github.io/example/standardmpc/