Refactor with StarAlgebras (#123)

* Refactor with StarAlgebras

* Fixes

* Fix format

* Fixes

Author: blegat

.github/workflows/ci.yml

@@ -35,6 +35,15 @@ jobs:
             ${{ runner.os }}-test-${{ env.cache-name }}-
             ${{ runner.os }}-test-
             ${{ runner.os }}-
+      - name: dev
+        shell: julia --project=@. {0}
+        run: |
+          using Pkg
+          Pkg.add([
+              PackageSpec(name="StarAlgebras", rev="mk/non_monomial_basis"),
+              PackageSpec(name="MultivariateBases", rev="master"),
+              PackageSpec(name="MultivariateMoments", rev="master"),
+          ])
       - uses: julia-actions/julia-buildpkg@v1
       - uses: julia-actions/julia-runtest@v1
       - uses: julia-actions/julia-processcoverage@v1

Project.toml

@@ -15,6 +15,7 @@ MultivariateMoments = "f4abf1af-0426-5881-a0da-e2f168889b5e"
 MultivariatePolynomials = "102ac46a-7ee4-5c85-9060-abc95bfdeaa3"
 MutableArithmetics = "d8a4904e-b15c-11e9-3269-09a3773c0cb0"
 SemialgebraicSets = "8e049039-38e8-557d-ae3a-bc521ccf6204"
+StarAlgebras = "0c0c59c1-dc5f-42e9-9a8b-b5dc384a6cd1"
 
 [compat]
 DataStructures = "0.18"

src/Bridges/Constraint/Constraint.jl

@@ -1,5 +1,6 @@
 module Constraint
 
+import StarAlgebras as SA
 import MultivariatePolynomials as MP
 import MultivariateBases as MB
 import SemialgebraicSets as SS

src/Bridges/Constraint/zero_polynomial.jl

@@ -1,23 +1,21 @@
 struct ZeroPolynomialBridge{
     T,
     F<:MOI.AbstractVectorFunction,
-    MT<:MP.AbstractMonomial,
-    MVT<:AbstractVector{MT},
+    B<:SA.ExplicitBasis,
 } <: MOI.Bridges.Constraint.AbstractBridge
     zero_constraint::MOI.ConstraintIndex{F,MOI.Zeros}
-    monomials::MVT
+    basis::B
 end
 
 function MOI.Bridges.Constraint.bridge_constraint(
-    ::Type{ZeroPolynomialBridge{T,F,MT,MVT}},
+    ::Type{ZeroPolynomialBridge{T,F,B}},
     model::MOI.ModelLike,
     f::MOI.AbstractVectorFunction,
-    s::PolyJuMP.ZeroPolynomialSet{SS.FullSpace,<:MB.MonomialBasis},
-) where {T,F,MT,MVT}
-    @assert MOI.output_dimension(f) == length(s.monomials)
-    zero_constraint =
-        MOI.add_constraint(model, f, MOI.Zeros(length(s.monomials)))
-    return ZeroPolynomialBridge{T,F,MT,MVT}(zero_constraint, s.monomials)
+    s::PolyJuMP.ZeroPolynomialSet{SS.FullSpace,<:MB.FullBasis{MB.Monomial},B},
+) where {T,F,B}
+    @assert MOI.output_dimension(f) == length(s.basis)
+    zero_constraint = MOI.add_constraint(model, f, MOI.Zeros(length(s.basis)))
+    return ZeroPolynomialBridge{T,F,typeof(s.basis)}(zero_constraint, s.basis)
 end
 
 function MOI.supports_constraint(
@@ -41,10 +39,14 @@ function MOI.Bridges.Constraint.concrete_bridge_type(
     ::Type{<:ZeroPolynomialBridge{T}},
     F::Type{<:MOI.AbstractVectorFunction},
     ::Type{
-        <:PolyJuMP.ZeroPolynomialSet{SS.FullSpace,<:MB.MonomialBasis,MT,MVT},
+        <:PolyJuMP.ZeroPolynomialSet{
+            SS.FullSpace,
+            <:MB.FullBasis{MB.Monomial},
+            B,
+        },
     },
-) where {T,MT,MVT}
-    return ZeroPolynomialBridge{T,F,MT,MVT}
+) where {T,B<:MB.SubBasis{MB.Monomial}}
+    return ZeroPolynomialBridge{T,F,B}
 end
 
 # Attributes, Bridge acting as an model
