Primal Dual Interior Point Linear Program Solver

An interior point method implementation for solving linear programming problems using a primal-dual central path following algorithm. See the full report on Report.pdf.

Overview

This solver implements a generalized primal-dual path following algorithm to solve linear programming problems from the LPnetlib library. The implementation is based on algorithms from Robert Vanderbei's textbook "Linear Programming, Foundations and Extensions" (Chapter 20).

Problem Formulation

The solver handles linear programs in the form:

$$\min c^T x$$

subject to:

• $Ax = b$
• $l \leq x \leq u$

where $l$ and $u$ are lower and upper bounds on the variables.

Features

• Primal-dual central path following algorithm
• Automatic rank reduction for redundant constraints
• Convergence tracking with barrier parameter $\mu$
• Support for box constraints without expanding the constraint matrix
• Efficient handling of sparse matrices

Requirements

• Julia (tested with recent versions)
• Required packages:
  • MatrixDepot
  • SparseArrays
  • LinearAlgebra
  • Printf
  • Plots (for convergence visualization in solver_plotting.jl)

Installation

Install the required packages:

using Pkg
Pkg.add("MatrixDepot")
Pkg.add("SparseArrays")
Pkg.add("LinearAlgebra")
Pkg.add("Printf")
Pkg.add("Plots")  # Optional, for plotting

Usage

Basic Usage

using MatrixDepot

# Load a problem from LPnetlib
md = mdopen("LPnetlib/lp_afiro")
problem = convert_matrixdepot(md)

# Solve the problem
tol = 1e-8
sol = iplp(problem, tol; maxit = 1000, verbose = true)

# Check results
if sol.flag
    println("Converged")
    println("Objective: ", sol.cs' * sol.xs)
else
    println("Did not converge")
end

Parameters

• tol: Tolerance for convergence (default: 1e-8)
• maxit: Maximum number of iterations (default: 1000)
• verbose: Print iteration details (default: false)
• name: Problem name for plot filename (solver_plotting.jl only)

Return Structure

The IplpSolution struct contains:

• x: Solution vector for the original problem variables
• flag: Convergence flag (true = converged, false = did not converge)
• cs, As, bs: Standard form problem components
• xs: Solution in standard form
• lam: Lagrange multipliers
• s: Dual slack variables

Algorithm

The solver uses a primal-dual interior point method that:

1. Initializes all primal and dual variables to positive values

2. Iteratively follows the central path by:

   • Computing residuals $\rho = b - Ax - f$ and $\sigma = c - A^T y - s + h$
   • Updating barrier parameter $\mu = \delta \frac{\gamma}{n+m}$ where $\gamma = f^T v + p^T q + t^T s + g^T h$
   • Solving Newton system for search directions
   • Computing step length to maintain positivity
   • Updating all variables

3. Terminates when $\mu < 10^{-13}$ (indicating convergence)

Hyperparameters

• $\delta = 0.1$: Controls barrier parameter reduction rate
• $r = 0.99$: Step length reduction factor
• $\mu_{\min} = 10^{-13}$: Minimum barrier parameter for convergence

Testing

Run the provided test cases:

# Test on class-provided problems
test_problem("LPnetlib/lp_afiro", true, -4.6475314e02)
test_problem("LPnetlib/lp_brandy", true, 1.5185099e03)
test_problem("LPnetlib/lp_fit1d", true, -9.1463780924e03)
test_problem("LPnetlib/lp_adlittle", true, 2.2549496e05)
test_problem("LPnetlib/lp_agg", true, -3.5991767287e07)
test_problem("LPnetlib/lp_stocfor1", true, -4.1131976219e04)
test_problem("LPnetlib/lp_25fv47", true, 5.5018458883e03)
test_problem("LPnetlib/lpi_chemcom", false, 0.0)

The test suite can be enabled by setting the appropriate if conditions to true at the end of the solver file.

Performance

• Most problems converge within 300 iterations
• Large problems (20,000+ variables/constraints) typically solve in 1-2 hours on consumer hardware
• The solver achieves high accuracy, often matching results to 7+ significant digits

Example Results

For lp_maros_r7 (144,848 non-zero entries):

• Converged in 224 iterations
• Final $\mu = 8.3 \times 10^{-14}$
• Objective: 1.4971851664796441e6
• Expected: 1.4971851665e6

Visualization

The solver_plotting.jl version includes convergence plotting that tracks:

• Objective value $c^T x$ over iterations
• Barrier parameter $\mu$ over iterations

Plots are automatically saved as {problem_name}.png.

Files

• Report.pdf: Full report for the CS 520 (Computational Methods in Optimization) course of this project.
• project2.jl: Main solver implementation with optimality checking
• plotter.jl: Enhanced version with convergence visualization
• Both files contain the same core algorithm with different output features

Known Issues

• lp_ganges: Produces a different (more optimal) solution than the reference, matching CPLEX results instead
• Very large or dense problems may require increased iteration limits
• Problems with $\mu > 1$ after several hundred iterations are likely infeasible or unbounded

References

Vanderbei, R. J. (2020). Linear Programming: Foundations and Extensions (5th ed.). Springer.

• Available at: https://link.springer.com/book/10.1007/978-3-030-39415-8
• Relevant sections: Part III (pages 295-385)