lp-mip-solver

How it is built

Architecture overview

An LP/MIP solver takes a model min cᵀx subject to L ≤ Ax ≤ U, l ≤ x ≤ u, with some variables required to be integer, and returns a provably optimal solution (or a proof that none exists). The model is stored as a sparse matrix, shrunk by presolve, and solved as a sequence of linear-programming relaxations. For a pure LP one relaxation is the answer, and the LP engine is chosen by problem size: dual simplex by default, an interior-point method or a first-order method for very large models. For a MIP the relaxation gives a bound; branch and cut tightens that bound with cutting planes, splits the problem on fractional integer variables, and discards every subtree whose bound cannot beat the best solution found so far. The one idea that makes it work: a MIP is a tree search in which every node is an LP, so the simplex method's ability to re-solve a slightly changed LP in a few pivots is the engine, and presolve, cuts and heuristics exist to make that tree small.

Parts:

• Presolve and postsolve — shrink and tighten the model, keep an undo stack — Presolve reductions
• Branch and cut with a warm-started LP engine — the tree, the bound, the cuts, one LP per node — Mixed-integer programming (branch and bound), Cutting planes
• Interior-point engine — fast from-scratch solve of large LPs — Interior-point method (primal-dual, Mehrotra), Sparse direct solver, Fill-reducing ordering
• First-order engine — matrix-vector-product-only LP solves, GPU friendly — primal-dual-hybrid-gradient (no page yet)
• Parallel tree search — independent nodes as tasks — Work-stealing deque

 model.mps / .lp ──► reader ──► sparse model ──► PRESOLVE ──► reduced model
                                                      │             │
                              postsolve stack ◄───────┘             ▼
                                                        ┌── root LP (simplex | IPM+crossover | PDLP)
 original-space solution ◄── postsolve ◄── incumbent ◄──┤        │ fractional x*
                                                        │        ▼
                                                        │   cut separation ──► cut pool ──► LP + cuts
                                                        │        │
                                                        │        ▼
                                                        └── B&B tree: node queue ──► dual simplex warm start
                                                                 ▲   │  branch on fractional var, prune by bound
                                                    heuristics ──┘   └──► new incumbent ⇒ tighter pruning

How the parts fit together

Model and storage

Everything downstream reads the constraint matrix, so it is held in compressed sparse form (Sparse matrices (CSR)): row access for presolve, activity bounds and cut generation, column access for pricing and the simplex ratio test, usually with both orientations kept. Industrial models are extremely sparse, so work is counted in nonzeros, not in rows × columns. HiGHS reads models from MPS or (CPLEX) LP files, and also lets callers build and modify a model through library calls (HiGHS docs, further features). MPS is a fixed, column-oriented punch-card-era text format with ROWS, COLUMNS, RHS and BOUNDS sections (lp_solve MPS description); the format itself has no standard direction of optimisation, which is a classic source of sign bugs when files move between solvers.

Presolve and postsolve

Presolve reductions run before any LP is solved and again at the MIP root: fixed variables are removed, singleton rows become bounds, redundant rows are dropped, bounds and big-M coefficients are tightened. Each reduction pushes an undo record so that postsolve maps the solution (and, for simplex, the basis) back to the original model. HiGHS presolves by default because the time saved is normally far greater than the cost; the exception is an LP solved by simplex from a supplied valid basis (HiGHS docs). Activity-based bound tightening needs rigorous bound arithmetic (interval-arithmetic, no page yet); finding duplicate or parallel rows is a hashing job (Fast hash functions). Contract: the reduced model has the same optimal value, and postsolve is cheap and deterministic. Gurobi developers describe a commercial presolver in Achterberg et al. (Sources).

Simplex engine

The Simplex method, in its dual form, is the default LP solver in HiGHS: faster and more robust in practice than the primal (HiGHS docs, solvers). It is the MIP's engine because a branch changes one variable bound and a cut adds one row, so the previous optimal basis is still dual feasible and the dual simplex reaches the new optimum in a handful of pivots (a hot start). An interior-point method would start over; HiGHS states that only simplex can be hot-started.

