Symmetry & SOS Optimization Compatibility Issues?

Hello all,

I am trying to implement rotational symmetry for an SOS problem. In particular, I followed the example problem here: Dihedral symmetry of the Robinson form · SumOfSquares.

My goal is to use a cyclic group instead for a simple rotation symmetry, I don’t need to implement the rotation piece, and this resulted in the following code.

using DynamicPolynomials
using SumOfSquares
using Mosek
using MosekTools
using PermutationGroups
using GroupsCore
using SymbolicWedderburn

include("/Users/elizabethc/Documents/Code/JULIA_SOS_Couette/CyclicGroup.jl")
struct CyclicAction <: Symmetry.OnMonomials end
SymbolicWedderburn.coeff_type(::CyclicAction) = Float64

@polyvar a[1:6]
poly3 = sum(a.^2) 

function SymbolicWedderburn.action(::CyclicAction, el::CyclicGroupElement, mono::AbstractMonomial)
       var_a3, var_a4 = a[4], a[3]; var_a5, var_a6 = a[6], a[5]
       sign_aodd = 1 <= el.residual <=2 ? -1 : 1
       sign_aeven = 2 <= el.residual ? -1 : 1
       return mono([a[1], a[2], a[3], a[4], a[5], a[6]] => [a[1], a[2], sign_aodd*var_a3, sign_aeven*var_a4, sign_aodd*var_a5, sign_aeven*var_a6])
       end
G = CyclicGroup(4)
pattern3 = Symmetry.Pattern(G, CyclicAction())

model3 = SOSModel(Mosek.Optimizer)
@variable(model3, α3)
con_ref3 = @constraint(model3, poly3 >= α3, symmetry = pattern3)
optimize!(model3)
@show termination_status(model3)
@show objective_value(model3)

This yields an infeasible result using the optimization toolbox Mosek, and an empty constraint result using the optimization toolbox CSDP, which makes no sense since polynomial poly is sum-of-squares by definition But, if I check that the symmetry is being correctly implemented using the following for loop

for g in G
    @show SymbolicWedderburn.action(CyclicAction(), g, poly)
end

it shows that the polynomial is invariant under the claimed rotation, as it should be. We can check the proper symmetry is being imposed with lin = sum(a).

To see if I understood some structure incorrectly, I went back to the dihedral Robinson from example from the JuMP/SumOfSquares documentation, and thought to implement the exact same method, but instead of the polynomial being the Robsinson form p(x,y), I generalized it so that we have the new polynomial P(x_1, x_2,y_1,y_2) = p(x_1,y_1) + p(x_2,y_2). The Wedderburn action should be the same but now adjusted to make sure that the rotations and reflections are happening to (x_i,y_i), a very easy adjustment, and easily checked again by the for loop. If I then attempt to run the optimization model, I receive the following error, with a corresponding stacktrace:

ERROR: No permutation can make [[0.0 1.0; -1.0 0.0], [0.0 -1.0; -1.0 0.0]] and [[0.0 1.0; -1.0 0.0], [0.0 1.0; 1.0 0.0]] match

If I change the polynomial to a sum-of-squares, it declares feasibility and then gives me an optimal value t=-7.740740..., which is completely incorrect it should be 0. I can’t make heads nor tails of these errors, and the only person I know in my group who has experience with this is hypothesizing that it is actually a symmetry package compatibility issue with SOS. Any insight on what the issue might be would be greatly appreciated.

1 Like

Welcome to the forum @ilovemath, this is a question for @blegat.

2 Likes

@ilovemath I’m behind the symmetry reduction; is there a minimal runnable example you could provide?
Is the CyclicGroup.jl the same as

?

@blegat Is this really an action OnMonomials? does transposition x → -x really fit into this?
Here:

I’ve rewritten the dihedral action as BySignedPermutations where
action(act, g, basis_elt) returns (new_basis_elt, ±1).

2 Likes

Thanks @abulak. Yes, that is the correct script. With this script, my existing example is close to minimal. Here is one that is even more minimal:

using DynamicPolynomials
using SumOfSquares
using Mosek
using MosekTools
using PermutationGroups
using GroupsCore
using SymbolicWedderburn
using LinearAlgebra

