{
"cells": [
{
"cell_type": "markdown",
"source": [
"# Sign symmetry"
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"**Adapted from**: Example 4 of [L09]\n",
"\n",
"[L09] Lofberg, Johan.\n",
"*Pre-and post-processing sum-of-squares programs in practice*.\n",
"IEEE transactions on automatic control 54, no. 5 (2009): 1007-1011."
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "(DynamicPolynomials.Variable{DynamicPolynomials.Commutative{DynamicPolynomials.CreationOrder}, MultivariatePolynomials.Graded{MultivariatePolynomials.LexOrder}}[x₁, x₂, x₃],)"
},
"metadata": {},
"execution_count": 1
}
],
"cell_type": "code",
"source": [
"using DynamicPolynomials\n",
"@polyvar x[1:3]"
],
"metadata": {},
"execution_count": 1
},
{
"cell_type": "markdown",
"source": [
"We would like to determine whether the following polynomial is a sum-of-squares."
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "1 + x₃² + x₁x₂ + x₂⁴ + x₁⁴",
"text/latex": "$$ 1 + x_{3}^{2} + x_{1}x_{2} + x_{2}^{4} + x_{1}^{4} $$"
},
"metadata": {},
"execution_count": 2
}
],
"cell_type": "code",
"source": [
"poly = 1 + x[1]^4 + x[1] * x[2] + x[2]^4 + x[3]^2"
],
"metadata": {},
"execution_count": 2
},
{
"cell_type": "markdown",
"source": [
"In order to do this, we can solve the following Sum-of-Squares program."
],
"metadata": {}
},
{
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Success: SDP solved\n",
"Primal objective value: 0.0000000e+00 \n",
"Dual objective value: 0.0000000e+00 \n",
"Relative primal infeasibility: 5.23e-17 \n",
"Relative dual infeasibility: 5.00e-11 \n",
"Real Relative Gap: 0.00e+00 \n",
"XZ Relative Gap: 1.00e-10 \n",
"DIMACS error measures: 7.85e-17 0.00e+00 5.00e-11 0.00e+00 0.00e+00 1.00e-10\n",
"CSDP 6.2.0\n",
"This is a pure primal feasibility problem.\n",
"Iter: 0 Ap: 0.00e+00 Pobj: 0.0000000e+00 Ad: 0.00e+00 Dobj: 0.0000000e+00 \n",
"Iter: 1 Ap: 8.10e-01 Pobj: 0.0000000e+00 Ad: 1.00e+00 Dobj: 4.0878274e+01 \n",
"Iter: 2 Ap: 1.00e+00 Pobj: 0.0000000e+00 Ad: 1.00e+00 Dobj: 2.7225538e+01 \n",
"* Solver : CSDP\n",
"\n",
"* Status\n",
" Result count : 1\n",
" Termination status : OPTIMAL\n",
" Message from the solver:\n",
" \"Problem solved to optimality.\"\n",
"\n",
"* Candidate solution (result #1)\n",
" Primal status : FEASIBLE_POINT\n",
" Dual status : FEASIBLE_POINT\n",
" Objective value : 0.00000e+00\n",
" Dual objective value : 0.00000e+00\n",
"\n",
"* Work counters\n",
" Solve time (sec) : 4.61102e-04\n",
"\n"
]
},
{
"output_type": "execute_result",
"data": {
"text/plain": "7-element DynamicPolynomials.MonomialVector{DynamicPolynomials.Commutative{DynamicPolynomials.CreationOrder}, MultivariatePolynomials.Graded{MultivariatePolynomials.LexOrder}}:\n 1\n x₃\n x₂\n x₁\n x₂²\n x₁x₂\n x₁²"
},
"metadata": {},
"execution_count": 3
}
],
"cell_type": "code",
"source": [
"import CSDP\n",
"solver = CSDP.Optimizer\n",
"using SumOfSquares\n",
"function sos_check(sparsity)\n",
" model = Model(solver)\n",
" con_ref = @constraint(model, poly in SOSCone(), sparsity = sparsity)\n",
" optimize!(model)\n",
" println(solution_summary(model))\n",
" return gram_matrix(con_ref)\n",
"end\n",
"\n",
"g = sos_check(Sparsity.NoPattern())\n",
"g.basis.monomials"
],
"metadata": {},
"execution_count": 3
},
{
"cell_type": "markdown",
"source": [
"As detailed in the Example 4 of [L09], we can exploit the *sign symmetry* of\n",
"the polynomial to decompose the large positive semidefinite matrix into smaller ones."
],
"metadata": {}
},
{
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Success: SDP solved\n",
"Primal objective value: 0.0000000e+00 \n",
"Dual objective value: 0.0000000e+00 \n",
"Relative primal infeasibility: 7.46e-17 \n",
"Relative dual infeasibility: 5.00e-11 \n",
"Real Relative Gap: 0.00e+00 \n",
"XZ Relative Gap: 1.21e-10 \n",
"DIMACS error measures: 1.21e-16 0.00e+00 5.00e-11 0.00e+00 0.00e+00 1.21e-10\n",
"CSDP 6.2.0\n",
"This is a pure primal feasibility problem.\n",
"Iter: 0 Ap: 0.00e+00 Pobj: 0.0000000e+00 Ad: 0.00e+00 Dobj: 0.0000000e+00 \n",
"Iter: 1 Ap: 7.15e-01 Pobj: 0.0000000e+00 Ad: 1.00e+00 Dobj: 4.7634216e+01 \n",
"Iter: 2 Ap: 1.00e+00 Pobj: 0.0000000e+00 Ad: 1.00e+00 Dobj: 3.5035802e+01 \n",
"* Solver : CSDP\n",
"\n",
"* Status\n",
" Result count : 1\n",
" Termination status : OPTIMAL\n",
" Message from the solver:\n",
" \"Problem solved to optimality.\"\n",
"\n",
"* Candidate solution (result #1)\n",
" Primal status : FEASIBLE_POINT\n",
" Dual status : FEASIBLE_POINT\n",
" Objective value : 0.00000e+00\n",
" Dual objective value : 0.00000e+00\n",
"\n",
"* Work counters\n",
" Solve time (sec) : 6.58989e-04\n",
"\n"
]
},
{
"output_type": "execute_result",
"data": {
"text/plain": "3-element Vector{DynamicPolynomials.MonomialVector{DynamicPolynomials.Commutative{DynamicPolynomials.CreationOrder}, MultivariatePolynomials.Graded{MultivariatePolynomials.LexOrder}}}:\n [x₂, x₁]\n [x₃]\n [1, x₂², x₁x₂, x₁²]"
},
"metadata": {},
"execution_count": 4
}
],
"cell_type": "code",
"source": [
"g = sos_check(Sparsity.SignSymmetry())\n",
"monos = [sub.basis.monomials for sub in g.blocks]"
],
"metadata": {},
"execution_count": 4
},
{
"cell_type": "markdown",
"source": [
"---\n",
"\n",
"*This notebook was generated using [Literate.jl](https://github.com/fredrikekre/Literate.jl).*"
],
"metadata": {}
}
],
"nbformat_minor": 3,
"metadata": {
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.11.3"
},
"kernelspec": {
"name": "julia-1.11",
"display_name": "Julia 1.11.3",
"language": "julia"
}
},
"nbformat": 4
}