Linear Matrix Inequality

jeanluc · August 4, 2021, 2:04pm

Hello everyone!

Here is a code I have tried to solve an LMI

But it failed.

How can I fix it? Will this code give the solution?

using Convex
using SCS
using LinearAlgebra
n = 2
P = Semidefinite(n)
Z12=zeros(Int8, 2, 2)
Z13=zeros(Int8, 2, 2)
Z14=zeros(Int8, 2, 2)
Z15=zeros(Int8, 2, 2)
Z16=zeros(Int8, 2, 2)
Z21=zeros(Int8, 2, 2)
Z23=zeros(Int8, 2, 2)
Z24=zeros(Int8, 2, 2)
Z25=zeros(Int8, 2, 2)
Z26=zeros(Int8, 2, 2)
Z31=zeros(Int8, 2, 2)
Z32=zeros(Int8, 2, 2)
Z34=zeros(Int8, 2, 2)
Z35=zeros(Int8, 2, 2)
Z36=zeros(Int8, 2, 2)
Z41=zeros(Int8, 2, 2)
Z42=zeros(Int8, 2, 2)
Z43=zeros(Int8, 2, 2)
Z45=zeros(Int8, 2, 2)
Z46=zeros(Int8, 2, 2)
Z51=zeros(Int8, 2, 2)
Z52=zeros(Int8, 2, 2)
Z53=zeros(Int8, 2, 2)
Z54=zeros(Int8, 2, 2)
Z56=zeros(Int8, 2, 2)
Z61=zeros(Int8, 2, 2)
Z62=zeros(Int8, 2, 2)
Z63=zeros(Int8, 2, 2)
Z64=zeros(Int8, 2, 2)
Z65=zeros(Int8, 2, 2)
P1=[1.1 0; 0 1.1]
P2=[2.1 0; 0 2.1]
P3=[3.1 0; 0 3.1]
P4=[4.1 0; 0 4.1]
P5=[5.1 0; 0 5.1]
P6=[6.1 0; 0 6.1]
Id=1*Matrix{Float64}(I, 12, 12)
C=1*Matrix{Float64}(I, 2, 2)
R=1*Matrix{Float64}(I, 2, 2)
S=0.3*Matrix{Float64}(I, 2, 2)
Q2=3*Matrix{Float64}(I, 2, 2)
Q1=2*Matrix{Float64}(I, 2, 2)
B=0.0001*Matrix{Float64}(I, 2, 2)
Q =[-2 1; 1 -3]
A = [0 1; 0 -0.5]
P1VARON=inv(P1)
P2VARON=inv(P2)
P3VARON=inv(P3)
P4VARON=inv(P4)
P5VARON=inv(P5)
P6VARON=inv(P6)
rho = Variable()
elementmatone=A'*P+P'A-C'*Q*C-rho*P+(P1+P2+P3+P4+P5)*1/2
R1=[elementmatone Z12 Z13 Z14 Z15 Z16]
R2=[Z21 2*P*P1VARON*P' Z13 Z14 Z15 Z16]
R3=[Z31   Z32  B'*P*P2VARON*P'*B  Z34 Z35 Z36]
R4=[Z41   Z42  Z43  2*P*P3VARON*P'  Z45 Z46]
R5=[Z51   Z52  Z53  Z54   B'*P4VARON*B Z56]
R6=[Z61   Z62  Z63  Z64    Z65 (2*P-C'*S)*P5VARON]
LMI1=[-R1;-R2;-R3;-R4;-R5;-R6]
p.objective = rho
p.constraints = [(LMI1>Id),(P>0)]
p = minimize(rho, constraints)
solve!(p)

Here is the error I got

julia> using Convex

julia> using SCS

julia> using LinearAlgebra

julia> n = 2
2

julia> P = Semidefinite(n)
Variable
size: (2, 2)
sign: real
vexity: affine
id: 106…234

julia> Z12=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z13=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z14=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z15=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z16=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z21=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z23=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z24=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z25=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z26=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z31=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z32=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z34=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z35=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z36=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z41=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z42=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z43=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z45=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z46=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z51=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z52=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z53=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z54=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z56=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z61=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z62=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z63=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z64=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> Z65=zeros(Int8, 2, 2)
2×2 Matrix{Int8}:
 0  0
 0  0

julia> P1=[1.1 0; 0 1.1]
2×2 Matrix{Float64}:
 1.1  0.0
 0.0  1.1

julia> P2=[2.1 0; 0 2.1]
2×2 Matrix{Float64}:
 2.1  0.0
 0.0  2.1

julia> P3=[3.1 0; 0 3.1]
2×2 Matrix{Float64}:
 3.1  0.0
 0.0  3.1

julia> P4=[4.1 0; 0 4.1]
2×2 Matrix{Float64}:
 4.1  0.0
 0.0  4.1

julia> P5=[5.1 0; 0 5.1]
2×2 Matrix{Float64}:
 5.1  0.0
 0.0  5.1

julia> P6=[6.1 0; 0 6.1]
2×2 Matrix{Float64}:
 6.1  0.0
 0.0  6.1

julia> Id=1*Matrix{Float64}(I, 12, 12)
12×12 Matrix{Float64}:
 1.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0

julia> C=1*Matrix{Float64}(I, 2, 2)
2×2 Matrix{Float64}:
 1.0  0.0
 0.0  1.0

julia> R=1*Matrix{Float64}(I, 2, 2)
2×2 Matrix{Float64}:
 1.0  0.0
 0.0  1.0

julia> S=0.3*Matrix{Float64}(I, 2, 2)
2×2 Matrix{Float64}:
 0.3  0.0
 0.0  0.3

julia> Q2=3*Matrix{Float64}(I, 2, 2)
2×2 Matrix{Float64}:
 3.0  0.0
 0.0  3.0

julia> Q1=2*Matrix{Float64}(I, 2, 2)
2×2 Matrix{Float64}:
 2.0  0.0
 0.0  2.0

julia> B=0.0001*Matrix{Float64}(I, 2, 2)
2×2 Matrix{Float64}:
 0.0001  0.0
 0.0     0.0001

julia> Q =[-2 1; 1 -3]
2×2 Matrix{Int64}:
 -2   1
  1  -3

julia> A = [0 1; 0 -0.5]
2×2 Matrix{Float64}:
 0.0   1.0
 0.0  -0.5

julia> P1VARON=inv(P1)
2×2 Matrix{Float64}:
 0.909091  0.0
 0.0       0.909091

julia> P2VARON=inv(P2)
2×2 Matrix{Float64}:
 0.47619  0.0
 0.0      0.47619

julia> P3VARON=inv(P3)
2×2 Matrix{Float64}:
 0.322581  0.0
 0.0       0.322581

julia> P4VARON=inv(P4)
2×2 Matrix{Float64}:
 0.243902  0.0
 0.0       0.243902

julia> P5VARON=inv(P5)
2×2 Matrix{Float64}:
 0.196078  0.0
 0.0       0.196078

julia> P6VARON=inv(P6)
2×2 Matrix{Float64}:
 0.163934  0.0
 0.0       0.163934

julia> rho = Variable()
Variable
size: (1, 1)
sign: real
vexity: affine
id: 243…140

julia> elementmatone=A'*P+P'A-C'*Q*C-rho*P+(P1+P2+P3+P4+P5)*1/2
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
+ (Convex.NotDcp; real)
├─ * (affine; real)
│  ├─ 2×2 adjoint(::Matrix{Float64}) with eltype Float64
│  └─ 2×2 real variable (id: 106…234)
├─ * (affine; real)
│  ├─ adjoint (affine; real)
│  │  └─ 2×2 real variable (id: 106…234)
│  └─ 2×2 Matrix{Float64}
├─ 2×2 Matrix{Float64}
├─ - (Convex.NotDcp; real)
│  └─ * (Convex.NotDcp; real)
│     ├─ real variable (id: 243…140)
│     └─ 2×2 real variable (id: 106…234)
└─ 2×2 Matrix{Float64}

julia> R1=[elementmatone Z12 Z13 Z14 Z15 Z16]
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
hcat (Convex.NotDcp; real)
├─ + (Convex.NotDcp; real)
│  ├─ * (affine; real)
│  │  ├─ 2×2 adjoint(::Matrix{Float64}) with eltype Float64
│  │  └─ 2×2 real variable (id: 106…234)
│  ├─ * (affine; real)
│  │  ├─ adjoint (affine; real)
│  │  │  └─ …
│  │  └─ 2×2 Matrix{Float64}
│  ├─ 2×2 Matrix{Float64}
│  ├─ - (Convex.NotDcp; real)
│  │  └─ * (Convex.NotDcp; real)
│  │     ├─ …
│  │     └─ …
│  └─ 2×2 Matrix{Float64}
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
└─ 2×2 Matrix{Int8}

julia> R2=[Z21 2*P*P1VARON*P' Z13 Z14 Z15 Z16]
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
hcat (Convex.NotDcp; real)
├─ 2×2 Matrix{Int8}
├─ * (Convex.NotDcp; real)
│  ├─ * (affine; real)
│  │  ├─ * (affine; real)
│  │  │  ├─ …
│  │  │  └─ …
│  │  └─ 2×2 Matrix{Float64}
│  └─ adjoint (affine; real)
│     └─ 2×2 real variable (id: 106…234)
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
└─ 2×2 Matrix{Int8}

julia> R3=[Z31   Z32  B'*P*P2VARON*P'*B  Z34 Z35 Z36]
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
hcat (Convex.NotDcp; real)
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
├─ * (Convex.NotDcp; real)
│  ├─ * (Convex.NotDcp; real)
│  │  ├─ * (affine; real)
│  │  │  ├─ …
│  │  │  └─ …
│  │  └─ adjoint (affine; real)
│  │     └─ …
│  └─ 2×2 Matrix{Float64}
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
└─ 2×2 Matrix{Int8}

julia> R4=[Z41   Z42  Z43  2*P*P3VARON*P'  Z45 Z46]
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
hcat (Convex.NotDcp; real)
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
├─ * (Convex.NotDcp; real)
│  ├─ * (affine; real)
│  │  ├─ * (affine; real)
│  │  │  ├─ …
│  │  │  └─ …
│  │  └─ 2×2 Matrix{Float64}
│  └─ adjoint (affine; real)
│     └─ 2×2 real variable (id: 106…234)
├─ 2×2 Matrix{Int8}
└─ 2×2 Matrix{Int8}

julia> R5=[Z51   Z52  Z53  Z54   B'*P4VARON*B Z56]
2×12 Matrix{Float64}:
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  2.43902e-9  0.0         0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0         2.43902e-9  0.0  0.0

julia> R6=[Z61   Z62  Z63  Z64    Z65 (2*P-C'*S)*P5VARON]
hcat (affine; real)
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
├─ 2×2 Matrix{Int8}
└─ * (affine; real)
   ├─ + (affine; real)
   │  ├─ * (affine; real)
   │  │  ├─ …
   │  │  └─ …
   │  └─ 2×2 Matrix{Float64}
   └─ 2×2 Matrix{Float64}

julia> LMI1=[-R1;-R2;-R3;-R4;-R5;-R6]
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
┌ Warning: Expression not DCP compliant. Trying to solve non-DCP compliant problems can lead to unexpected behavior.
└ @ Convex C:\Users\user\.julia\packages\Convex\UNspE\src\dcp.jl:25
transpose (Convex.NotDcp; real)
└─ hcat (Convex.NotDcp; real)
   ├─ transpose (Convex.NotDcp; real)
   │  └─ - (Convex.NotDcp; real)
   │     └─ …
   ├─ transpose (Convex.NotDcp; real)
   │  └─ - (Convex.NotDcp; real)
   │     └─ …
   ├─ transpose (Convex.NotDcp; real)
   │  └─ - (Convex.NotDcp; real)
   │     └─ …
   ├─ transpose (Convex.NotDcp; real)
   │  └─ - (Convex.NotDcp; real)
   │     └─ …
   ├─ transpose (constant; negative)
   │  └─ 2×12 Matrix{Float64}
   └─ transpose (affine; real)
      └─ - (affine; real)
         └─ …

julia> p.objective = rho
ERROR: UndefVarError: p not defined
Stacktrace:
 [1] top-level scope
   @ REPL[66]:1

julia> p.constraints = [(LMI1>Id),(P>0)]
ERROR: UndefVarError: p not defined
Stacktrace:
 [1] top-level scope
   @ REPL[67]:1

julia> p = minimize(rho, constraints)
ERROR: MethodError: no method matching minimize(::Variable, ::typeof(constraints))
Closest candidates are:
  minimize(::Convex.AbstractExpr) at C:\Users\user\.julia\packages\Convex\UNspE\src\problems.jl:74
  minimize(::Convex.AbstractExpr, ::Constraint...; numeric_type) at C:\Users\user\.julia\packages\Convex\UNspE\src\problems.jl:72
  minimize(::Convex.AbstractExpr, ::Array{var"#s39", N} where {var"#s39"<:Constraint, N}; numeric_type) at C:\Users\user\.julia\packages\Convex\UNspE\src\problems.jl:74
Stacktrace:
 [1] top-level scope
   @ REPL[68]:1

julia> solve!(p)
ERROR: UndefVarError: p not defined
Stacktrace:
 [1] top-level scope
   @ REPL[69]:1

Any assistance you can provide would be greatly appreciated.

zdenek_hurak · August 4, 2021, 2:58pm

This is not an LMI, is it? The terms R2 throug R5 contain (weighted) squares of the variable P. Do you have a mathematical description of your problem? It might be useful to accompany your code with unambiguous description of the problem you want to solve.

jeanluc · August 4, 2021, 3:53pm

Thanks a lot. How can I post the LMI here?

zdenek_hurak · August 4, 2021, 4:11pm

While writing the response, you can insert LaTeX code pretty much in the same way you include it while writing a report or paper. Inline \cos\alpha^2+\sin\alpha^2=1 by inserting the LaTeX math between two $ symbols or displayed

\cos\alpha^2+\sin\alpha^2=1

by inserting the LaTeX math between two pairs $$ of symbols.

jeanluc · August 4, 2021, 5:31pm

Thnx a lot

\begin{bmatrix} \Xi & 0 & 0 & 0& 0 & 0 \\ 0 &2P\Pi_1^{-1}P^T & 0 & 0& 0& 0 \\ 0 & 0 & B^TP\Pi_2^{-1}P^TB & 0 & 0 & 0\\ 0 & 0 & 0 & 2P\Pi_3^{-1}P^T& 0& 0 \\ 0 & 0 & 0 & 0 & B^T\Pi_4^{-1}B& 0 \\ 0 & 0 & 0 & 0 & 0&(2P-C^TS)\Pi_5^{-1} \\ \end{bmatrix} < 0

in which

\begin{split} \Xi=&A^TP+P^TA+(\Pi_1+\Pi_2+\Pi_3+\Pi_4+\Pi_5)\frac{\epsilon}{2}-\rho P-C^TQC\\ \end{split}

I have used the R_i in the code as rows to make a block Matrix because it is a block LMI.

zdenek_hurak · August 4, 2021, 6:04pm

A few comments:

You do not have to build the huge matrix (function) by yourself - leave it up to the Convex.jl parser. Just create 6 constraints and then stack them into a vector of constraints.
My remark from the original response still holds, you have not reacted to it. If P is a (matrix) variable, the second, third and fourth entries on the diagonal contain a product of P, hence your constraints are not linear.
The fifth entry on the diagonal is irrelevant, it is just a constant matrix.
But there is another problem with it: it is a product of P and P' scaled by a positive constant. No way to make it negative definite. Are you sure your problem statement is correct? In fact, this comment applies to some other components on the diagonal.
Your matrices P1 through P5 in your code (or ∏₁ through ∏₅ in the above typeset math) are just scaled identity matrices. You can reduce the complexity of your formulas quite a bit by exploiting this knowledge about the structure of these matrices.

added later

Semidefinite matrix variables (the matrix P in your case) are assumed symmetric in Convex.jl (and elsewhere, I bet). No need to perform tranposition then.

jeanluc · August 4, 2021, 6:19pm

Dear Zdeněk
Yes, you are right I need to recheck the problem.
The nonlinear terms and the constant term
Matrix P is also symmetric as you mentioned I need to reconsider it again.
For about the construction of block LMI via Convex.jl are there some examples?

Merci

zdenek_hurak · August 4, 2021, 7:32pm

In fact, stacking together several LMIs into a single LMI does not differ from how other constraints are handled in Convex – you can just create several constraints and attach them to your problem. The parser will do what is needed. But if an example is needed, here comes the LMI Hello world example

julia> P = Semidefinite(2)
Variable
size: (2, 2)
sign: real
vexity: affine
id: 174…604

julia> using LinearAlgebra

julia> I₂ = Matrix{Float64}(I,2,2)
2×2 Matrix{Float64}:
 1.0  0.0
 0.0  1.0

julia> c1 = A'*P+P*A <= I₂
<= constraint (affine)
├─ + (affine; real)
│  ├─ * (affine; real)
│  │  ├─ 2×2 adjoint(::Matrix{Float64}) with eltype Float64
│  │  └─ 2×2 real variable (id: 174…604)
│  └─ * (affine; real)
│     ├─ 2×2 real variable (id: 174…604)
│     └─ 2×2 Matrix{Float64}
└─ 2×2 Matrix{Float64}

julia> c2 = P >= I₂
>= constraint (affine)
├─ 2×2 real variable (id: 174…604)
└─ 2×2 Matrix{Float64}

julia> constraints = [c1,c2]
2-element Vector{Constraint}:
 <= constraint (affine)
├─ + (affine; real)
│  ├─ * (affine; real)
│  │  ├─ 2×2 adjoint(::Matrix{Float64}) with eltype Float64
│  │  └─ 2×2 real variable (id: 174…604)
│  └─ * (affine; real)
│     ├─ 2×2 real variable (id: 174…604)
│     └─ 2×2 Matrix{Float64}
└─ 2×2 Matrix{Float64}
 >= constraint (affine)
├─ 2×2 real variable (id: 174…604)
└─ 2×2 Matrix{Float64}

julia> problem = satisfy(constraints)
minimize
└─ 0
subject to
├─ <= constraint (affine)
│  ├─ + (affine; real)
│  │  ├─ * (affine; real)
│  │  │  ├─ …
│  │  │  └─ …
│  │  └─ * (affine; real)
│  │     ├─ …
│  │     └─ …
│  └─ 2×2 Matrix{Float64}
└─ >= constraint (affine)
   ├─ 2×2 real variable (id: 174…604)
   └─ 2×2 Matrix{Float64}

status: `solve!` not called yet

julia> solve!(problem,SCS.Optimizer)
----------------------------------------------------------------------------
    SCS v2.1.3 - Splitting Conic Solver
    (c) Brendan O'Donoghue, Stanford University, 2012
----------------------------------------------------------------------------
Lin-sys: sparse-indirect, nnz in A = 22, CG tol ~ 1/iter^(2.00)
eps = 1.00e-05, alpha = 1.50, max_iters = 5000, normalize = 1, scale = 1.00
acceleration_lookback = 10, rho_x = 1.00e-03
Variables n = 5, constraints m = 13
Cones:  primal zero / dual free vars: 2
    linear vars: 8
    sd vars: 3, sd blks: 1
Setup time: 1.01e-04s
----------------------------------------------------------------------------
 Iter | pri res | dua res | rel gap | pri obj | dua obj | kap/tau | time (s)
----------------------------------------------------------------------------
     0| 1.11e+19  1.16e+19  1.00e+00 -1.96e+19  2.82e+19  2.52e+19  1.05e-04 
    16| 5.19e+08  7.80e+08  1.00e+00  3.60e+08  2.61e+19  2.61e+19  9.32e-03 
----------------------------------------------------------------------------
Status: Infeasible
Timing: Solve time: 9.34e-03s
    Lin-sys: avg # CG iterations: 4.00, avg solve time: 2.33e-06s
    Cones: avg projection time: 3.84e-07s
    Acceleration: avg step time: 5.34e-04s
----------------------------------------------------------------------------
Certificate of primal infeasibility:
dist(y, K*) = 0.0000e+00
|A'y|_2 * |b|_2 = 1.1973e-10
b'y = -1.0000
============================================================================
┌ Warning: Problem status INFEASIBLE; solution may be inaccurate.
└ @ Convex ~/.julia/packages/Convex/UNspE/src/solution.jl:250

Obviously, the Lyapunov inequality does not admit a positive definite solution, hence the matrix A is not stable. This can be verified by

julia> eigvals(A)
2-element Vector{Float64}:
 -0.7102362843766895
  1.1093314611610228

jeanluc · August 4, 2021, 8:19pm

I have checked the problem the nonlinear terms and the constant term were exactly the same as what I have written here.

Does this mean that it cannot be solved?

I really appreciate your kind helps.

zdenek_hurak · August 4, 2021, 9:19pm

Look at, say, the second entry on the diagonal. Taking into consideration that the P1 matrix (as used in the code, or Π1 in the formula) is just an identity matrix multiplied by a positive constant, we can simplify the corresponding inequality to P*P<0. Now, do you assume that P is a real matrix or complex? If real, there is clearly no feasible solution.

Similarly – but without any degree of freedom – have a look a the fifth entry on the diagonal. There are no variables there (the only variables in your problem are P and ρ ). Since both B and Π4 is again just a diagonal matrices scaled by a positive number, the fifth entry on the diagonal is a positive definite matrix. No way to make it comply with the <0 inequality.

zdenek_hurak · August 5, 2021, 10:33am

I don’t immediately see how this vector inequality transforms into (any part of) your original matrix inequality. Anyway, this inequality is not linear either.

jeanluc · August 5, 2021, 11:32am

Dear Zdeněk

Thank you so much for your great notes. It seems it needs more investigation.