Term sparsity

Adapted from: Example 3.5 of [WML20b]

[WML20a] Wang, Jie, Victor Magron, and Jean-Bernard Lasserre. TSSOS: A Moment-SOS hierarchy that exploits term sparsity. arXiv preprint arXiv:1912.08899 (2020).

[WML20b] Wang, Jie, Victor Magron, and Jean-Bernard Lasserre. Chordal-TSSOS: a moment-SOS hierarchy that exploits term sparsity with chordal extension. arXiv preprint arXiv:2003.03210 (2020).

using DynamicPolynomials
@polyvar x[1:3]

(DynamicPolynomials.PolyVar{true}[x₁, x₂, x₃],)

We would like to find the minimum value of the polynomial

poly = x[1]^2 - 2x[1]*x[2] + 3x[2]^2 - 2x[1]^2*x[2] + 2x[1]^2*x[2]^2 - 2x[2]*x[3] + 6x[3]^2 + 18x[2]^2*x[3] - 54x[2]*x[3]^2 + 142x[2]^2*x[3]^2

\[ 2x_{1}^{2}x_{2}^{2} + 142x_{2}^{2}x_{3}^{2} - 2x_{1}^{2}x_{2} + 18x_{2}^{2}x_{3} - 54x_{2}x_{3}^{2} + x_{1}^{2} - 2x_{1}x_{2} + 3x_{2}^{2} - 2x_{2}x_{3} + 6x_{3}^{2} \]

The minimum value of the polynomial can be found to be zero.

import CSDP
solver = CSDP.Optimizer
using SumOfSquares
function sos_min(sparsity)
    model = Model(solver)
    @variable(model, t)
    @objective(model, Max, t)
    con_ref = @constraint(model, poly - t in SOSCone(), sparsity = sparsity)
    optimize!(model)
    return value(t), moment_matrix(con_ref)
end

bound, ν = sos_min(NoSparsity())
bound

3 . 0 4 5 1 9 6 8 7 2 4 4 2 0 9 1 e - 1 0

We find the corresponding minimizer (0, 0, 0) by matching the moments of the moment matrix with a dirac measure centered at this minimizer.

extractatoms(ν, 1e-6)

Atomic measure on the variables x[1], x[2], x[3] with 1 atoms:
 at [2.0703586132122684e-6, 1.2367899143611297e-6, 2.2362354270007176e-8] with weight 1.0000000014903563

We can see below that the basis contained 6 monomials hence we needed to use 6x6 PSD matrix variables.

ν.basis

MonomialBasis{DynamicPolynomials.Monomial{true},DynamicPolynomials.MonomialVector{true}}(DynamicPolynomials.Monomial{true}[x₁x₂, x₂x₃, x₁, x₂, x₃, 1])

Using the monomial/term sparsity method of [WML20a] based on cluster completion, we find the same bound.

bound, ν = sos_min(MonomialSparsity())
bound

3 . 0 4 5 1 9 6 8 7 2 4 4 2 0 9 1 e - 1 0

Which is not suprising as no sparsity reduction could be performed.

[sub.basis for sub in ν.sub_moment_matrices]

1-element Array{MonomialBasis{DynamicPolynomials.Monomial{true},DynamicPolynomials.MonomialVector{true}},1}:
 MonomialBasis{DynamicPolynomials.Monomial{true},DynamicPolynomials.MonomialVector{true}}(DynamicPolynomials.Monomial{true}[x₁x₂, x₂x₃, x₁, x₂, x₃, 1])

Using the monomial/term sparsity method of [WML20b] based on chordal completion, the lower bound is smaller than 0.

bound, ν = sos_min(MonomialSparsity(ChordalCompletion()))
bound

0 . 0 0 3 5 5 1 2 0 0 9 9 0 4 9 8 6 5 9 7

However, this bound was obtained with an SDP with 4 matrices of size 3x3.

[sub.basis for sub in ν.sub_moment_matrices]

4-element Array{MonomialBasis{DynamicPolynomials.Monomial{true},DynamicPolynomials.MonomialVector{true}},1}:
 MonomialBasis{DynamicPolynomials.Monomial{true},DynamicPolynomials.MonomialVector{true}}(DynamicPolynomials.Monomial{true}[x₂x₃, x₂, x₃])
 MonomialBasis{DynamicPolynomials.Monomial{true},DynamicPolynomials.MonomialVector{true}}(DynamicPolynomials.Monomial{true}[x₂x₃, x₂, 1])
 MonomialBasis{DynamicPolynomials.Monomial{true},DynamicPolynomials.MonomialVector{true}}(DynamicPolynomials.Monomial{true}[x₁, x₂, 1])
 MonomialBasis{DynamicPolynomials.Monomial{true},DynamicPolynomials.MonomialVector{true}}(DynamicPolynomials.Monomial{true}[x₁x₂, x₁, 1])

This page was generated using Literate.jl.