Empirical Asymptotic Runtime Analysis of Linear Programming Algorithms

Abstract: This paper takes an empirical look at asymptotic runtime growth rates for the most widely used algorithms for solving linear programming (LP) problems across a set of six optimization application areas that are known to produce large and difficult LP models. On the algorithm side, we consider the simplex method, interior-point methods, and PDHG. On the model side, we use a LLM to create families of instances in different application areas, allowing us to study model types and sizes that are simultaneously synthetic and realistic. The results indicate that simple regression models typically predict observed runtimes quite well within a model class, and that asymptotic behavior can vary significantly between the different algorithms. This may have a significant impact on which algorithms will be most effective for solving large LP models in the future.

• The paper presents an extensive empirical analysis of LP algorithms, comparing simplex, interior-point, and PDHG methods.
• It employs LLM-generated synthetic benchmarks and power-law OLS regression to quantify the runtime scaling of solver components.
• The study reveals tradeoffs between scalability and solution accuracy, emphasizing crossover bottlenecks and the need for improved refinement methods.

Empirical Asymptotic Runtime Analysis of Linear Programming Algorithms

Introduction

This paper presents an extensive empirical investigation into the asymptotic runtime behavior of the leading algorithms for solving linear programming (LP) problems: simplex, interior-point methods, and the Primal-Dual Hybrid Gradient (PDHG) method. The study is unique in its use of LLMs to automatically generate synthetic, yet structurally realistic, parametric LP instance families in six high-impact application domains. Through ordinary least squares (OLS) regression, with a focus on power-law fits, the work elucidates the scaling dynamics of solver runtimes as a function of constraint matrix nonzeros, contributing key insights relevant to both algorithmic selection and expectations for future hardware and model size trends.

Algorithms and Their Architectural Implications

The simplex method, interior-point methods (IPMs), and PDHG represent substantially different computational architectures and tradeoffs between iteration cost, parallelizability, and solution accuracy:

• Simplex: Basis-driven, with cycling through adjacent vertices; strong solution accuracy (mean absolute constraint violation ≈ $10^{-10}$) and superior support for reoptimization, but fundamentally sequential and difficult to parallelize at scale.
• Interior-Point: Employs repeated solutions of weighted normal equations, leveraging sparse Cholesky factorization and fill-reducing orderings. IPMs exhibit local quadratic convergence, are markedly amenable to multi-core parallelism in the factorizations, but depend on potentially non-parallelizable preprocessing (ordering) steps, and often require post-processing (crossover) to obtain basic feasible solutions. Achievable accuracy is typically $10^{-8}$ (without crossover).
• PDHG: First-order, matrix-vector-multiplication-dominated method highly suitable for GPU acceleration. Each iteration is computationally cheap, but convergence (and thus solution accuracy) is slow; for practical tolerances, moderate ($10^{-4}$ with academic norms, or $10^{-6}$ for commercial solvers). Crossover can mitigate accuracy deficits but may be bottlenecked by low-quality intermediate solutions.

The study enforces commercially-relevant convergence criteria across algorithms. A notable aspect is the detailed quantification of solution accuracy gaps between the approaches, emphasizing practical implications for downstream applications.

Synthetic Model Generation and Benchmarking Protocol

To enable controlled, scalable experimentation, the author utilizes ChatGPT 5.4 to synthesize benchmark-generating scripts in six application domains: Airline Fleet Assignment, Gas Network Optimization, Production Planning, Supply-Chain Network Design, Telecommunication Network Design, and Electrical Power Unit Commitment. For each domain, model families range over parameter values ("benchmark ladders") yielding instances from trivially solvable sizes up to the largest feasible within hardware constraints (up to tens of millions of constraints/nonzeros after presolve).

Each algorithm is benchmarked across these model families. The core performance metric is wall-clock runtime, regressed against constraint matrix nonzeros using the power law form $runtime = a \times (\text{nnz})^b$, with $b$ and its confidence bounds extracted via log-log OLS. $R^2$ (in log space) quantifies fit quality. The study further decomposes total runtime across major algorithmic subcomponents (presolve, ordering, iterations, and crossover).

