Tuesday, October 4, 2011

Solving a large sparse system of equations with an LP solver

When we want to solve a large, sparse system of linear equations

image

using an LP solver we expect this to be reasonably fast. After all, solving systems of equations is what an LP solver does all the time. So we introduce a zero objective and solve. With Gurobi a sample problem gives:

Starting Gurobi...
Optimize a model with 9801 rows, 9801 columns and 48609 nonzeros
Presolve removed 4599 rows and 4599 columns
Presolve time: 0.06s
Presolved: 5202 rows, 5202 columns, 47242 nonzeros

Iteration    Objective       Primal Inf.    Dual Inf.      Time
       0      handle free variables                          0s
    3958      handle free variables                          5s
    4975      handle free variables                         11s
    5199    0.0000000e+00   0.000000e+00   0.000000e+00     25s
    5199    0.0000000e+00   0.000000e+00   0.000000e+00     25s

Solved in 5199 iterations and 25.19 seconds

Actually we can much improve on this. By specifying an advanced basis. We know the equations (slacks) will be non-basic in the solution and all the variables will be basic. If we set up an advanced basis for this we see for the same problem:

Starting Gurobi...
Optimize a model with 9801 rows, 9801 columns and 48609 nonzeros
Iteration    Objective       Primal Inf.    Dual Inf.      Time
       0    0.0000000e+00   0.000000e+00   0.000000e+00      0s

Solved in 0 iterations and 0.45 seconds

Indeed no simplex iterations were needed: just an invert.

 

Update

In response to a comment about Octave: that is still faster:

octave-3.2.4.exe:3> tmp=load("-ascii","mat.txt");
octave-3.2.4.exe:4> b=load("-ascii","rhs.txt");
octave-3.2.4.exe:5> a=sparse(tmp(:,1),tmp(:,2),tmp(:,3));
octave-3.2.4.exe:6> nnz(a)
ans =  48609
octave-3.2.4.exe:7> tic();x=a\b;toc()
Elapsed time is 0.057 seconds.
octave-3.2.4.exe:8>

Notes:

  • Not completely fair of course: the LP solver does still a lot of extra things (scaling, pricing etc.).
  • Note that I used a sparse matrix in Octave. Not sure what method they use to solve this system.
  • This was a structured test matrix: see: http://www.public.iastate.edu/~akmitra/aero361/design_web/Laplace.pdf. For more unstructured matrices Gurobi may do much better relative to Octave.