How to plot the feasibility domain of LMI

[Wail Hamdi]

hello sir Johan

please can you help me how to plot the feasibility domain of different LMIs like this one

I want to compare two theorems and see which one has the largest feasibility domain. are those pictures plotted using Yalmip or something else?

[Johan Löfberg]

An LMI defines a convex set, and that is definitely not what you show ( it is discrete)

To plot a set, you use the command plot

>> help lmi/plot
 plot  Plots the feasible region of a set of constraints
 
  p = plot(C,x,c,n,options)
 
  Note that only convex sets are allowed, or union of convex sets
  represented using binary variables (either defined explictly or
  introduced by YALMIP when modelling, e.g., mixed integer linear
  programming representable operators)
 
  C:  Constraint object
  x:  Plot variables [At most three variables]
  c:  color [double] ([r g b] format) or char from 'rymcgbk'
  n:  #vertices [double ] (default 100 in 2D and 300 otherwise)
  options: options structure from sdpsettings

Various examples all over the Wiki, as in

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

[Wail Hamdi]

I'm sorry I didn't explain it properly
for example, I have the following matrices:
A1=[-0.5 5; a -1];
A2=[-0.5 5; a+b -1];
they are a function of 'a' and 'b'
I want to plot according to 'a' and 'b' the points where an LMI is  feasible
for example the following LMI
Ai'*P+P*Ai<0

then I want to plot the points where another LMI is  feasible, for example:

Ai'*Pi+Pi*Ai<0

I have seen similar maps in many articles where the authors demonstrated that an LMI is better than another because it is applicable for more couples (a,b)

[Johan Löfberg]

then you have to do what they do, basically grid over a and b and solve and check if feasible.

quickly done with an optimizer

sdpvar a b
P = sdpvar(2);

A1=[-0.5 5; a -1];
A2=[-0.5 5; a+b -1];

Model = [A1'*P+P*A1 <= 0,
         A2'*P+P*A2 <= 0,
         P >= eye(2)];

S = optimizer(Model,[],sdpsettings('solver','mosek'),[a;b],P);
[A,B] = meshgrid(-3:0.25:3,-1:.1:1);
[~,infeas] = S([A(:)';B(:)']);
plot(A(find(~infeas)),B(find(~infeas)),'bo');
hold on
plot(A(find(infeas)),B(find(infeas)),'ro');

[Wail Hamdi]

yes I have done this before but I found that it depends on the Solver, for example, Lmilab succeed to solve the problem where MOsek did not
so there may be points feasible but the solver could not solve them, so this will give me imprecision in mapping.

[Johan Löfberg]

There is no other way.

If you used lmilab and yalmip, then you cannot trust any result as lmilab does not tell yalmip if the problem was solved or not.