Each pivot needs the basis inverse only through solves, so the solver keeps a sparse LU factorisation of the basis (LU and Cholesky factorization) and updates it after each basis change (product-form or Forrest-Tomlin style, Forrest and Tomlin 1972) instead of refactoring; it refactors after a bounded number of updates (HiGHS: simplex_update_limit). Pricing uses edge weights (simplex_dual_edge_weight_strategy). Floating-point error bites here: ratio tests need tolerances, scaling matters, and a drifting factorisation shows up as an infeasible basis, so values are recomputed after refactoring, with cost perturbation as the usual defence against degeneracy. Parallelising one dual simplex is hard; Huangfu and Hall (2018) do it, but the current docs say later serial improvements were not carried over, the speed-up is typically bounded by memory channels, and independent instances per worker are the better use of cores (HiGHS docs, parallel).

Interior-point engine

For large LPs the Interior-point method (primal-dual, Mehrotra) needs, as a rule of thumb, a few tens of iterations (order of), each solving a Newton system (Newton root finding) on the perturbed optimality conditions: one sparse factorisation of the normal equations or augmented system (Sparse direct solver). Fill-in decides the cost, so the matrix is reordered first (Fill-reducing ordering); HiGHS' direct IPM, HiPO, offers nested dissection (METIS, a Graph partitioning application), approximate minimum degree and reverse Cuthill-McKee, is parallel, and per its docs can beat the other HiGHS IPM by more than an order of magnitude as problems grow (HiGHS docs, solvers). That other IPM, IPX, solves the Newton systems iteratively with preconditioned Conjugate gradient (Schork and Gondzio 2020): memory-light, but the iteration count depends on the preconditioner.

The seam to know: an interior point ends near the middle of the optimal face, not at a vertex. A MIP needs a basis to warm-start simplex, so IPM output goes through crossover; HiGHS runs it by default and can skip it or run it only when the IPM did not certify optimality. HiGHS also lets the user choose the MIP root LP solver, and uses an IPM to compute the analytic centre for its centralised-rounding heuristic.

First-order engine

primal-dual-hybrid-gradient (no page yet) is the PDLP method: restarted PDHG with preconditioning, using only sparse matrix-vector products, so it needs little memory beyond the model and maps onto a GPU (Applegate et al. 2021; GPU compute shaders covers the launch model). HiGHS ships cuPDLP-C and a native HiPDLP for NVIDIA (and, for HiPDLP, AMD) GPUs, and its docs are candid: on a CPU they are unlikely to compete, there is no crossover to a basic solution, and first-order LP solvers are "very much in their infancy, and are not robust" (HiGHS docs, parallel). Use it for models too big to factorise when moderate accuracy is enough: PDHG converges slowly to tight tolerances.

Branch and cut

Mixed-integer programming (branch and bound) is the controlling loop: pick an open node, solve its LP, prune if it is infeasible or its bound cannot beat the incumbent, otherwise branch. The decisions that matter are the branching rule (reliability branching combines pseudocosts with strong branching on unreliable candidates; Achterberg, Koch and Martin 2005), node selection, and how many cuts to add. The open-node queue is a priority queue keyed by bound or estimate (d-ary heap).

Cutting planes are separated at the root and selectively in the tree: a valid inequality that removes the fractional LP optimum without removing an integer-feasible point is appended as a row (hence the hot start above). A cut pool with ageing drops inactive cuts, and generators filter badly scaled cuts, since huge coefficient ranges can cut off a true solution (Floating-point error). HiGHS exposes pool ageing, cut separation at nodes, symmetry detection, clique tables and heuristics including RENS, RINS, feasibility jump, shifting and ZI-round (options). Primal heuristics (diving, the feasibility pump of Fischetti, Glover and Lodi 2005, RENS of Berthold 2014) are not covered by a block page here; a good incumbent early lets every pruning test fire sooner.

Parallelism and reproducibility

Open nodes are independent tasks, so a work-stealing scheduler (Work-stealing deque) keeps workers busy, and shared structures (incumbent, cut pool, clique table) need disciplined publication (Memory ordering). The price is determinism: thread timing changes the order incumbents are found, hence which nodes are pruned. HiGHS' development docs are not consistent here (the solvers page calls its MIP solver largely single-threaded with a multithreaded prototype; the parallel page describes multi-threaded tree exploration when parallel is not off) and are marked work in progress. SCIP's suite has a separate parallel framework, UG (SCIP 9.0 report). Randomised tie-breaking is seeded (PRNGs (PCG, xoshiro)) so runs repeat and benchmark claims can be tested across seeds.

