Creation date: 2020-06-17
A lot of things happened on this blog lately as I've changed the layout completely. If you haven't seen the blog about the blog :D Check out Using Franklin.jl
Now that I don't have two layouts one for Patrons and one for others I can use all the features which only exist in the new layout for everyone and spend more time implementing new ones ;)
Please comment here or reach me otherwise if you have improvements for the new Blog.
Okay enough housekeeping! This blog is again about the constraint solver.
If you haven't read about it before please check out the first post.
I'm basically trying to build a full constraint solver from scratch in Julia, which isn't an easy task. In this series you can read about my journey.
This post is about the indicator constraint again. You might want to check out Indicator Part 1 before continue reading.
This post is about extending JuMP.jl to support constraints of the form:
@constraint(m, b => {[x, y, z] in CS.AllDifferentSet()})
whereas before only
@constraint(m, b => { x + y + z <= 10)})
so affine inner constraints were allowed.
Let's just assume we are a user and want to use that constraint:
using JuMP
using ConstraintSolver
const CS = ConstraintSolver
m = Model(CS.Optimizer)
@variable(m, 0 <= x <= 1, Int)
@variable(m, 0 <= y <= 1, Int)
@variable(m, a, Bin)
@constraint(m, a => {[x,y] in CS.AllDifferentSet()})
@objective(m, Max, a)
optimize!(m)
I got this error:
julia> indicator()
ERROR: MethodError: no method matching
_build_indicator_constraint(::JuMP.var"#_error#68"{Symbol},
::VariableRef,
::VectorConstraint{VariableRef,ConstraintSolver.AllDifferentSetInternal,VectorShape},
::Type{MathOptInterface.IndicatorSet{MathOptInterface.ACTIVATE_ON_ONE,S}
where S<:MathOptInterface.AbstractScalarSet})
Closest candidates are:
_build_indicator_constraint(::Function, ::AbstractVariableRef, ::ScalarConstraint, ::Type{MathOptInterface.IndicatorSet{A,S} where S<:MathOptInterface.AbstractScalarSet}) where A at /home/ole/.julia/dev/JuMP/src/indicator.jl:13
For a constraint which JuMP supports but the solver doesn't the output looks like:
julia> second_order_cone()
ERROR: Constraints of type MathOptInterface.VectorOfVariables-in-MathOptInterface.SecondOrderCone are not supported by the solver and there are no bridges that can reformulate it into supported constraints.
Then one would need to write a MOI.supports_constraint
function for that type and actually implement it.
In our case we need to work on the parsing first and change only little on the solver side as we mostly implement the indicator constraint already. In our case there is no real difference of what the inner constraint looks like.
The error message shows us that we need to create a _build_indicator_constraint
function and also points us to where that function exists in JuMP for ScalarConstraint
.
function _build_indicator_constraint(
_error::Function, variable::AbstractVariableRef,
constraint::ScalarConstraint, ::Type{MOI.IndicatorSet{A}}) where A
set = MOI.IndicatorSet{A}(moi_set(constraint))
return VectorConstraint([variable, jump_function(constraint)], set)
end
We need to have
"""
Support for indicator constraints with a set constraint as the right hand side
"""
function JuMP._build_indicator_constraint(
_error::Function, variable::JuMP.AbstractVariableRef,
constraint::JuMP.VectorConstraint, ::Type{MOI.IndicatorSet{A}}) where A
set = CS.IndicatorSet{A}(variable, MOI.VectorOfVariables(constraint.func), constraint.set, 1+length(constraint.func))
vov = VariableRef[variable]
append!(vov, constraint.func)
return JuMP.VectorConstraint(vov, set)
end
So we need constraint::VectorConstraint
here to not only support the all different constraint but all kind of in Set
constraints.
We can then extract the variable indices of the inner constraint and create a new set. The new set is needed because of:
struct IndicatorSet{A, S <: AbstractScalarSet} <: AbstractVectorSet
set::S
IndicatorSet{A}(set::S) where {A, S <: AbstractScalarSet} = new{A,S}(set)
end
Here S
only allows AbstractScalarSet
. The set is defined in MOI and you can jump to the code directly when clicking on the link below the code ;) Oh yeah a new feature...
I've just added my own set for this:
struct IndicatorSet{A} <: MOI.AbstractVectorSet
ind_var :: JuMP.VariableRef
func :: MOI.VectorOfVariables
set :: MOI.AbstractVectorSet
dimension::Int
end
I'm not only saving func
and set
but also ind_var
which is the activation variable. Currently I'm not using the variable so I could remove that. The func
and set
are used for dispatching later and set
should be called set
because it is the name in the original MOI.IndicatorSet
. The dimension
which we set to 1+length(constraint.func)
is important as JuMP and MOI not work as one might naively think :D I'm still quite sure that I don't know how it works exactly but this is a very rough outline of what I think. Please correct me if I'm wrong in the comments.
One might think (at least I did)
Writing
m = Model(CS.Optimizer)
@variable(m, 0 <= x <= 2, Int)
calls add_variable
and add_constraint
for lower and upper bound or Interval
such that there is a one to one correspondents to my MOI.add_variable
and MOI.add_constraint
code but this is not how it works.
There are several reasons why this is not the case. One obvious one is that JuMP allows bridges such that when interval is not supported (which actually doesn't get called even if it is supported in this case) it will be reformulated to 0 <= x
and x <= 2
. Some more of these bridges exist and some are more complicated than this trivial example i.e I don't have to support ACTIVATE_ON_ZERO
in the indicator constraint which the user would write like
@constraint(m, !b => { ... })
Because it can be reformulated by adding an extra variable i.e. y
such that y is the opposite of b
(an extra constraint) and then y => { ... }
.
I think what JuMP and MOI roughly do is creating an internal model which gets all variables and constraints and when calling optimize!
it gets copied over to the actual optimizer. In this step the order is actually different than the order one would normally think of:
variables
constraints
just like the user modeled it...
The so called VectorOfVariables
constraints like a in Set()
are added first and those variables are added first then the constraint and then the rest of the variables. I'm not too sure why this is the case. I think it has something to do with the fact that some solvers which are not written in Julia support it only in this way. Again please correct me if I'm talking bullshit. ;)
Nevertheless I had the wrong dimension first (somehow I've set it to 0
) probably because I didn't know why I need it :/ and that was the reason why variables weren't added to the model because MOI "thought": "Oh this constraint doesn't have variables so I don't need to add them".
After adding JuMP._build_indicator_constraint
we now have to support that constraint:
function MOI.supports_constraint(
::Optimizer,
func::Type{MOI.VectorOfVariables},
set::Type{IS},
) where {A, IS<:CS.IndicatorSet{A}}
return A == MOI.ACTIVATE_ON_ONE || A == MOI.ACTIVATE_ON_ZERO
end
So here I used the same structure as in the other supports_constraint
but use CS.IndicatorSet
here.
Additionally of course:
function MOI.add_constraint(
model::Optimizer,
vars::MOI.VectorOfVariables,
set::IS,
) where {A, T<:Real,IS<:CS.IndicatorSet{A}}
In it I have a normal:
internals = ConstraintInternals(
length(com.constraints) + 1,
vars,
set,
Int[v.value for v in vars.variables]
)
as well as for the inner constraint:
inner_internals = ConstraintInternals(
0,
MOI.VectorOfVariables(vars.variables[2:end]),
set.set,
Int[v.value for v in vars.variables[2:end]]
)
inner_constraint = init_constraint_struct(typeof(set.set), inner_internals)
Here the 2:end
is that we don't have the first variable which is the indicator variable be part of the inner constraint. The last line creates the correct struct depending on what inner constraint we have.
i.e I've defined:
function init_constraint_struct(::Type{AllDifferentSetInternal}, internals)
AllDifferentConstraint(
internals,
Int[], # pval_mapping will be filled later
Int[], # vertex_mapping => later
Int[], # vertex_mapping_bw => later
Int[], # di_ei => later
Int[], # di_ej => later
MatchingInit(),
Int[]
)
end
The last thing we have to do is to update the prune_constraint!
and other _constraint
functions such that they support CS.IndicatorSet
and MOI.IndicatorSet
.
That looks like this:
function prune_constraint!(
com::CS.CoM,
constraint::IndicatorConstraint,
- fct::VAF{T},
+ fct::Union{MOI.VectorOfVariables, VAF{T}},
set::IS;
logs = true,
-) where {A, T<:Real, ASS<:MOI.AbstractScalarSet, IS<:MOI.IndicatorSet{A, ASS}}
+) where {A, T<:Real, ASS<:MOI.AbstractScalarSet, IS<:Union{IndicatorSet{A}, MOI.IndicatorSet{A, ASS}}}
The rest stayed basically the same. This is available in v0.2.0 as I had some bugs in v0.1.8 :D
Hope you found this helpful and also enjoyed the new layout and functionality like diffs and links to the GitHub code snippets.
Thanks for reading and special thanks to my six patrons!
List of patrons paying more than 4$ per month:
Site Wang
Gurvesh Sanghera
Szymon Bęczkowski
Currently I get 19$ per Month using Patreon and PayPal.
For a donation of a single dollar per month you get early access to the posts. Try it out at the start of a month and if you don't enjoy it just cancel your subscription before pay day (end of month).
I'll keep you updated on Twitter OpenSourcES as well as my more personal one: Twitter Wikunia_de