@@ -68,14 +70,14 @@ end
 
 # Attributes, Bridge acting as a constraint
 function MOI.get(
-    model::MOI.ModelLike,
+    ::MOI.ModelLike,
     ::MOI.ConstraintSet,
-    bridge::ZeroPolynomialBridge{T,F,MT,MVT},
-) where {T,F,MT,MVT}
+    bridge::ZeroPolynomialBridge{T,F,MB.SubBasis{MB.Monomial,M,V}},
+) where {T,F,M,V}
     return PolyJuMP.ZeroPolynomialSet(
         SS.FullSpace(),
-        MB.MonomialBasis{MT,MVT},
-        bridge.monomials,
+        MB.FullBasis{MB.Monomial,M}(),
+        bridge.basis,
     )
 end
 function MOI.get(
@@ -91,5 +93,5 @@ function MOI.get(
     bridge::ZeroPolynomialBridge,
 )
     values = MOI.get(model, MOI.ConstraintDual(attr.result_index), bridge)
-    return MM.measure(values, bridge.monomials)
+    return MM.measure(values, bridge.basis)
 end

src/Bridges/Constraint/zero_polynomial_in_algebraic_set.jl

@@ -3,26 +3,34 @@ using LinearAlgebra
 struct ZeroPolynomialInAlgebraicSetBridge{
     T,
     F<:MOI.AbstractVectorFunction,
-    BT<:MB.AbstractPolynomialBasis,
+    Z<:SA.AbstractBasis,
     DT<:SS.AbstractAlgebraicSet,
     MT<:MP.AbstractMonomial,
     MVT<:AbstractVector{MT},
 } <: MOI.Bridges.Constraint.AbstractBridge
     zero_constraint::MOI.ConstraintIndex{
         F,
-        PolyJuMP.ZeroPolynomialSet{SS.FullSpace,BT,MT,MVT},
+        PolyJuMP.ZeroPolynomialSet{
+            SS.FullSpace,
+            Z,
+            MB.SubBasis{MB.Monomial,MT,MVT},
+        },
     }
     domain::DT
-    monomials::MVT
+    basis::MB.SubBasis{MB.Monomial,MT,MVT}
 end
 
 function MOI.Bridges.Constraint.bridge_constraint(
-    ::Type{ZeroPolynomialInAlgebraicSetBridge{T,F,BT,DT,MT,MVT}},
+    ::Type{ZeroPolynomialInAlgebraicSetBridge{T,F,Z,DT,MT,MVT}},
     model::MOI.ModelLike,
     f::MOI.AbstractVectorFunction,
-    s::PolyJuMP.ZeroPolynomialSet{<:SS.AbstractAlgebraicSet},
-) where {T,F,BT,DT,MT,MVT}
-    p = MP.polynomial(MOI.Utilities.scalarize(f), s.monomials)
+    s::PolyJuMP.ZeroPolynomialSet{
+        <:SS.AbstractAlgebraicSet,
+        Z,
+        MB.SubBasis{MB.Monomial,MT,MVT},
+    },
+) where {T,F,Z,DT,MT,MVT}
+    p = MP.polynomial(MOI.Utilities.scalarize(f), s.basis)
     # As `*(::MOI.ScalarAffineFunction{T}, ::S)` is only defined if `S == T`, we
     # need to call `similar`. This is critical since `T` is
     # `Float64` when used with JuMP and the coefficient type is often `Int` with