Data and control flow

One MIP solve.

1. Read the MPS or LP file into a sparse model (linear in file size).
2. Presolve to a fixpoint, pushing undo records; collect clique and implication data for cuts and propagation.
3. Solve the root LP (dual simplex; or IPM plus crossover; or PDLP for huge models on a GPU). Usually the largest single LP of the run.
4. Root loop: separate cuts, add the best as rows, re-solve with the dual simplex, run root heuristics, stop when the bound stalls.
5. Branch: push two children; each node warm-starts from its parent's basis, usually a few pivots.
6. On each improved incumbent update the cutoff, prune the queue, optionally start a neighbourhood-search heuristic.
7. Stop at the gap tolerance (mip_rel_gap, mip_abs_gap in HiGHS) or a node, time or leaf limit.
8. Postsolve the incumbent to the original variables; report status and bound.

One large LP. Presolve, interior point (each iteration one sparse factorisation or CG solve), crossover if a basic solution is needed, postsolve. PDLP skips crossover and spends its time in matrix-vector products.

One hot-started re-solve. A caller edits a bound or coefficient, the stored basis is repaired (surplus basic variables dropped or slacks added), and simplex continues; this is the pattern behind column generation and rolling-horizon planning.

Scaling

1. Model size: presolve first (cheapest, pays on generated models); then the engine that fits the shape: dual simplex for MIP and moderate LPs, IPM for large LPs, PDLP on a GPU beyond that.
2. Factorisation memory: nested dissection or minimum degree, multifrontal factorisation for large structured models, or an iterative IPM to trade time for memory.
3. Tree size: cuts, branching, symmetry handling and early incumbents shrink the tree by far more than a faster node; a better bound beats a faster LP.
4. Wall-clock: tree search is memory-bound (HiGHS: unlikely to help beyond about eight threads); for many instances run one per worker.
5. Beyond one machine: decomposition (Benders, column generation) or a distributed framework such as UG changes the formulation or scheduler, not the engine.

Failure modes

Real implementations

The open solvers' relative speed is tracked in Hans Mittelmann's public benchmarks (LP, MILP) at plato.asu.edu/bench.html; the standard test set is MIPLIB 2017. Neither is reproduced here, since the tables change with every release.

Sources

• Q. Huangfu, J. A. J. Hall, Parallelizing the dual revised simplex method, Math. Prog. Comp. 10, 2018. paper
• L. Schork, J. Gondzio, Implementation of an interior point method with basis preconditioning, Math. Prog. Comp. 12, 2020. paper
• D. Applegate et al., Practical Large-Scale Linear Programming using Primal-Dual Hybrid Gradient, NeurIPS 2021. paper
• J. J. H. Forrest, J. A. Tomlin, Updated triangular factors of the basis to maintain sparsity in the product form simplex method, Math. Programming 2, 1972. paper
• S. Mehrotra, On the implementation of a primal-dual interior point method, SIAM J. Optim. 2, 1992. paper
• R. E. Bixby, Solving real-world linear programs: a decade and more of progress, Operations Research 50, 2002. paper
• T. Achterberg, T. Koch, A. Martin, Branching rules revisited, Oper. Res. Lett. 33, 2005. paper
• T. Achterberg, R. E. Bixby, Z. Gu, E. Rothberg, D. Weninger, Presolve reductions in mixed integer programming, INFORMS J. Computing 32, 2020. paper
• T. Achterberg, SCIP: solving constraint integer programs, Math. Prog. Comp. 1, 2009. paper
• M. Fischetti, F. Glover, A. Lodi, The feasibility pump, Math. Programming 104, 2005. paper
• T. Berthold, RENS, Math. Prog. Comp. 6, 2014. paper
• A. Gleixner et al., MIPLIB 2017: data-driven compilation of the 6th mixed-integer programming library, Math. Prog. Comp. 13, 2021. paper
• S. Bolusani et al., The SCIP Optimization Suite 9.0, 2024. report
• HiGHS documentation: solvers, parallel, further features, options
• lp_solve project, MPS file format. page