# Hypercube

Given a Sum-of-Squares constraint on an algebraic set:

\[g_1(x) = , \ldots, g_m(x) = 0 \Rightarrow p(x) \ge 0.\]

We can either use the certificate:

\[p(x) = s(x) + \lambda_1(x) g_1(x) + \cdots + \lambda_m(x) g_m(x), s_0(x) \text{ is SOS},\]

or

\[p(x) \equiv s(x) \pmod{\la g_1(x), \ldots, g_m(x) \ra}, s_0(x) \text{ is SOS}.\]

the second one leads to a *simpler* SDP but needs to compute a *Gr\"obner* basis:

- SemialgebraicSets implements Buchberger's algorithm.
- The
`@set`

macro recognizes variable fix, e.g.,`x = 1`

and provides shortcut. - If you know a \alert{better} way to take modulo, better create your \alert{own} type of algebraic set!

We illustrate this in this example.

```
using DynamicPolynomials
@polyvar x[1:3]
p = sum(x)^2
using SumOfSquares
S = algebraicset([xi^2 - 1 for xi in x])
```

```
Algebraic Set defined by 3 equalities
-1//1 + x[1]^2 = 0
-1//1 + x[2]^2 = 0
-1//1 + x[3]^2 = 0
```

We will now search for the minimum of `x`

over `S`

using Sum of Squares Programming. We first need to pick an SDP solver, see here for a list of the available choices.

```
import CSDP
solver = optimizer_with_attributes(CSDP.Optimizer, MOI.Silent() => true)
function min_algebraic(S)
model = SOSModel(solver)
@variable(model, α)
@objective(model, Max, α)
@constraint(model, c, p >= α, domain = S)
optimize!(model)
@show termination_status(model)
@show objective_value(model)
end
min_algebraic(S)
```

`-4.404406006575101e-10`

Note that the minimum is in fact `1`

. Indeed, since each variables is odd (it is either `-1`

or `1`

) and there is an odd number of variables, their sum is odd. Therefore it cannot be zero!

We can see that the Gröbner basis of `S`

was computed

```
@show S.I.gröbner_basis
S.I.algo
```

`Buchberger(1.4901161193847656e-8, SemialgebraicSets.presort!, SemialgebraicSets.normal_selection)`

The Gröbner basis is simple to compute in this case as the vector of `xi^2 - 1`

is already a Gröbner basis. However, we still need to divide polynomials by the Gröbner basis which can be simplified in this case.

```
const MP = MultivariatePolynomials
const SS = SemialgebraicSets
struct HypercubeIdeal{V} <: SS.AbstractPolynomialIdeal
variables::Vector{V}
end
struct HypercubeSet{V} <: SS.AbstractAlgebraicSet
ideal::HypercubeIdeal{V}
end
MP.variables(set::HypercubeSet) = MP.variables(set.ideal)
MP.variables(ideal::HypercubeIdeal) = ideal.variables
Base.similar(set::HypercubeSet, ::Type) = set
SS.ideal(set::HypercubeSet) = set.ideal
function Base.rem(p, set::HypercubeIdeal)
return MP.polynomial(map(MP.terms(p)) do term
mono = MP.monomial(term)
new_mono = one(mono)
for (var, exp) in powers(mono)
if var in set.variables
exp = rem(exp, 2)
end
new_mono *= var^exp
end
MP.coefficient(term) * new_mono
end)
end
H = HypercubeSet(HypercubeIdeal(x))
min_algebraic(H)
```

`-4.404406006575101e-10`

Let's now try to find the correct lower bound:

```
function min_algebraic_rational(S, d)
model = SOSModel(solver)
@variable(model, q, SOSPoly(MP.monomials(x, 0:d)))
deno = q + 1
@constraint(model, c, deno * p >= deno, domain = S)
optimize!(model)
@show termination_status(model)
end
```

`min_algebraic_rational (generic function with 1 method)`

With `d = 0`

, it's the same as previously

`min_algebraic_rational(H, 0)`

`INFEASIBLE::TerminationStatusCode = 2`

But with `d = 1`

, we can find the correct lower bound

`min_algebraic_rational(H, 1)`

`OPTIMAL::TerminationStatusCode = 1`

*This page was generated using Literate.jl.*