Type unstability in SymbolicWedderburn

[sebastiendesignolle]

Given an $m\times m$ matrix $M$, I'm interested in solving relaxations of the problem

$$\max\left\{\sum_{x,y=1}^mM_{xy}a_xb_y~|~a_1,\ldots,a_m,b_1,\ldots,b_m=\pm1\right\},$$

in particular when the matrix $M$ is highly symmetric.

Here is the way I'm coding it at the moment, but this leads to a lot of unstability in the symmetry reduction, see the attached image.
Am I missing something simple?

import COSMO
import Dualization
import DynamicPolynomials as DP
import JuMP
import LinearAlgebra as LA
import PermutationGroups as PG # v0.4.3
import Polymake
import SemialgebraicSets
import SumOfSquares as SOS # v0.7.3
# import SymbolicWedderburn # v0.3.7

function test(
    M::Matrix,
    gens::Vector{Vector{Int}} = Vector{Int}[];
    solver = COSMO.Optimizer,
    level = 1
)
    m1, m2 = size(M)
    DP.@polyvar x[1:m1+m2]
    a = x[1:m1]
    b = x[m1+1:m1+m2]
    f = LA.dot(a, M, b)
    S = SemialgebraicSets.intersect([SemialgebraicSets.@set x[i]^2 == 1 for i in 1:m1+m2]...)

    model = SOS.SOSModel(Dualization.dual_optimizer(solver))
    JuMP.@variable(model, t)
    JuMP.@objective(model, Min, t)

    if isempty(gens)
        pattern = nothing
    else
        @assert m1 == m2
        m = m1
        generators = vcat([PG.Perm{UInt16}([g; g .+ m]) for g in gens], [PG.Perm{UInt16}([m+1:2m; 1:m])])
        G = PG.PermGroup(generators)
        pattern = SOS.Symmetry.Pattern(G, SOS.Symmetry.VariablePermutation())
    end

    JuMP.@constraint(model, c, f ≤ t; domain = S, symmetry = pattern, maxdegree = 2level)
    JuMP.set_silent(model)
    JuMP.optimize!(model)
    return JuMP.objective_value(model)
end

# 24-cell as a simple example
pol = Polymake.polytope.wythoff("D4", 1)
# vertices
v = Matrix{Int}(pol.VERTICES)[:, 2:end]
# Gram matrix
p = v * v'
# generators of the action on the vertices, i.e., on the lines/columns of p
gens = Polymake.to_one_based_indexing(pol.GROUP.VERTICES_ACTION.GENERATORS)

# @profview of the following line gives the image below
test(p, gens; level = 1)

Also, I'm just starting working with SOS.jl so any comments on my syntax is welcome.
Thanks in advance!

[kalmarek]

this is sure a pretty autumn mountains landscape!

[kalmarek]

here's a smaller mwe:

G, S = let m1
    m1 = 3
    M = ones(Int, m1, m1)

    DP.@polyvar x[1:m1]
    a = x[1:m1]
    f = LA.dot(a, M, a)
    S = [SemialgebraicSets.@set x[i]^2 == 1 for i in 1:m1]
    generators = [PG.perm"(1,2,3)", PG.perm"(1,2)"]

    let p = PG.gens(G, 1), v = MP.variables(S[1])
        subs = v => SOS.Symmetry._map_idx(i -> v[i^p], v)
        @code_warntype MP.pair_zip(MP._monomial_vector_to_variable_tuple(subs))
    end

    G, S
end

warntype returns Any here

[blegat]

@kalmarek thanks for reducing it, the issue is clear now. This issue seems to have been introduced by JuliaAlgebra/MultivariatePolynomials.jl#268

[blegat]

@sebastiendesignolle should be fixed by DynamicPolynomials v0.6.1, let me know if that's not the case

[sebastiendesignolle]

Just run it, and the profile looks way less colourful, thanks a lot!