The Programming Paradigm to Find One Solution Among 8,080,104 Candidates
mlabonne.github.io
mlabonne.github.io
:- use_module(library(clpz)).
solution(N) :-
N in 0..10_000,
N #= 13*_,
N #= 19*_,
N #= 37*_.
Note how seamlessly constraint programming blends into Prolog: We simply state these constraints as relations that must hold for logic variables, just as always in Prolog! In fact, plain Prolog can be regarded as a particular instance of constraint logic programming (CLP), namely as constraint logic programming over Herbrand terms, abbreviated as CLP(H).Sample query and answer:
?- solution(N), indomain(N).
N = 0
; N = 9139.
We see from this answer that there are in fact two solutions to the puzzle: 0 and 9139. Prolog reports all solutions on backtracking.A nice feature of the Scryer Prolog toplevel is that we can simply press "a" to enumerate all solutions. For example, if we increase the upper bound to 100_000 as mentioned in the course of the article, we get:
?- N in 0..100_000,
N #= 13*_, N #= 19*_, N #= 37*_,
indomain(N).
N = 0
; N = 9139
; N = 18278
; N = 27417
; N = 36556
; N = 45695
; N = 54834
; N = 63973
; N = 73112
; N = 82251
; N = 91390.
So, 11 solutions (not 10).Hence, the sample code shown in the article should probably change:
model.NewIntVar(1, ...)
to army = model.NewIntVar(0, ...)
in order to correctly model the equation shown before the code.The knapsack problem shown in the article can also be easily expressed and solved:
?- Vs = [Bread,Meat,Beer],
Vs ins 0..sup,
Bread + 3*Meat + 7*Beer #=< 19,
Popularity #= 3*Bread + 10*Meat + 26*Beer,
labeling([max(Popularity)], Vs).
yielding solutions in decreasing order of popularity: Vs = [2,1,2], Bread = 2, Meat = 1, Beer = 2, Popularity = 68
; Vs = [5,0,2], Bread = 5, Meat = 0, Beer = 2, Popularity = 67
; Vs = [0,4,1], Bread = 0, Meat = 4, Beer = 1, Popularity = 66
; Vs = [1,1,2], Bread = 1, Meat = 1, Beer = 2, Popularity = 65
; ... .
This answer shows that the optimal solution with popularity 68 is unique in this case. The program output in the article shows "3712 popularity", yet the text also mentions 68: The number 3712 stems from using 1_000 as the capacity.It is notable how compactly and conveniently such tasks can be expressed and solved with Prolog.
The equation should be 1 ≤ army ≤ 10000, since you cannot have an army of zero soldiers.
I've been working with some young devs who only learned Python and web stuff (?!). One thing I've noticed about these devs is that they have some really bad intuitions about how long it takes to brute force problems like this. Brute forcing the hardest problem in this article requires less than 3 seconds for a C program but requires 3 minutes for a Python program.
Bumping the capacity to 10,000 (so from 8 million options to one trillion options) takes about 40 minutes in C, but would probably take all day in Python. Once you get much past 10,000 the effort of pulling in a library starts to make sense because, even though C brute force is still tractable, you need to start parallelizing things and that's more effort than an import statement.
Also, note to author: there are several mistakes in this article. There are bonds errors in your army examples and you used capacity = 19 instead of capacity = 1000 in one of the questions for the bread/meat/beer problem.
you may not necessarily be timing what you think you are timing
Well, yes, there are all sorts of tricks, including the use of a solver :)
In my earlier post, read "brute-force" as "full enumeration". The comment was about dev's intuitions re: the time required to fully enumerate a space; ie when is it worth thinking about anything other than the stupidest solution possible.
but
real world applied optimisation problems (not simple illustrative toy problems as in the article) often have an exponentially large number of candidate or feasible solutions, attempting to brute force the search space in C is not going to be a fruitful exercise. e.g. you might be trying to solve a facility location problem to place `m` depots on nodes in a graph with `n` vertices, where each node can have 0 or 1 depot, so the space of feasible solutions is `2^n`, the powerset of the set of nodes. for a smallish applied problem, n might be `n=5000` . there might be a bunch of set-cover-like constraints that require a bunch of demand nodes on the graph are close enough to a depot node to be "covered", and objective function costs for each depot constructed or the distance from each demand node to the closest depot.
commercial mixed integer solvers such as Gurobi have a whole bag of tricks to reduce the problem into a smaller problem during a presolve phase before the real MIP solver even gets woken up. a presolve phase employs a bunch of heuristics to detect common kinds of problem structure or substructure -- for example in the toy "optimisation and beer" problem it might immediately recognise that the only constraint is a knapsack constraint and palm the entire problem to a specialised knapsack solver which trivially solves the problem. more realistic applied problems are often weirder and messier -- more generally a presolve might be able to eliminate some but not all of the decision variables, or use substructure-specific techniques harvested from the literature to prove and inject additional constraints that make the real MIP solver's search space smaller.
My observation was intended primarily as a pedagogic critique. The C solutions to these problems don't require any cleverness and are transliterations of the Python equivalents.
I agree re: optimizers of course. The article would be better if the final example was something where even really cleverly optimized "is this really correct?" enumeration techniques under-performed the optimizer.
If you need 100 clock cycles per brute-force check, you can brute-force 40 million cases per second.
8-million is actually very small. Modern computers (GPUs really) are on the order of reasonably brute-forcing 2^40-ish search spaces these days (or one-in-a-quadrillion, depending on the complexity per check)
~1GHz clock on GPUs x 4096 shaders == 4-trillion gpu-shader core cycles per second. If you're willing to wait 1-hour (3600 seconds), that's 14-quadrillion gpu-shader clock cycles.
Part II only takes 73 µs for the small example, and 6.77s for the large one (the one in the title).
¹ah, I re-read the article and "The lower bound is 1 since we know there's an army" … but that's not ever established in the math parts. Do we know there's an army? Our scouts seem … eccentric. Of course, while the prose claims a lower bound, the math immediately below that paragraph contradicts it with "0 ≤ army ≤ 10000"… (given the waffling in the article on this bound, I suspect the author's initial implementation of a constraint solver probably returned 0, at which point they updates some parts of the article…)
a brute force triple loop would check a bunch of possibilities that aren't solutions either the very naive 1000×1000×1000 (for one billion) or slightly less naive 1000×333×142 (for 47.3 million).
Not that a billion would take that long either on a modern computer, but if you look at the OR-Tools solution, it's really not any longer or more complicated than loops would be. As long as you have good tools (coughprologcough) it doesn't take very long for constraint solving to come just as naturally as loops.
So the N=1k takes about 6.5ms on my desktop, and the N=10k takes a mere 400ms. Only when I try 100k do I need to actually wait, for an egregious 40s.
And, a more trifling note, their variable bounds are too generous...
meat = model.NewIntVar(0, 1 + capacity // 3, 'meat')
beer = model.NewIntVar(0, 1 + capacity // 7, 'beer')
reduces the triple-loop estimate to a mere 47M cases.1) Constraint Logic Programming — quite cool as long as you have experience in Prolog, but choice of underlying solvers is limited: Prolog runtime has built-in solver.
2) Answer Set Programming - syntactically similar to Prolog, but later compiled (grounded) to elementary logic solvable with SAT techniques — very efficient.
3) Constraint Programming via Solver API — there are many solvers: Choco, Oscar, OR-tools that are usable programmatically. This approach is advertised in the article.
4) common API - there are libraries, e.g. Numberjack, that generalize API over several solvers
5) Constraint Programming via Modelling Language: there exist languages designed just to define Constraint Programming models. MiniZinc is at the moment the most mature and has the most extensive support from industry.
6) related techniques: SAT, SMT
Normally I start with 5) as it allows me to quickly prototype a solution for the problem. Only when it's not enough, I switch to 3).
$ python3 -c 'import math; print(math.lcm(13, 19, 37))'
9139