Empirical Results

Iterative Phase Scaling

For sufficiently large models, iteration-dominated runtime is well-described by power laws in nearly all cases ($R^2 > 0.9$). Key empirical findings include:

• Dual simplex: Exhibits the highest scaling exponents (often $b\in [1.5,2.1]$), consistent with prior theoretical and empirical reports of poor parallel scaling, especially at large instance sizes.
• Interior-point methods: Typical exponents are in $[1.5,1.9]$ for iterations, with fill-reducing ordering showing even higher exponents in some cases (notably for the Fleet and GasNet sets), potentially offsetting the parallel performance of factorizations.
• PDHG: Delivers the lowest exponents (commonly $10^{-8}$0 for well-fitted cases), reflecting superior scaling, but at the cost of lower post-iteration accuracy. Exceptionally noisy iteration times in GasNet suggest model-dependent sensitivity to data conditioning.

Crossover and Subcomponent Bottlenecks

Crossover, essential for transforming interior or PDHG solutions to basic feasible solutions, exhibits exponents up to $10^{-8}$1, occasionally exceeding those for iteration phases. This is problematic for PDHG, where low-accuracy final iterates can lead to pathological crossover performance—observed empirically as failure to produce basic solutions for large GasNet and several Fleet instances.

Ordering and presolve exhibit exponents between ~1 and 2.4, with fit quality universally high except for TelecomND ordering, likely due to model/data idiosyncrasies.

Absolute (Finite-Size) Performance

On state-of-the-art hardware (EPYC 9575F CPUs, RTX Pro 6000 GPUs), absolute performance (wall-clock time) depends strongly on the algorithm-model pair and presence/absence of crossover. Without crossover, PDHG (on GPU) is faster than interior-point (on CPU) for GasNet and ProdPlan; simplex/interior-point is sometimes preferable on Fleet, SCND, TelecomND, and UnitCommit.

With crossover enabled, interior-point plus crossover is more robust—PDHG often fails crossover or is vastly slower due to initial accuracy issues. In many cases, crossover time dominates total runtime (up to 95% for PDHG-SCND), making it a significant limiting factor as models scale further.

Implications and Theoretical Reflections

This comprehensive study demonstrates a persistently sharp tradeoff between asymptotic scalability and solution accuracy in commercial LP solvers. PDHG shows consistently favorable scaling, but its practical advantage is compromised by weak solution accuracy and severe crossover bottlenecks, particularly as model size and complexity increase. Interior-point methods strike a practical balance but are susceptible to poorly scalable subcomponents (ordering, crossover).

These findings reinforce several theoretical and practical implications:

• For ultra-large LPs where moderate accuracy is sufficient, first-order GPU-friendly methods like PDHG may become dominant, provided improved crossover or solution polishing techniques can be developed.
• For applications demanding high-precision basic solutions (e.g., for MIP branch-and-bound, sensitivity analysis), simplex/interior-point with robust crossover remains essential, but hardware and algorithmic bottlenecks in sequential steps must be addressed.
• Use of LLM-generated synthetic models is a promising approach for scaling benchmarking but raises questions about data realism and the transferability of results to specific industry applications.

The study highlights open questions on developing high-quality low-accuracy-to-basic-solution procedures, improving parallelizability in presolve and crossover, and establishing better realism metrics for synthetic LP data.

Conclusion

The empirical analysis establishes that different LP algorithms demonstrate distinct asymptotic scaling profiles when subjected to realistic, LLM-generated large-scale benchmarks. PDHG offers superior theoretical scaling at the expense of solution accuracy, with crossover posing a critical bottleneck. Simplex and interior-point methods provide higher accuracy, but their scalability is fundamentally limited by sequential and combinatorial subroutines. Achieving both high accuracy and superior scalability will require new algorithmic innovations, particularly around low-accuracy solution refinement and high-performance crossover.

This work serves as a foundation for both practitioners seeking to select LP algorithms for large-scale applications and researchers aiming to design next-generation solvers optimized for emerging hardware and unprecedented problem sizes.