1

I need to solve a linear program (LP) with a mix of binary and continuous variables (MILP).
The data I use as constraints and objectives come from pandas data frame manipulations, so they are in matrix (or I should probably say numpy array) format, and the variables to solve for are sometimes 100's / 1000's.
I know how to use matrices and vectors to setup and solve MILP's in R, but as I coded part of the data handling in python, I wanted to include the LP part in the python script.

So I looked for MILP solvers for python, and I found this post:

Python Mixed Integer Linear Programming

and this one:

https://realpython.com/linear-programming-python/#what-is-mixed-integer-linear-programming

Apart from many links no longer working in the first post, my impression, after browsing through the docs of some of the suggested solvers, was that the LP functions in scipy, which take matrices and vectors as input, cannot solve MILP's (see also this post); on the other hand, interfaces like pulp etc., which do solve MILP's, seem to have their own very special input where one has to specify one by one the equations and name each variable. That would be very inconvenient, as I already have a simple numpy array where each column is a variable and each row is a vector of constraint coefficients.

Does anyone know of a python interface to a MILP solver that takes matrices and vectors as input, instead of requiring to write everything as explicit text equations?

Thanks!

PS I should probably mention that I am looking for interfaces to non-commercial solvers.

EDIT adding my current conclusion, after following the advice from users, for future record and in case it helps someone else

In short: as pulp does the job and comes with a (by all accounts) very good LP solver (CBC), I gave in and accepted to rewrite my problem as a set of equations and inequalities, instead of using the matrix format.
It seems to work, it's indeed fast, at least for the problems I tried so far.
My advice to anyone trying to use it on large problems with many variables would be: look up LpAffineExpression, it's very useful to put together complicated constraints. But first you need to define your variables using LpVariable. A variable can also be a list of variables (that's what stumped me initially).

Before trying pulp, I tried cvxopt, which is easy to install, but then it's true, like the users below said, it has its own bizarre format, and GLPK unfortunately is not a very fast solver.

When I still tried to stick to my matrix input format, I tried cvxpy. On paper this is very good, and in fact would allow me to write the problem as a sum of squares min, which is in fact what my original problem was (I turned it into a LP by replacing sum of squares with sum of abs values, and converting that to LP by known techniques).
I got all the way to writing my objective and constraints in matrix format, and then the solver failed... in fact even the example problem published in the cvxpy documentation failed with the same error. Big waste of time.

Anyway, lesson learned: no point fighting: sometimes you just have to do what python wants you to do ;)

Improve this question
17
  • Can you update a small example and the expected value of each variable, please? Commented Sep 7, 2021 at 9:16
  • Just write yourself a wrapper: matrix-form -> whatever-form-you need. This is trivial and keeps full freedom about picking your library. This is valuable, as all libraries have their downsides. Commented Sep 7, 2021 at 9:16
  • @sascha : thanks; thing is, I am quite new to python, as might have transpired from my post; so manipulating objects in complex ways to 'wrap' a matrix into whichever form the solver likes is not going to be 'trivial' at all for me :( In the meantime I found that lpSolve could be a solution, BUT then it only works with Python 2, apparently, and even there the installation is quite a palaver (and frankly from my previous work I also found lpSolve much slower than other solvers). I think for now I will have to stick to R. lpsolve.sourceforge.net/5.5/Python.htm Commented Sep 7, 2021 at 9:35
  • 1
    Symphony is (imho) a toy-solver compared to Cbc. (There is a reason it's not part of the benchmarks i linked) But well... if it works it works... You can use the same matrix-based interfaces in Python: e.g. Cylp and cvxpy (cp.Minimize(cp.sum_squares(A*x - b))). But you will hit much harder roadblocks there not related to the APIs. The former: good luck building it if you are on Windows. The latter: good luck with open-source solver support (besides GLPK). Commented Sep 7, 2021 at 13:27
  • 1
    Not sure if this is useful for you, but CVXR/CVXPY has a matrix interface (more flexible than just a single A,b,c matrix/vector representation). The R and Python versions are close in how you would set up a model. The main disadvantage (for me) is that it cannot handle multi-dimensional variables (like X(i,j,k)) but that does not seem a problem for your model. (Sidenote: for many practical models min cx, Ax=b is not a very useful representation. Hence most serious modeling tools allow some kind of algebraic input). Commented Sep 8, 2021 at 8:16

2 Answers 2

Reset to default
1

There are some specifics missing and no real definition for how close to exactly what you have as input must be accepted, but Google's OR tools mention a method to define problems for the various solvers they wrap through arrays. This may be workable for you.

For 100s/1000s of variables, this is probably not an option, but I have added an answer to the scipy.optimize question with a toy implementation of branch-and-bound on top of linprog which takes it's args in arrays. If your problem is highly constrained, it might give some answers in a day or two :).

Improve this answer
Sign up to request clarification or add additional context in comments.

1 Comment

Thank you, this is interesting! In the meantime I 'surrendered' to the syntax required by pulp, and I have to say it's indeed easier than a matrix format, especially when sparse. As for the execution time, I found that pulp's CBC accepts a time limit parameter, and from some tests it emerged that the CBC solution gets very close to a very satisfactory objective in much less time than it takes to reach absolute optimality, so we are always using that. I am sure further developments are possible, but for now this fits our needs.
0

For future reference: Make sure you have at least pip install scipy>=1.9.3 and use from scipy.optimize import milp as explained here using the HiGHS solvers via method='highs'

Improve this answer

Comments

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.