Why I'm betting on heuristic solvers for combinatorial optimization

In combinatorial optimization, the word “solver” almost automatically evokes the world of mixed-integer linear programming (MILP). It’s a Pavlovian reflex most of us share, and I blame it on our common academic background: the first solver many of us ever saw was solving a MILP, most likely written with an algebraic modeling library.

Over the past few years I’ve noticed, at least in my network, growing attention towards heuristic-based solvers, e.g. Timefold, Hexaly, Seeker. Not that MILP solvers don’t explore the solution space directly (via the so-called primal heuristics), but they do so in support of an exact algorithm. In heuristic solvers, by contrast, search in the primal space is the main engine, often built on much simpler paradigms such as local search.

I’m writing this article to try to organize what I’ve learned about the differences in design and mindset between traditional solvers and heuristic solvers, and to explain why I’ve become a big fan of the latter.

There are two main reasons for my optimism. The first concerns the present: these solvers considerably widen the range of problems that can be attacked with a general-purpose tool. The solving engine asks almost nothing of the model: it doesn’t need linearity or any specific algebraic structure, it only needs to be able to evaluate the effect of a change to the solution. In practice, most problems can be written down “as they are”, and the modeler’s time goes into the business problem rather than into contortions to please the solver. The second concerns the future: I expect that within a few years a growing share of problems that are today the exclusive domain of hand-written algorithms will become general-purpose solver territory.

MILP solvers are an extraordinary technology. They have reached a remarkable level of sophistication thanks to decades of research and development, the modeling libraries are extremely mature, and any operations research specialist can use them effectively to model and solve the most diverse problems. They are the state of the art for network flow, unit commitment, supply chain network design, and so on: in general, every class of problems where the LP relaxation is strong.

The modeling approach of a MILP solver is algebraic in nature: everything that depends on a decision becomes a variable, and the mutual consistency of the variables is enforced through constraints. Objective function and constraints must be linear expressions. Several common families of nonlinearity (logical predicates, piecewise linear functions) can be modeled, but model performance often suffers significantly, typically because the LP relaxation gets much weaker. In practice, it’s very easy to ruin an efficient model by adding a single class of ugly constraints. Other families of nonlinearity, e.g. generic continuous or black-box functions, simply cannot be modeled at all.

Consider, as an example, a classic CVRP. In a MILP formulation, the combinatorial structure of the problem, namely the fact that a solution is a set of sequences, is scattered across a constellation of binaries that the solver has to painstakingly piece back together. To be fair, the state of the art, using advanced branch-and-price-and-cut techniques, can solve instances of a few hundred customers to proven optimality, which is no small feat. My point is a different one: the moment the client’s problem comes with its umpteen real-world constraints, that road requires writing an ad hoc algorithm, which is exactly what in a previous post I called the last resort, because of its high development and maintenance costs.1

At this point, it’s reasonable to ask whether the combinatorial structures that show up most often could be declared rather than reconstructed. This approach has long been practiced in the world of constraint programming (CP). A global constraint such as AllDifferent, Cumulative, NoOverlap or Circuit is an abstraction that carries with it an algorithm specialized in handling the combinatorial structure of that constraint. In the AllDifferent example, the modeler declares that those variables must all be different, and the solver activates a specialized propagator that prunes the search tree far more efficiently than a linear reformulation would. In general, the more of the problem’s combinatorial structure is captured by global constraints, the more effective propagation is, and with it the search.2

There is, however, a second line of research in constraint programming that, in my opinion, doesn’t get all the credit it deserves. Around the turn of the millennium, Pascal Van Hentenryck and Laurent Michel asked themselves a question that only sounds obvious in hindsight: if constraint programming made tree search declarative, why not do the same for local search? Out of that came first Localizer, then the Comet language, and finally the book “Constraint-Based Local Search”, where the whole line of work is systematized.3

The two central ideas of Constraint-Based Local Search (CBLS) are exactly what we find today under the hood of heuristic solvers. The first is invariants, one-way constraints of the form “this quantity is by definition the sum of those other ones”, which the system commits to keeping up to date at every change to the solution. The second is differentiable objects: constraints and objective functions capable of telling how much a move would cost before the move is actually performed.4 Put together, they make it possible to evaluate hundreds of thousands of moves per second instead of recomputing the solution’s evaluation from scratch every time.