@@ -32,12 +40,16 @@ function MOI.Bridges.Constraint.bridge_constraint(
     zero_constraint = MOI.add_constraint(
         model,
         MOI.Utilities.vectorize(MP.coefficients(r)),
-        PolyJuMP.ZeroPolynomialSet(SS.FullSpace(), s.basis, MP.monomials(r)),
+        PolyJuMP.ZeroPolynomialSet(
+            SS.FullSpace(),
+            s.zero_basis,
+            MB.SubBasis{MB.Monomial}(MP.monomials(r)),
+        ),
     )
-    return ZeroPolynomialInAlgebraicSetBridge{T,F,BT,DT,MT,MVT}(
+    return ZeroPolynomialInAlgebraicSetBridge{T,F,Z,DT,MT,MVT}(
         zero_constraint,
         s.domain,
-        s.monomials,
+        s.basis,
     )
 end
 
@@ -48,42 +60,60 @@ function MOI.supports_constraint(
 ) where {T}
     return true
 end
+
 function MOI.Bridges.added_constrained_variable_types(
     ::Type{<:ZeroPolynomialInAlgebraicSetBridge},
 )
     return Tuple{Type}[]
 end
+
 function MOI.Bridges.added_constraint_types(
-    ::Type{<:ZeroPolynomialInAlgebraicSetBridge{T,F,BT,DT,MT,MVT}},
-) where {T,F,BT,DT,MT,MVT}
-    return [(F, PolyJuMP.ZeroPolynomialSet{SS.FullSpace,BT,MT,MVT})]
+    ::Type{<:ZeroPolynomialInAlgebraicSetBridge{T,F,Z,DT,MT,MVT}},
+) where {T,F,Z,DT,MT,MVT}
+    return [(
+        F,
+        PolyJuMP.ZeroPolynomialSet{
+            SS.FullSpace,
+            Z,
+            MB.SubBasis{MB.Monomial,MT,MVT},
+        },
+    )]
 end
+
 function MOI.Bridges.Constraint.concrete_bridge_type(
     ::Type{<:ZeroPolynomialInAlgebraicSetBridge{T}},
     F::Type{<:MOI.AbstractVectorFunction},
-    ::Type{<:PolyJuMP.ZeroPolynomialSet{DT,BT,MT,MVT}},
-) where {T,BT,DT<:SS.AbstractAlgebraicSet,MT,MVT}
+    ::Type{<:PolyJuMP.ZeroPolynomialSet{DT,Z,MB.SubBasis{MB.Monomial,MT,MVT}}},
+) where {T,Z,DT<:SS.AbstractAlgebraicSet,MT,MVT}
     G = MOI.Utilities.promote_operation(-, T, F, F)
-    return ZeroPolynomialInAlgebraicSetBridge{T,G,BT,DT,MT,MVT}
+    return ZeroPolynomialInAlgebraicSetBridge{T,G,Z,DT,MT,MVT}
 end
 
 # Attributes, Bridge acting as an model
 function MOI.get(
-    ::ZeroPolynomialInAlgebraicSetBridge{T,F,BT,DT,MT,MVT},
+    ::ZeroPolynomialInAlgebraicSetBridge{T,F,Z,DT,MT,MVT},
     ::MOI.NumberOfConstraints{
         F,
-        PolyJuMP.ZeroPolynomialSet{SS.FullSpace,BT,MT,MVT},
+        PolyJuMP.ZeroPolynomialSet{
+            SS.FullSpace,
+            Z,
+            MB.SubBasis{MB.Monomial,MT,MVT},
+        },
     },
-) where {T,F,BT,DT,MT,MVT}
+) where {T,F,Z,DT,MT,MVT}
     return 1
 end
 function MOI.get(
-    b::ZeroPolynomialInAlgebraicSetBridge{T,F,BT,DT,MT,MVT},
+    b::ZeroPolynomialInAlgebraicSetBridge{T,F,Z,DT,MT,MVT},
     ::MOI.ListOfConstraintIndices{
         F,
-        PolyJuMP.ZeroPolynomialSet{SS.FullSpace,BT,MT,MVT},
+        PolyJuMP.ZeroPolynomialSet{
+            SS.FullSpace,
+            Z,
+            MB.SubBasis{MB.Monomial,MT,MVT},
+        },
     },
-) where {T,F,BT,DT,MT,MVT}
+) where {T,F,Z,DT,MT,MVT}
     return [b.zero_constraint]
 end
 