include("/Users/elizabethc/Documents/Code/JULIA_SOS_Couette/CyclicGroup.jl")
struct CyclicAction <: Symmetry.OnMonomials end
SymbolicWedderburn.coeff_type(::CyclicAction) = Float64

@polyvar a[1:2]
function SymbolicWedderburn.action(::CyclicAction, el::CyclicGroupElement, mono::AbstractMonomial)
                     var_a1, var_a2 = a[2], a[1]
                     sign_aodd = 1 <= el.residual <=2 ? -1 : 1
                     sign_aeven = 2 <= el.residual ? -1 : 1
                     return mono([a[1], a[2]] => [sign_aodd*var_a1, sign_aeven*var_a2])
end
G = CyclicGroup(4)
pattern = Symmetry.Pattern(G, CyclicAction())
poly = a[1]^2 + a[2]^2
model = SOSModel(Mosek.Optimizer)
@variable(model, t)
@objective(model, Max, t)
con_ref = @constraint(model, poly -t >=0, symmetry = pattern)
optimize!(model)

This code yields an infeasible result. Not a clue why. (It is late, there could be a typo, will double check in the morning.)

As to the dihedral example, I went through trying to understand the exact differences here for the Robinson example. With the adjustment to use BySignedPermutation as the generator for the pattern for the SOS optimization problem for the Robinson form example, it throws the same error as I mentioned in my post when I generalize to p(x_1, y_1) + p(x_2, y_2) (same notation as in og post). So enforcing the Symmetry.pattern either with DihedralAction or DihedralActionSP does not seem to matter, which I’m not entirely surprised by.

P.S. The first struct for DihedralAction will throw an error if OnMonomials is not written as Symmetry.OnMonomials. I’m not including the files in lines 7-9 so unless it is specified in there that the prefix is unnecessary this is a typo.

@ilovemath it seems to me that the provided action function does not define a group action.
Could you make sure that the following tests pass?

G = CyclicGroup(4)
act = CyclicAction()
basis = monomials(a, 0:2)

using Test
@testset "Action axioms" begin
    SW = SymbolicWedderburn
    for b in basis
        # the group identity acts as the identity function
        @test SW.action(act, one(G), b) == b
    end
    for b in basis
        for g in G
            for h in G
                # the action is right-associative, i.e. b^(g*h) = (b^g)^h
                @test SW.action(act, g * h, b) ==
                      SW.action(act, h, SW.action(act, g, b))
            end
        end
    end
end
2 Likes

Hi @ilovemath, thanks for reporting this issue. We are adding automatic checks to make sure you get a nice error if your action is not a group homomorphism in add checks for group action correctness by kalmarek · Pull Request #56 · kalmarek/SymbolicWedderburn.jl · GitHub

1 Like

Thanks all, this was helpful. The problem ended up being extremely minor and a bit dumb. For future interested parties, here is the correct version of the minimal example:

using DynamicPolynomials
using SumOfSquares
using Mosek
using MosekTools
using PermutationGroups
using GroupsCore
using SymbolicWedderburn
using LinearAlgebra

include("/Users/elizabethc/Documents/Code/JULIA_SOS_Couette/CyclicGroup.jl")
struct CyclicAction <: Symmetry.OnMonomials end
SymbolicWedderburn.coeff_type(::CyclicAction) = Float64

@polyvar a[1:2]
function SymbolicWedderburn.action(::CyclicAction, el::CyclicGroupElement, mono::AbstractMonomial)
                     if isodd(el.residual)
                         var_a1, var_a2 = a[2], a[1]
                     else
                         var_a1, var_a2 = a[1], a[2]
                     end
                     sign_aodd = 1 <= el.residual <=2 ? -1 : 1
                     sign_aeven = 2 <= el.residual ? -1 : 1
                     return mono([a[1], a[2]] => [sign_aodd*var_a1, sign_aeven*var_a2])
end
G = CyclicGroup(4)
pattern = Symmetry.Pattern(G, CyclicAction())
poly = a[1]^2 + a[2]^2
model = SOSModel(Mosek.Optimizer)
@variable(model, t)
@objective(model, Max, t)
con_ref = @constraint(model, poly -t >=0, symmetry = pattern)
optimize!(model)