Hardware-agnostic linear programming on the GPU

Guillaume Dalle, Michael Klamkin, Simeon Schaub

2026-08-13

Slides

gdalle.github.io/JuliaCon2026

Motivation

Breakfast optimization

A bakery can produce croissants (1.2€) or pains au chocolat (1.5€):

	Butter	Sugar	Flour
Croissant	35g	20g	15g
Pain au chocolat	25g	30g	20g

The ingredients have unit costs and supply limits:

	Butter	Sugar	Flour
Unit cost	2€/kg	1€/kg	0.5€/kg
Max supply	19kg	14kg	12kg

How much of each should they bake to maximize profits?

Mathematical formulation

• Decision variables:
  • Quantities \(c, p\)
  • Ingredient weights \(b, s, f\)
• Objective:
  • Revenue
  • Cost
• Constraints:
  • Max supply
  • Ingredients

\[ \begin{align*} \max_{c, p, b, s, f} \quad & +(1.2c + 1.5p) \\ \phantom{\max_{c, p, b, s, f}} & - (2b + 1s + 0.5f) \\ \text{s.t.} \quad & 0 \leq c, p \\ & 0 \leq b \leq 19 \\ & 0 \leq s \leq 14 \\ \phantom{\max_{c, p, b, s, f}} & 0 \leq f \leq 12 \\ & 35c + 25p = 1000b \\ & 20c + 30p = 1000s \\ & 15c + 20p = 1000f \end{align*} \]

Linear programming

Optimize a linear objective subject to affine constraints:

\[\text{minimize} \quad c^\top x \quad \text{subject to} \quad A x \leq b\]

More useful than it seems:

• For convex problems (subroutine of Frank-Wolfe)
• For integer problems (subroutine of Branch & Bound)

Parallelization opportunities & challenges

Why parallelize?

• Because we have GPUs!
• To solve large instances
• To solve several instances

Why not parallelize?

• Best CPU solvers: simplex & interior points
• Sequential++ (sparse matrix factorizations)
• Applications require high precision

Visualization with lpviz.net by (Grand and Klamkin 2026).

Primal-dual algorithms

GPUs are best leveraged with first-order methods (Lu and Yang 2025a):

• Same intuition as deep learning: gradient only
• Core primitive: sparse matrix-vector product (SpMV)
• Recent proof of concept… coded in Julia (Lu and Yang 2025b)

CoolPDLP design

State of the art

Package	Hardware	Type-generic	Batchable
cuPDLPx	NVIDIA, AMD	No	No
HPR-LP-C	NVIDIA	No	No
HiGHS	NVIDIA	No	No
cuopt	NVIDIA	No	Yes
MPAX	agnostic (via JAX)	No	Yes
BatchPDLP.jl	NVIDIA	No	Yes
CoolPDLP.jl	Any	Yes	Yes

CoolPDLP is hardware-agnostic

• Coded with GPUArrays.jl and KernelAbstractions.jl.
• Runs on any GPU or CPU.
• Can use optimized linear algebra with e.g. cuSPARSE

CoolPDLP is type-generic

• The constraint matrix \(A\) in \(Ax \leq b\) is arbitrary
  • Dense
  • Sparse
  • Structured (!!!)
• The floating point precision is arbitrary
  • Float32
  • Float64
  • exotic types

CoolPDLP is batchable

Solves a bunch of problems in lockstep:

\[\text{minimize} \quad [c_1, ..., c_B]^\top x \quad \text{subject to} \quad A x \leq [b_1, ..., b_B]\]

The core primitive becomes sparse matrix-matrix product (SpMM).

Package demo

Set up

Useful imports

using CoolPDLP
import GPUSelect
import HiGHS, JuMP, MathOptInterface as MOI
import MathOptBenchmarkInstances as MOBI
using JLArrays, SparseArrays
ENV["DATADEPS_ALWAYS_ACCEPT"] = true

Precompiling packages...
   1150.7 ms  ✓ QuartoNotebookWorkerJSONExt (serial)
  1 dependency successfully precompiled in 1 seconds

true

