Simd
Simd

Reputation: 21254

How to read an MPS file for scipy milp

Scipy 1.9 now supports milp. This is great news but to use it I need to read in an MPS file specifying the problem. Scipy doesn't seem to have a function to read in MPS files.

How can I do that?

Upvotes: 1

Views: 1119

Answers (1)

joni
joni

Reputation: 7157

As mentioned in the comments, this feature is not implemented yet. As a workaround, you could read the .mps file with a modelling framework like PuLP or python-mip and then extract all the data of the underlying MILP.

Here's a naïve (and thus slow) implementation of this approach by means of python-mip:

import numpy as np
from mip import Model
from scipy.sparse import csr_matrix

def read_mps(file: str):
    """
    Reads a .mps and saves all the data of the MILP:

    min c^T * x

    s.t. b_l <= A*x <= b_u
          lb <=   x <= ub
                x_i integer if integrality[i] = 1
    """
    mdl = Model(solver_name="CBC")
    mdl.read(file)

    # model parameters
    num_vars = len(mdl.vars)
    num_cons = len(mdl.constrs)

    # variable types and bounds
    lb = np.zeros(num_vars)
    ub = np.inf*np.ones(num_vars)
    integrality = np.zeros(num_vars)
    for i, var in enumerate(mdl.vars):
        lb[i] = var.lb
        ub[i] = var.ub
        if var.var_type != "C":
            integrality[i] = 1

    # objective
    c = np.zeros(num_vars)
    for i, var in enumerate(mdl.vars):
        if var in mdl.objective.expr:
            c[i] = mdl.objective.expr[var]
    if mdl.sense != "MIN":
        c *= -1.0

    # constraint coefficient matrix
    b_l = -np.inf*np.ones((num_cons))
    b_u = np.inf*np.ones((num_cons))
    row_ind = []
    col_ind = []
    data    = []
    for i, con in enumerate(mdl.constrs):
        if con.expr.sense == "=":
            b_l[i] = con.rhs
            b_u[i] = con.rhs
        elif con.expr.sense == "<":
            b_u[i] = con.rhs
        elif con.expr.sense == ">":
            b_l[i] = con.rhs
        for j, var in enumerate(mdl.vars):
            if var in (expr := con.expr.expr):
                coeff = expr[var]
                row_ind.append(i)
                col_ind.append(j)
                data.append(coeff)
    A = csr_matrix((data, (row_ind, col_ind)), shape=(num_cons, num_vars))
    return c, b_l, A, b_u, lb, ub, integrality

Then, you can read all the data and pass it to milp like this:

from scipy.optimize import milp, Bounds, LinearConstraint

c, b_l, A, b_u, lb, ub, integrality = read_mps("example.mps")
bounds = Bounds(lb, ub)
constraints = LinearConstraint(A, b_l, b_u)

res = milp(c=c, constraints=constraints, bounds=bounds, integrality=integrality)

Having said this, the above workaround is only useful for small to mid-sized model files.

By the way, it's worth mentioning that scipy.optimize.milp uses the HiGHS solver under the hood, for which an interface was also recently built into PuLP and there are plans to do the same for python-mip. Both modelling frameworks can read .lp and .mps files and pass them to the interfaced solvers, so there's no urgent need to wait for this feature in scipy. However, to be fair, PuLP only interfaces the executables of some solvers instead of shipping them with the python package and further assumes the end-user has them already installed. IMO, this isn't really a convenient approach.

Upvotes: 1

Related Questions