Solving Riccati equation using LMI algorithm

I am tring to solve Riccati equation using LMI algorithm, but I can't find a solver that can support quadratic semidefinite constraints. The goal is to minimize P then calculate the optimal gain K.

I tried: 'mosek', 'scs', 'sdpt3'

I run the following code and get an error:

Warning: Solver not applicable (scs does not support quadratic semidefinite constraints)

LMI problem infeasible or failed

%% Vehicle Parameters
m=1573;
Iz=2873;
l1=1.58;
l2=1.1;
lh=7;
Ca1=80000;
Ca2=80000;
vx_max=25;
vx_min=15;

%% Matrix Variables
A_22=-2*(Ca2+Ca1)/(m*vx_max);
A_23= 2*(Ca2+Ca1)/m;
A_24=-2*(Ca2*l2-Ca1*l1)/(m*vx_max);
A_42=-2*(Ca2*l2-Ca1*l1)/(Iz*vx_max);
A_43= 2*(Ca2*l2-Ca1*l1)/Iz;
A_44=-2*(-Ca2*((l2)^2)+Ca1*((l1)^2))/(Iz*vx_max);
A_1 = [0 1 0 0; 0 A_22 A_23 A_24; 0 0 0 1; 0 A_42 A_43 A_44];

%% B&Q&R
B = [0; 2*Ca2/m; 0; 2*Ca2*l2/Iz];
q11=1;
q12=0;
q21=0;
q22=0;
Q=diag([q11,q12,q21,q22]);
R = 1;

%% Riccati Equations
% Care_max = A_1'*P + P*A_1 - P*B*inv(R)*B'*P + Q <= 0;
% Care_min = A_2'*P + P*A_2 - P*B*inv(R)*B'*P + Q <= 0;

%% Define LMI variables
P = sdpvar(4);
L = sdpvar(1, 4);

%% Define LMI constraints
% F = [A_1'*P + P*A_1 - P*B*inv(R)*B'*P + Q <= 0;     P >= 0;     [1; 0; 0; 0]*L*B == 1;
%      L >= 0;
%      (A_1 + B*L)'*P + P*(A_1 + B*L) <= 0];

F = [[A_1'*P+P*A_1+Q P*B; B'*P -eye(1,1)] <= 0;     P >= 0;     [1; 0; 0; 0]*L*B == 1;
     L >= 0;
     (A_1 + B*L)'*P + P*(A_1 + B*L) <= 0];

%% Define objective function
obj = trace(Q*P);

%% Define LMI options
opts = sdpsettings('solver', 'scs', 'verbose', 0);

%% Solve LMI
sol = optimize(F, obj, opts);

%% Check feasibility and extract solution
if sol.problem ~= 0
    error('LMI problem infeasible or failed');
else
    P_sol = value(P);
    L_sol = value(L);
end

%% Compute optimal gain
K = inv(R)*B'*P_sol;

%% Display results
disp('Optimal Gain =');
disp(K);

Answer 1

[Samyuktha]

Hi Gilad,

You can choose any of the three algorithms that Quadratic Programming Algorithms ('quadprog') supports:

• 'interior-point-convex' (default)
• 'trust-region-reflective'
• 'active-set'

Since you have a positive semidefinite problem with a large number of linear constraints and not a large number of variables, try 'active-set'. Please refer to the following documentation links for more information:

Choosing the Algorithm, When the Solver Fails

Hope it helps!