Creating bridges with JuMP.jl and MOI.jl

Happy new year everyone and welcome back to my blog!

In this post you'll learn something about the bridge system in JuMP/MOI and how to create your own simple bridges.

Disclaimer: Basically all of this was created by trial and error which means there are probably better ways to do this. This post might be some guidance for doing the same for your project until the proper way is formally documented 😉

I'll try to keep this up to date when I learn new stuff about the system myself.

This post is part of my long long series about how to write a constraint solver from almost scratch. Feel free to start at the beginning and come back to this post in a year 😄 or just read this as it's relatively independent of the rest.

Additionally feel free to checkout the sidebar of this blog to see all posts in this series as well as other posts that might be of interest to you.

What are bridges?

If you are interested in the full details please checkout this paper: "MathOptInterface: a data structure for mathematical optimization problems". Mathieu Besançon also wrote a nice high level overview about bridges and JuMP about a year ago on his blog Bridges as an extended dispatch system.

In short form bridges allow the user of a solver to write down the problem in a natural form like:

\[ \operatorname{min} 2x + 7 \\ 5x - 17 \geq 10 \]

Even though the solver might not be able to directly support these types of constraints and objectives. Often a solver supports a so called standard form which doesn't really have to be THE standard.

Let's assume the solver itself only supports maximizing a single variable and affine equality constraints.

Then the above problem can be expressed as:

\[ \begin{aligned} \operatorname{max} & -y \\ 2x + 7 &= y \\ -5x + 17 + z &= -10 \\ z &\geq 0 \\ \end{aligned} \]

These reformulations are always identical such that it doesn't make sense to do these on the solver level. It's also not always such an easy task as it might seem.

Furthermore, it's obvious that the modeler doesn't want to do the reformulations as it makes it harder to write and read.

This gap between the model and the solver is bridged by JuMP/MOI with bridges 😉

There are standard bridges which are always ready to be used like

\[ \geq \quad \underset{\cdot -1}{\Longleftrightarrow} \quad \leq \]

some others might make it too complicated and should only be used if explicitly requested.

A bit more detail

There are some critical parts to make this work

• The solver package must specify which constraints/variables/objectives are supported

• JuMP/MOI must decide which bridge to use

The first one is rather easier and you might have seen this in an earlier posts where I explain how I made ConstraintSolver compatible with JuMP: How to create a JuMP solver.

As an example:

const = MOI.ScalarAffineFunction{Float64}
MOI.supports_constraint(
    ::Optimizer, 
    ::Type{SAF}, 
    ::Type{MOI.EqualTo{Float64}}) = true

The supports_constraint function is called to determine whether a constraint is supported or not. The developer of the solver needs to define methods with the correct types to specify what it allows as an input.

A harder part is done by JuMP/MOI itself as there are normally several bridges possible. This is in contrast to multiple dispatch where methods need to be unambiguous. Internally a graph is created which represents all the possible ways of creating a supported model, using bridges, as paths in the graph. Afterwards one of the shortest paths is used as the bridge.

Here shortest paths basically means that we want to create as few constraints/variables as possible. Again more detail can be read in the paper.

Which bridge do I want to create?

Before I want to show you how to create your own bridges let me start with explaining what I personally want to do with them.

My solver does implement == and <= but not >=. Now as mentioned before there is a bridge which converts from >= to <= such that this is normally not a problem. Unfortunately this currently doesn't work when >= is used as the inner part of an indicator constraint i.e:

@constraint(m, b => { 2x-y >= 7 })

is supported in my ConstraintSolver by using a hack in all versions \(\leq\) 0.6.0 namely by bridging it myself when the indicator constraint is created.

if ASS isa Type{MOI.GreaterThan{T}}	
    inner_terms = [	
        MOI.ScalarAffineTerm(-v.scalar_term.coefficient, v.scalar_term.variable_index) for v in func.terms if v.output_index == 2	
    ]	
    inner_constant = -inner_constant	
    inner_set = MOI.LessThan{T}(-set.set.lower)	
end

where ASS the AbstractScalarSet.

