MPS-import for scipy.linprog benchmarking

mps_import.py

""" Hacky import & solving of MPS-files (only tested for netlib)
    - source of mps-files: http://www.netlib.org/lp/data/
        - instances need to be uncompressed -> see readme in link-folder
    - uses cvxopt's mps-parsing
        - internal problem-modification / canonicalization by cvxopt "_inmatrixform"
    - only implemented / tested for problems without explicit bound-/range-constraints
        - see problem table within readme in link-folder for instance-type table
    - probably ugly/unnecessary dense/sparse-stuff (code is old and hacky)
    - just to be careful: cvxopt's license is not compatible with scipy!

    This code was used privately for benchmarking an own IPM-solver. It's
    probably a very bad idea to call linprog for some of these not-that-small
    problems because linprog seems to be using dense-matrices.

    Alternative mps-reading approach (idea only; no implementation):
        https://github.com/biosustain/swiglpk (python-bindings for GLPK)
        see also (usage-issue / question): https://github.com/biosustain/swiglpk/issues/19

        This might be an interesting route, because there seems to be less
        internal problem-modification's and GLPKs MPS-reader is probably
        more mature too.

    Sascha-Dominic Schnug 2017
"""


import numpy as np
import scipy.sparse as sp
import scipy.optimize as spopt
from cvxopt.base import matrix, spmatrix
from cvxopt.modeling import op
from cvxopt.solvers import lp


def read_mps_preprocess(filepath):
    problem = op()
    problem.fromfile(filepath)
    mat_form = problem._inmatrixform(format='dense')
    format = 'dense'
    assert mat_form
    lp, vmap, mmap = mat_form

    variables = lp.variables()
    x = variables[0]
    c = lp.objective._linear._coeff[x]
    inequalities = lp._inequalities
    G = inequalities[0]._f._linear._coeff[x]
    h = -inequalities[0]._f._constant
    equalities = lp._equalities
    A, b = None, None
    if equalities:
        A = equalities[0]._f._linear._coeff[x]
        b = -equalities[0]._f._constant
    elif format == 'dense':
        A = matrix(0.0, (0,len(x)))
        b = matrix(0.0, (0,1))
    else:
        A = spmatrix(0.0, [], [],  (0,len(x)))  # CRITICAL
        b = matrix(0.0, (0,1))

    c = np.array(c).flatten()
    G = np.array(G)
    h = np.array(h).flatten()
    A = np.array(A)
    b = np.array(b).flatten()

    return c, G, h, A, b

""" LINPROG """
def solve_with_linprog(problem):
    c, G, h, A, b = problem
    G = sp.csc_matrix(G)
    A = sp.csc_matrix(A)
    result = spopt.linprog(c, G.todense(), h, A.todense(), b)
    return result

""" CVXOPT """
def scipy_sparse_to_spmatrix(A):
    A = sp.csc_matrix(A)  # hack
    coo = A.tocoo()
    SP = spmatrix(coo.data.tolist(), coo.row.tolist(), coo.col.tolist(), size=A.shape)
    return SP

def solve_cvxopt(problem, solver='conelp'):
    # solver='mosek' works too if mosek is available
    c, G, h, A, b = problem
    c = matrix(c)
    G = scipy_sparse_to_spmatrix(G)
    h = matrix(h)
    A = scipy_sparse_to_spmatrix(A)
    b = matrix(b)
    result = lp(c, G, h, A, b, solver=solver)
    return result['primal objective']

""" RUN """
afiro = read_mps_preprocess('afiro_.mps')
print('linprog: ', solve_with_linprog(afiro))
print('cvxopt conelp: ', solve_cvxopt(afiro))