A Practical Combinatorial Interior Point Method for Minimum-Cost Flow

1 Introduction

• The authors report that practical IPM implementations had not emerged, while recursive data structures, depth, and constants obstructed a faithful implementation of.They frame this as a gap between the asymptotics of the best known algorithms and those of algorithms that can be run.Source: §1 · pp. 2–3
• The paper gives the first practical implementation of the almost-linear-time minimum-cost flow algorithm of [CKL+25].It discretizes the implementation to address numerical precision and removes the need for oracle access to the optimal cost.Source: §1 · p. 3
• The implementation estimates progress more carefully and uses empirically motivated rules to take larger steps and make the IPM feasible in practice.The authors report orders-of-magnitude improvement in progress per step.Source: §1 · p. 3
• The solver is evaluated against LEMON’s cost scaling and network simplex solvers on 12 instance families from.The paper notes network simplex’s speed on small networks and cost scaling’s better scaling on larger inputs; both exhibit super-linear scaling on most benchmarks.Source: §1 · p. 3
• On grid wide, scaling is m^1.04 versus LEMON’s m^1.48, and the solver scales better than both LEMON solvers on vision rnd and grid wide.The comparison is notable despite the new solver’s larger constant factors and its status as a first, not-constant-optimized implementation.Source: §1 · p. 3

2 Almost-Linear-Time Min-Cost Flow in Theory

• The method starts from a strictly feasible flow and updates it by a circulation chosen to reduce the potential while keeping the flow strictly within capacity bounds.The potential measures closeness to optimal cost and feasibility; if the optimal cost is unknown, it can be found by binary search.Source: §2 · p. 4
• An approximate min-ratio cycle reduces the potential by Ω(1/m^o(1)), yielding convergence in m^1+o(1) iterations to a solution that can be rounded to minimum-cost flow.The update changes the flow by a circulation, preserving demand constraints, and approximations to gradients and lengths can be maintained because their multiplicative changes are bounded.Source: §2 · p. 4
• A low-stretch spanning tree gives, in expectation, a fundamental cycle approximating the min-ratio cycle, and O(log n) independent trees yield a good cycle with high probability.The argument combines expected stretch, Markov’s inequality, and a Chernoff bound.Source: §2 · p. 4
• As flow updates and lengths change, the data structure identifies off-tree edges whose cycles remain approximate min-ratio cycles or swaps edges to restore cycle quality.This addresses degradation of the spanning trees during the sequence of updates.Source: §2 · p. 4

3 Transfer to Practice

• Flow atoms preserve exact arithmetic for routed amounts, and the largest atom size that still yields the exact optimum is chosen to limit update cost.Discretization loses strict feasibility, so the implementation sets the threshold at the hard bounds to keep the exact optimum reachable.Source: §3 · pp. 4–5
• A cost scale λ replaces dependence on the unknown optimal cost, with a shifted barrier potential and periodic increases to λ.The horizontal barrier shift handles loss of strict feasibility, while the vertical shift makes the barrier minimum zero.Source: §3 · p. 5
• More accurate progress bounds split cycle progress into cost and forward- and reverse-barrier changes, enabling more carefully calculated steps and larger updates.The cost term remains exact; the forward barrier uses a second-order upper bound, and the reverse barriers are bounded by a matched, offset barrier term.Source: §3 · p. 5
• The implementation uses one shortest-path tree as a practical substitute for intricate low-stretch trees and evaluates its fundamental cycles for progress.It gathers cycle aggregates with depth-first search and modified union-find, then maintains the tree using a link-cut data structure.Source: §3 · p. 6
• A queue prioritizes off-tree edges by estimated progress, rechecking aggregates and step sizes when a candidate is popped because estimates can become stale.The progress bound uses six cycle aggregates, including cost, barrier derivatives, and minimum forward residual capacity.Source: §3 · pp. 5–6
• The implementation repairs its tree by edge swaps as lengths change, avoiding full recomputation; a phase refreshes all aggregates and queue estimates.The authors report that repairs were efficient enough that they never needed to redraw the tree from scratch.Source: §3 · p. 6
• Preprocessing splits the residual graph into strongly connected components so the implementation can maintain a root-reachability invariant for routing flow away from the root.The authors say this reduces unroutable cycles and leaves the optimum unchanged because flow cannot be routed between components.Source: §3 · p. 6
• The implementation doubles λ when many cycles make no progress, then terminates by rounding when doubling has little effect; it also uses heuristics to detect livelocks.For choosing ε, the authors report that a setting retaining at least 65% of progress for half of cycles that progress at ε = 1 works well.Source: §3 · pp. 6–7