@@ -102,8 +132,8 @@ function MOI.get(
     set = MOI.get(model, attr, bridge.zero_constraint)
     return PolyJuMP.ZeroPolynomialSet(
         bridge.domain,
-        set.basis,
-        bridge.monomials,
+        set.zero_basis,
+        bridge.basis,
     )
 end
 # TODO ConstraintPrimal
@@ -121,10 +151,11 @@ function MOI.get(
 )
     dual = MOI.get(model, attr, bridge.zero_constraint)
     set = MOI.get(model, MOI.ConstraintSet(), bridge.zero_constraint)
-    μ = MM.measure(dual, set.monomials)
+    μ = MM.measure(dual, set.basis)
     I = SS.ideal(bridge.domain)
-    return [dot(rem(mono, I), μ) for mono in bridge.monomials]
+    return [dot(rem(mono, I), μ) for mono in bridge.basis.monomials]
 end
+
 function MOI.get(
     model::MOI.ModelLike,
     attr::PolyJuMP.MomentsAttribute,

src/PolyJuMP.jl

@@ -1,6 +1,7 @@
 module PolyJuMP
 
 import MutableArithmetics as MA
+import StarAlgebras as SA
 
 # MultivariatePolynomials extension
 

src/SAGE/SAGE.jl

@@ -1,6 +1,7 @@
 module SAGE
 
 import MutableArithmetics as MA
+import MultivariateBases as MB
 import MultivariatePolynomials as MP
 import MathOptInterface as MOI
 import JuMP
@@ -40,9 +41,10 @@ struct Signomials{M<:Union{Nothing,Int,MP.AbstractMonomial}} <:
 end
 Signomials() = Signomials(nothing)
 _index(_, ::Nothing) = nothing
-_index(monos, mono::MP.AbstractMonomial) = findfirst(isequal(mono), monos)::Int
-function JuMP.moi_set(c::Signomials, monos)
-    return Cone(Signomials(_index(monos, c.monomial)), _exponents_matrix(monos))
+_index(basis, mono::MP.AbstractMonomial) = MB.monomial_index(basis, mono)::Int
+function JuMP.moi_set(c::Signomials, basis::MB.SubBasis{MB.Monomial})
+    monos = basis.monomials
+    return Cone(Signomials(_index(basis, c.monomial)), _exponents_matrix(monos))
 end
 
 """
@@ -55,9 +57,10 @@ struct Polynomials{M<:Union{Nothing,Int,MP.AbstractMonomial}} <:
     monomial::M
 end
 Polynomials() = Polynomials(nothing)
-function JuMP.moi_set(c::Polynomials, monos)
+function JuMP.moi_set(c::Polynomials, basis::MB.SubBasis{MB.Monomial})
+    monos = basis.monomials
     return Cone(
-        Polynomials(_index(monos, c.monomial)),
+        Polynomials(_index(basis, c.monomial)),
         _exponents_matrix(monos),
     )
 end
@@ -77,9 +80,9 @@ function JuMP.build_constraint(
         _error("unsupported keyword argument `$key`.")
     end
     coefs = PolyJuMP.non_constant_coefficients(p)
-    monos = MP.monomials(p)
-    cone = JuMP.moi_set(set, monos)
-    shape = PolyJuMP.PolynomialShape(monos)
+    basis = MB.SubBasis{MB.Monomial}(MP.monomials(p))
+    cone = JuMP.moi_set(set, basis)
+    shape = PolyJuMP.PolynomialShape(basis)
     return PolyJuMP.bridgeable(
         JuMP.VectorConstraint(coefs, cone, shape),
         JuMP.moi_function_type(typeof(coefs)),
@@ -136,7 +139,7 @@ function decomposition(
     tol::Real,
     result_index::Int = 1,
 )
-    monos = con_ref.shape.monomials
+    monos = con_ref.shape.basis.monomials
     attr = DecompositionAttribute(tol, result_index)
     return Decomposition([
         MP.polynomial(a, monos) for

src/attributes.jl

@@ -1,5 +1,5 @@
 """
-    MomentsAttribute(N)
+    MomentsAttribute(result_index)
     MomentsAttribute()
 
 A constraint attribute for the vector of moments corresponding to the