Let’s turn, then, to the architecture of modern heuristic solvers. The main ingredients I’ve identified are the following: high-level modeling abstractions beyond the classic integer and continuous variables, search moves specialized for those abstractions, handling of some classes of constraints through penalties in the objective function, and an incremental recomputation engine for the state of the solution.

High-level decision abstractions are a relatively recent addition, and in my opinion the most interesting development on the heuristic solver front. For instance, both Timefold and Hexaly support variable-length lists as decision variables. This abstraction is ideal for efficiently modeling problems whose decisions are permutations. In the CVRP example, we can model the visits performed by a vehicle with a single list, and the capacity constraints become a sum computed by walking the list with a for loop. The solver knows that object is a sequence, and can therefore apply operators specialized for that data structure, e.g. reversing a subsequence or moving a segment from one list to another, rather than just flipping a binary (and hoping something sensible comes out). Similar things can be said about unordered sets and intervals.5 The result isn’t just expressive elegance: by leveraging these abstractions, heuristic solvers are reporting numbers that until recently were out of the question for a general-purpose solver.6

What follows is a modeling approach I’d call computation-oriented, rather than algebraic. You declare a limited set of decision variables, and every quantity that depends on them is computed through procedures familiar to any programmer: loops, conditionals, array indexing. The result is not only a much smaller set of decision variables (and therefore a more efficient search), but also that the modeler is forced to think about what the real degrees of freedom of the business problem are, which is an extremely healthy exercise. In my experience, the question “what are the decisions?” leads to far more useful conversations with domain experts than the question “what are the variables?”.

Search moves specialized for certain classes of data structures and constraints are, in my opinion, the analogue of the global constraints of constraint programming. It’s the same principle, put to a different use. In constraint programming, encoding the combinatorial structure serves to prune the search tree efficiently, whereas in local search it serves to move from one solution to another without breaking feasibility.

Finally, the incremental computation engine and the handling of constraints through penalties are a direct inheritance from Constraint-Based Local Search. In particular, for those classes of constraints that make it hard to move from one solution to another without breaking feasibility, it’s far more efficient to formulate a high-priority cost function and let the solver use that gradient to guide the search. However, for this approach to be effective, part of the problem’s combinatorial structure has to be encoded through the appropriate structures and global constraints. This sentence from Hexaly’s documentation sums up the idea nicely:

Turning too many constraints into objectives will lead to a model without any structure, whereas one of the strengths of Hexaly Optimizer is to perform moves preserving the feasibility and thus “making sense”.

I’d also like to say a few words about scalability with respect to the available hardware. This is a topic on which I changed my mind several times while writing this article, so I hope what I’m about to say comes out sensibly (also because whatever I write on the subject, I feel like I’m doing someone an injustice). Let’s start from the fact that MILP solvers are notoriously hard to parallelize, and in my experience, on generic workloads, one observes speedups close to the logarithm of the number of cores.7 But local search doesn’t scale arbitrarily either: independent portfolio-style searches reduce the variance of the outcome more than the time to solution, and parallel move evaluation runs into the fact that the search itself remains sequential. The big difference, in my opinion, lies in where the bottleneck ends up. In heuristic solvers it’s the throughput of incremental evaluations, which is a matter of data locality, memory layout, and so on, whereas in branch-and-cut it’s the logical structure of the search, where hardware helps less. In other words, my very high-level impression is that there is a lot more room for engineering and performance improvement in heuristic solvers.

At this point, after singing the praises of heuristic solvers for twenty-odd paragraphs, let’s turn to the dual side (heh heh) of the matter.

One limitation of heuristic solvers is that they generally provide no dual bound, and cannot prove that a problem is infeasible. The optimality gap is useful as a well-defined termination criterion, certainly more understandable than an arbitrary time limit. I find it much less useful as a certificate of solution quality, because it certifies optimality with respect to a model of reality, while in industrial problems modeling approximations and data errors dominate by several orders of magnitude. In particular, a piece of data can be stale, approximate, estimated, or simply dirty, and against imperfect data a 0% gap doesn’t bring much value. Moreover, the baseline we compare against is almost never the optimum, but the plan that is currently made by hand or by a system of simple static rules, and the delta against that is the metric the client is actually buying. Proving infeasibility, on the other hand, is far more useful, because more often than not it points to a problem in the input data that would otherwise have been hard to spot.