Baselines and Instances

• The comparison uses LEMON 1.3.1 on five DIMACS core instance families and adds six sizes to grid wide.The family starts at α = 1, halves α until exact optimality, and warm-starts each subsequent size from the prior instance’s α.Source: §4 · p. 7

Plots

• The plots show edge count against runtime in seconds on logarithmic axes, with separate lines for the IPM and both LEMON solvers.Each family-size pair is measured once on a single seed.Source: §4 · p. 7

Failures

• On the largest road flow instance, smaller α values cause livelock, while grid long is retired after early termination persists on the next instance.The road flow run terminates only at its smallest terminating α and returns a feasible but non-optimal flow.Source: §4 · p. 7

Scaling

• Our exponent is competitive on five named families, wins on vision rnd and grid wide, and loses on netgen 8 and goto sr.For grid wide, the reported exponents are m^1.04 for the IPM and m^1.48 for the best baseline.Source: §4 · p. 8

Extending vision rnd

• On the extended LB07-bunny-lrg vision rnd family, the IPM has the best exponent but does not beat both LEMON solvers.LEMON cost scaling times out at 48 hours on the largest datapoint, while network simplex has a lower exponent than in the overview.Source: §4 · p. 8

Other Approaches

• The submitted design was selected over alternatives because it was more stable and scaled competitively across families and sizes.Alternatives included extending the barrier, different LSST designs, cost-scale and convergence rules, forest routing, and dynamic precision.Source: §5 · p. 8

Coupling Between Cost Scale and Precision

• Cost scale λ and precision α are coupled: exact-convergence strategies at either extreme showed poor scaling, making their joint steering important.Carefully increasing λ can converge at larger α but scales poorly; a constant λ with corresponding α also showed poor scaling in first experiments.Source: §5 · p. 8

Future Work

• Understanding why the IPM scales well on grid wide but struggles on other families could guide better steering rules.The authors identify this as future work toward practical almost-linear minimum-cost flow on more instance families.Source: §5 · p. 8

• AI was used for spelling and grammar checks, readability feedback, sentence polishing, literature discovery, and plot generation.Source: A AI Use · p. 13
• Agentic AI was used for coding and evaluation in some parts of the implementation.Source: A AI Use · p. 13
• AI generated the Appendix B proof, which two authors independently verified.Source: A AI Use · p. 13

• The proof applies Bernoulli’s inequality for negative exponents and Jensen’s inequality, with equality only when the arguments coincide.The equality condition in the Jensen step is x = 0 or all b_i equal.Source: B Proof of Upper Bound · pp. 13–14
• For γ > 0, the theorem gives unique positive a and b whose fitted curve matches the sum’s first two derivatives at zero.The fit is formed for S as a sum of terms f_a_i,b_i.Source: B Proof of Upper Bound · p. 13
• The fit lies strictly above S(x) for x > 0 unless all b_i are equal, when the fit equals S.Source: B Proof of Upper Bound · pp. 13–14
• Matching derivatives at zero yields b = M/W and then a = M b^(1+γ), establishing existence, uniqueness, and positivity.Source: B Proof of Upper Bound · pp. 13–14
• The result holds for every positive real γ, not only integer exponents.Source: B Proof of Upper Bound · p. 14