That kinda works fine but I wanted to be more fancy in writing a bridge for this.

At the current stage the bridge isn't really that fancy as it is still very much specialized for my problem but I think I'm getting there.

How to create a bridge?

First of all we need to create a struct which holds the constraint index (CI) and therefore also the type of the constraint. The type of a constraint consists of two types the function type and the set type: i.e \(2x+7 <= 7\) can be separated into \(2x+7\) as the function of type ScalarAffineFunction and the set is <= 7 which is of type GreaterThan.

struct IndicatorGreaterToLessBridge{
    T,
    F<:MOI.VectorAffineFunction,
    G<:MOI.VectorAffineFunction,
    A
} <: MOIBC.FlipSignBridge{T,MOI.IndicatorSet{A},MOI.IndicatorSet{A},F,G}
    constraint::CI{F,MOI.IndicatorSet{A,MOI.LessThan{T}}}
end

This is the new struct for our bridge which is mostly copied from the GreaterToLessBridge that already exists in MOI.

Next up is the part where we communicate to MOI that we support a new type of constraint namely an indicator constraint of the form b => { ... >= .}.

function MOI.supports_constraint(
    ::Type{<:IndicatorGreaterToLessBridge{T}},
    ::Type{<:MOI.VectorAffineFunction},
    ::Type{IS},
) where {T,A,IS<:MOI.IndicatorSet{A,MOI.GreaterThan{T}}}
    return true
end

There are some more steps we need to do to actually make this work as we currently support that constraint without a bridge.

First of all we need to disallow that constraint without the bridge by changing the corresponding supports_constraint method:

function MOI.supports_constraint(
    ::Optimizer,
    func::Type{VAF{T}},
    set::Type{IS},
) where {A,T<:Real,ASS<:MOI.AbstractScalarSet,IS<:MOI.IndicatorSet{A,ASS}}
    if ASS <: MOI.GreaterThan
        return false
    end
    return A == MOI.ACTIVATE_ON_ONE || A == MOI.ACTIVATE_ON_ZERO
end

This returns false now for indicator constraints with a GreaterThan set.

Finally we need to add the bridge to be considered. I do this in the Optimizer function which creates the optimizer itself. It now returns a LazyBridgeOptimizer (lbo).

lbo = MOIB.full_bridge_optimizer(optimizer, options.solution_type)
MOIB.add_bridge(lbo, CS.ReifiedGreaterToLessBridge{options.solution_type})
MOIB.add_bridge(lbo, CS.IndicatorGreaterToLessBridge{options.solution_type})
return lbo

For convenience I've defined:

const MOIB = MathOptInterface.Bridges
const MOIBC = MathOptInterface.Bridges.Constraint

Few more things to do:

• Tell MOI which new variables and/or constraints are created by this bridge

  • This is important for finding the best bridge as mentioned earlier

• Define which bridge constraint type the bridge returns

• Define the bridge transformation itself

function MOIB.added_constrained_variable_types(::Type{<:IndicatorGreaterToLessBridge})
    return []
end
function MOIB.added_constraint_types(
    ::Type{<:IndicatorGreaterToLessBridge{T,F,G,A}},
) where {T,F,G,A}
    return [(F, MOI.IndicatorSet{A,MOI.LessThan{T}})]
end

These two methoids define which new constraint/variable we create. It's in an array as some bridges create several new constraints. We don't need any additional variables -> we can return an empty array in added_constrained_variable_types.

Now on to the concrete type of the bridge we want to return. This is needed as our struct from before is a parametric struct.

function MOIBC.concrete_bridge_type(
    ::Type{<:IndicatorGreaterToLessBridge{T}},
    ::Type{<:MOI.VectorAffineFunction},
    ::Type{IS},
) where {T,A,IS<:MOI.IndicatorSet{A,MOI.GreaterThan{T}}}
    return IndicatorGreaterToLessBridge{T,MOI.VectorAffineFunction{T},MOI.VectorAffineFunction{T},A}
end

Finally we want to specify the bridge itself by transforming the function and the set.