Another aspect not to be underestimated is the mental shift required to use this kind of solver effectively. MILP solvers have enormous market penetration, and any operations research specialist can get their hands on them. None of this holds (yet) for heuristic solvers. The modeler has to learn to think in terms of decisions and derived quantities rather than variables and constraints, has to decide which constraints are structure (and should be declared) and which are penalties (and should go into the objective), and needs a high-level idea of which moves the solver is able to make. The literature is thinner and models are not easily portable from one solver to another. It’s not an insurmountable barrier, but it is a barrier, and for some companies it can be a significant obstacle to adoption.

Finally, there are classes of problems for which local search isn’t necessarily the right tool. The most obvious case I’ve already mentioned at the beginning: problems where the LP relaxation is strong, where a MILP solver reaches proven optimality in reasonable time. Less obvious is the case of heavily constrained problems where finding any feasible solution at all is the hard part, and the penalty-in-the-objective approach doesn’t work well. I’d also add mostly continuous problems, where simplex and interior point methods are hard to beat.

Some of these limitations are being mitigated by an increasingly deep hybridization of exact and heuristic methods. Hexaly, for instance, activates specialized solution strategies (MILP, interior point, and so on) when it recognizes certain structures in the model, and is therefore able to provide dual bounds or prove infeasibility. It’s interesting to note how this hybridization is the mirror image of what happens in MILP solvers. In a MILP solver the exact algorithm is the main engine and heuristics work in its service. In Hexaly it’s the other way around: local search is the main engine, and exact algorithms are only invoked when the solver recognizes a structure in the model they are good at. I wouldn’t be surprised to see this kind of hybridization applied to other heuristic solvers in the near future.

To conclude, let me try to sum things up. The ideas behind heuristic solvers have been around for decades and rest on a solid literature, but it’s only in recent years that the engineering of the products has truly made a leap in quality. The result is that classes of problems that until recently were out of reach for a generic solver are now tractable, with more than satisfactory results. Every signal I see suggests this trajectory will continue. Accordingly, I expect the job of the operations researcher in industry to shift more and more from designing specialized algorithms to deeply understanding business problems. As far as I’m concerned, that’s great news! It means I’ll get to spend more time reasoning about the client’s problem and iterating faster. Indeed, most of the value of an optimization project is created in the ongoing dialogue with domain experts, in the continuous refinement of the decision process: every hour the solver saves me is an hour I can spend understanding the domain better.

Let me know about your experiences with solvers, whether you agree with my impressions, or whether I’ve left something out.

Footnotes

  1. Trust me: no matter how complete you think your optimization engine’s constraint portfolio is, you will run into the case where the client has one or two quirks you hadn’t thought of, and they will ruin your day. ↩

  2. I’ve always found it curious that this approach is the exact opposite of MILP modeling best practices. In particular, dense constraints are generally poison for a MILP solver. ↩

  3. Pascal Van Hentenryck and Laurent Michel, Constraint-Based Local Search, MIT Press, 2005. My dear friend and former boss Jean-François Baffier introduced me to this book during my internship at IIJ Lab. I’m also told that a paper copy is rather hard to come by these days. ↩

  4. Modern heuristic solvers don’t ask the user to evaluate the deltas of local moves by hand; they infer them automatically from the structure of the problem. For instance, Timefold uses an engine called Bavet to process “Constraint Streams” through a declarative, functional-style API. The engine represents the relationships between variables and constraints as a network of nodes; when local search proposes a change, the engine automatically propagates the change through the network of operators. I advice taking a good look at the source code, it’s really instructive. ↩

  5. I’m very curious to see whether new modeling abstractions will be introduced. I wouldn’t be surprised to see maps appear as a decision construct in the near future, which would make one-to-one and many-to-one assignment problems natural to model. ↩

  6. For instance, Hexaly publishes benchmarks on public CVRP instances with three thousand to thirty thousand customers, with average gaps to the best known solutions below 3% within one minute. The result is all the more striking when you consider it comes from a generic solver rather than a specialized algorithm. ↩

  7. On the subject, I recommend Gurobi’s talk How to Exploit Parallelism in Linear Programming and Mixed-Integer Programming. ↩

← All posts