# Sum-Of-Squares SOC constraint in JuMP

**URL:** <https://discourse.julialang.org/t/sum-of-squares-soc-constraint-in-jump/50091>\
**Category:** Optimization (Mathematical)\
**Tags:** jump, optimization\
**Created:** [November 13, 2020, 12:54pm UTC](https://discourse.julialang.org/t/sum-of-squares-soc-constraint-in-jump/50091 "2020-11-13T12:54:52Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![andrewrosemberg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/andrewrosemberg/32/9902_2.png) [@andrewrosemberg](https://discourse.julialang.org/u/andrewrosemberg)\
**Post date:** [November 13, 2020, 12:54pm UTC](https://discourse.julialang.org/t/sum-of-squares-soc-constraint-in-jump/50091/1 "2020-11-13T12:54:52Z")

</div>

I’m trying to add a constraint of type `(x'Ax <= b)` in JuMP to be interpreted as SOC by the solver.  
A is Positive definite.  
Seems to be the same error as presented in [this issue](https://github.com/jump-dev/MosekTools.jl/issues/2#issuecomment-471686494), but `A+A'` has positive eigenvalues:

```julia
eigvals(cov(d).data + cov(d).data')
10-element Array{Float64,1}:
  2.0014222400461272
  2.0605866115342453
  2.2799706986926673
  2.5317489320096396
  2.7815912421610607
  3.459123499910021
  4.181503965703311
  4.970808766218942
  6.5318231430488245
 53.903463261426204

```

`A` also has positive eignenvalues:

```julia
eigvals(cov(d).data)
10-element Array{Float64,1}:
  1.0007111200230636
  1.0302933057671226
  1.1399853493463337
  1.2658744660048198
  1.3907956210805303
  1.7295617499550104
  2.0907519828516556
  2.485404383109471
  3.2659115715244122
 26.951731630713102

```

From stacktrace, I can see that this happens at the call of `cholesky` ([ref](https://github.com/jump-dev/MathOptInterface.jl/blob/7b24c6d3e8fd27c5004ce2aef2986d4a2de0f3ad/src/Bridges/Constraint/quad_to_soc.jl#L75)) and seems that the “error” was already reported: [https://github.com/JuliaLang/julia/issues/30831](https://github.com/JuliaLang/julia/issues/30831)

Could the solution be to call `cholesky(Hermitian(matrix))` or the same with `Symmetric ` in MOI?

---

<div class="post-metadata">

**Author:** ![andrewrosemberg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/andrewrosemberg/32/9902_2.png) [@andrewrosemberg](https://discourse.julialang.org/u/andrewrosemberg)\
**Post date:** [November 13, 2020, 4:30pm UTC](https://discourse.julialang.org/t/sum-of-squares-soc-constraint-in-jump/50091/2 "2020-11-13T16:30:58Z")

</div>

> [@andrewrosemberg](#):
>
> x’Ax \<= b

Minimal working case:

```julia
using Random
using JuMP
using ECOS

n = 2
rng = MersenneTwister(1)
means = randn(rng, n)
sqrt_cov = rand(rng, n, n) ./ 2.0
Σ = sqrt_cov' * sqrt_cov

model = Model(ECOS.Optimizer)
@variable(model, x[i=1:n] >= 0);
@variable(model, y[i=1:n] <= 0);
v = x + y
@objective(model, Min, 0);
@constraint(model, v'Σ*v / 321.0^2 <= 2);

optimize!(model)

```

---

<div class="post-metadata">

**Author:** ![andrewrosemberg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/andrewrosemberg/32/9902_2.png) [@andrewrosemberg](https://discourse.julialang.org/u/andrewrosemberg)\
**Post date:** [November 13, 2020, 4:33pm UTC](https://discourse.julialang.org/t/sum-of-squares-soc-constraint-in-jump/50091/3 "2020-11-13T16:33:26Z")

</div>

Curious enough, the following doesn’t error:

```julia
using Random
using JuMP
using ECOS

n = 2
rng = MersenneTwister(1)
means = randn(rng, n)
sqrt_cov = rand(rng, n, n) ./ 2.0
Σ = sqrt_cov' * sqrt_cov

model = Model(ECOS.Optimizer)
@variable(model, x[i=1:n] >= 0);
@objective(model, Min, 0);
@constraint(model, x'Σ*x / 321.0^2 <= 2);

optimize!(model)

```

---

<div class="post-metadata">

**Author:** ![andrewrosemberg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/andrewrosemberg/32/9902_2.png) [@andrewrosemberg](https://discourse.julialang.org/u/andrewrosemberg)\
**Post date:** [November 13, 2020, 4:49pm UTC](https://discourse.julialang.org/t/sum-of-squares-soc-constraint-in-jump/50091/4 "2020-11-13T16:49:30Z")

</div>

and even more curious, the following doesn’t error:

```julia
using Random
using JuMP
using ECOS

n = 2
rng = MersenneTwister(1)
means = randn(rng, n)
sqrt_cov = rand(rng, n, n) ./ 2.0
Σ = sqrt_cov' * sqrt_cov

model = Model(ECOS.Optimizer)
@variable(model, x[i=1:n] >= 0);
@variable(model, y[i=1:n] <= 0);
@variable(model, v[i=1:n])
@constraint(model, v .== x + y)
@objective(model, Min, 0);
@constraint(model, v'Σ*v / 321.0^2 <= 2);

optimize!(model)

```

---

<div class="post-metadata">

**Author:** ![odow](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/odow/32/28685_2.png) [@odow](https://discourse.julialang.org/u/odow)\
**Post date:** [November 15, 2020, 9:35pm UTC](https://discourse.julialang.org/t/sum-of-squares-soc-constraint-in-jump/50091/5 "2020-11-15T21:35:31Z")

</div>

Just use the `SecondOrderCone` syntax:

```nohighlight
@constraint(model, [t; x] in SecondOrderCone())

```

---

<div class="post-metadata">

**Author:** ![andrewrosemberg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/andrewrosemberg/32/9902_2.png) [@andrewrosemberg](https://discourse.julialang.org/u/andrewrosemberg)\
**Post date:** [November 16, 2020, 8:52am UTC](https://discourse.julialang.org/t/sum-of-squares-soc-constraint-in-jump/50091/6 "2020-11-16T08:52:04Z")

</div>

yeah I did that, I was just curious why this didn’t work. However, Joaquim just showed me that indeed when `v` is not a variable, the constraint uses a matrix that is Positive-Semidefinite instead of Positive-Definite.

```julia
v'Σ*v = (x + y)'Σ*(x + y) = x'Σ*(x + y) + y'Σ*(x + y) = x'Σ*x + x'Σ*y + y'Σ*x + y'Σ*y
Since x'Σ*y == y'Σ*x, then:
v'Σ*v = (x + y)'Σ*(x + y) = x'Σ*(x + y) + y'Σ*(x + y) = x'Σ*x + x'Σ*y + y'Σ*x + y'Σ*y = x'Σ*x + 2x'Σ*y + y'Σ*y
And:
x'Σ*x + 2x'Σ*y + y'Σ*y == [x; y]' [Σ Σ;Σ Σ] [x; y] 
At last:
[Σ Σ;Σ Σ] >= 0 and not > 0

```