Let's start with the set as it's simpler:

function MOIBC.map_set(::Type{<:IndicatorGreaterToLessBridge}, set::MOI.IndicatorSet{A,<:MOI.GreaterThan}) where A
    inner_set = set.set
    return MOI.IndicatorSet{A}(MOI.LessThan(-inner_set.lower))
end
function MOIBC.inverse_map_set(::Type{<:IndicatorGreaterToLessBridge}, set::MOI.IndicatorSet{A,<:MOI.LessThan}) where A
    inner_set = set.set
    return MOI.IndicatorSet{A}(MOI.GreaterThan(-inner_set.upper))
end

The first one transforms the set from GreaterThan to LessThan by having a look at the set of the set (what I define as inner_set) and then multiply it with \(-1\). The inverse_map_set is the same but the other way around. This is needed for things like ConstraintPrimal which I haven't used yet but it's always helpful to have the function implemented already.

This will now work in the sense that there will be no methods missing and the code will run through as expected. This is due to the fact that the FlipSignBridge that we extend already has a fallback function for mapping the function.

In our case however this will transform the VectorAffineFunction for both the indicator variable and the affine constraint.

A vector affine function consists of the following:

help?> MOI.VectorAffineFunction
  VectorAffineFunction{T}(terms, constants)

  The vector-valued affine function A x + b, where:

    •    A is a sparse matrix specified by a list of VectorAffineTerm objects.

    •    b is a vector specified by constants

  Duplicate indices in the A are accepted, and the corresponding coefficients are summed together.

The terms have an index and a ScalarAffineTerm which get summed together if they have the same index.

Something like: \(b => { 2x + y + 7 >= 10}\) is then:

MathOptInterface.VectorAffineFunction{Float64}(
    MathOptInterface.VectorAffineTerm{Float64}[
        MathOptInterface.VectorAffineTerm{Float64}(1, MathOptInterface.ScalarAffineTerm{Float64}(1.0, MathOptInterface.VariableIndex(1))), 
        MathOptInterface.VectorAffineTerm{Float64}(2, MathOptInterface.ScalarAffineTerm{Float64}(2.0, MathOptInterface.VariableIndex(2))), 
        MathOptInterface.VectorAffineTerm{Float64}(2, MathOptInterface.ScalarAffineTerm{Float64}(1.0, MathOptInterface.VariableIndex(3)))
    ], [0.0, 0.0])

Now we want to multiply the ScalarAffineTerms with \(-1\) if the index is 2. I wrote some little uninteresting code snippet about it therefore, I'll not show it here. Just wanted to mention that this might be a gotcha one can run into.

Next up

In this next up section I often try to give a view into a possible future. That future often is part of parallel universe where I do exactly what I have planned so maybe another you can read it. 😉

You might want to join the JuliaLang Zulip and especially my ConstraintSolver stream to follow the list of priorities.

1. Have a strictly less than constraint.

   • This will most likely work with a tolerance and then use a MIP solver to get a better rhs

2. Ability to anti prune in a reified constraint which kinda needs 1.

   • Issue #222

3. Combining constraints with && and || inside indicator and reified constraints on both sides

   • Issue #205

4. Continue with 1D const element constraints

   • PR #213

   • This needs some form of 3. when used inside a reified or indicator constraint

Thanks everyone for reading and see you in the next couple of weeks!

If you enjoyed this post please consider supporting me and help me reach my goal of 20 patrons in 2021. 🎉

Special thanks to my 10 patrons!

Special special thanks to my >4$ patrons. The ones I thought couldn't be found 😄

• Site Wang

• Gurvesh Sanghera

• Szymon Bęczkowski

For a donation of a single dollar per month you get early access to these posts. Your support will increase the time I can spend on working on this blog.

There is also a special tier if you want to get some help for your own project. You can checkout my mentoring post if you're interested in that and feel free to write me an E-mail if you have questions: o.kroeger <at> opensourc.es

I'll keep you updated on Twitter OpenSourcES.

Want to be updated? Consider subscribing on Patreon for free