Computing bounds

This is an ongoing series about: How to build a constraint solver?
It’s part 18 of the series so we have came a long way. In the last weeks you might have seen that I implemented a lot of stuff in between.
Since a few days ago I have a bigger project now: Visualizing the data we have about the new coronavirus COVID-19. I’ve made a blog post about it Visualizing COVID-19 but there are a lot of things to do to make it really useful. It was just a fast way of trying stuff out. I want to make it interactive and scale by population/population density/area as well as getting an idea of the likeliness to get infected in each country by combining those data.

Nevertheless I want to explain a major change in the ConstraintSolver in this post.

The problem

Last time I implemented the option to optimize on linear objectives. The problem is that I had only a very simple way of computing bounds. That means if we have the following problem:

m = Model(CS.Optimizer) 
@variable(m, 0 <= x[1:10] <= 15, Int)
@constraint(m, sum(x) <=  15)
@objective(m, Max, sum(x))

and run optimize!(m). It takes a crazy amount of time to proof that the optimal solution is 15. What I did so far was that the bound computation only works with the variables. That means we have a look at the maximum value of each variable and add them up to obtain a bound. Here it will be 150. With the additional constraint sum(x) <= 15. We can immediately reduce that to 15 but well the ConstraintSolver isn’t that intelligent.

My first idea was to use some knapsack heuristic for this and started programming and realized way too late that it’s not working for negative coefficients or negative values.

Then I tried to only use it for constraints where it actually works like the one above and someone opened an issue issue #83. It’s quite a good problem to see the next issue my approach had.

I’ll show the problem in a minute. First I want to explain the approach I had in a bit more detail. Every constraint had a function which had as an input the objective function and then I computed with the knapsack heuristic a bound based on each constraint and objective separately and used the tightest bound.

The problem is of course that I used each constraint individually instead of combining them:

Code from the mentioned issue:

using JuMP
using ConstraintSolver
const CS = ConstraintSolver

model = Model()

# Variables
@variable(model, 0 <= inclusion[h = 1:3] <= 1, Int);
@variable(model, 0 <= allocations[h = 1:3, a = 1:3] <= 1, Int);
@variable(model, 0 <= days[h = 1:3, a = 1:3] <= 5, Int);

# Constraints
@constraint(model, must_include[h = 1:3],
    sum(allocations[h,a] for a in 1:3) <= inclusion[h]
);
# at least n
@constraint(model, min_hospitals,
    sum(inclusion[h] for h in 1:3) >= 3
);
# every h must be...