GPU environment setup

backend = if parse(Bool, get(ENV, "CI", "false"))
    JLBackend()
else
    # GPUSelect.auto_install!()
    GPUSelect.Backend()
end

JLBackend(false)

Setting up a problem

dataset = MOBI.Netlib
name = MOBI.list_instances(dataset)[1]
qps, path = MOBI.read_instance(dataset, name);
milp = MILP(qps; dataset, name, path)

 Downloading artifact: netliblp

MILP instance 25fv47 from dataset Netlib:
- types:
  - values Float64
  - vectors Vector{Float64}
  - matrices SparseArrays.SparseMatrixCSC{Float64, Int64}
- variables: 1571
  - 1571 continuous
  - 0 integer
- constraints: 821
  - 305 inequalities
  - 516 equalities
- nonzeros: 10400

Solving on the CPU

algo = PDLP(
    Float32, Int32, SparseMatrixCSC;
    termination_reltol = 1.0f-1, time_limit = 30.0,
)
solve(milp, algo);  # compile
sol, stats = solve(milp, algo);
stats

Convergence stats with termination status OPTIMAL:
- KKT relative errors: primal 2.036e-02, dual 4.446e-02, gap 3.684e-02
- time elapsed: 0.006 seconds
- KKT passes: 100

Solving on the GPU

algo_gpu = PDLP(
    Float32, Int32, GPUSparseMatrixCOO; backend,
    termination_reltol = 1.0f-1, time_limit = 30.0,
)
solve(milp, algo_gpu);  # compile
sol_gpu, stats_gpu = solve(milp, algo_gpu);
stats_gpu

Convergence stats with termination status OPTIMAL:
- KKT relative errors: primal 2.036e-02, dual 4.446e-02, gap 3.684e-02
- time elapsed: 0.454 seconds
- KKT passes: 100

JuMP interface

model = JuMP.read_from_file(path; format=MOI.FileFormats.FORMAT_MPS)
JuMP.relax_integrality(model)
JuMP.set_optimizer(model, CoolPDLP.Optimizer)
JuMP.set_silent(model)
JuMP.set_attribute(model, "float_type", Float32)
JuMP.set_attribute(model, "int_type", Int32)
JuMP.set_attribute(model, "matrix_type", GPUSparseMatrixCOO)
JuMP.set_attribute(model, "backend", backend)
JuMP.set_attribute(model, "termination_reltol", 1.0e-1)
JuMP.set_attribute(model, "time_limit", 30.0)
JuMP.optimize!(model)
JuMP.is_solved_and_feasible(model)

┌ Warning: Got mismatched float type: solving in Float32 but returning the solution in Float64.
└ @ CoolPDLP ~/.julia/packages/CoolPDLP/fN3I2/src/MOI_wrapper.jl:225

true

Conclusion

Perspectives

• More functionalities: presolving, feasibility, polishing, infeasibility detection, crossover
• Kernel fusion with Reactant.jl
• Extension to nonlinear problems (quadratic first)
• Applications to decision-focused learning (Mandi et al. 2024)

References

Grand, Evan, and Michael Klamkin. 2026. “Lpviz: Interactive Linear Programming Visualization.” Pre-published April 30. https://doi.org/10.48550/arXiv.2604.27518.

Lu, Haihao, and Jinwen Yang. 2025a. “An Overview of GPU-based First-Order Methods for Linear Programming and Extensions.” Pre-published June 2. https://doi.org/10.48550/arXiv.2506.02174.

Lu, Haihao, and Jinwen Yang. 2025b. “cuPDLP.jl: A GPU Implementation of Restarted Primal-Dual Hybrid Gradient for Linear Programming in Julia.” Operations Research 73 (6): 3440–52. https://doi.org/10.1287/opre.2024.1069.

Mandi, Jayanta, James Kotary, Senne Berden, et al. 2024. “Decision-Focused Learning: Foundations, State of the Art, Benchmark and Future Opportunities.” Journal of Artificial Intelligence Research 80 (August): 1623–701. https://doi.org/10.1613/jair.1.15320.