search

coin-or / python-mip

Python-MIP: collection of Python tools for the modeling and solution of Mixed-Integer Linear programs
Eclipse Public License 2.0
518 stars 92 forks source link

Speeding up problem set-up #146

Open hadware opened 3 years ago

hadware commented 3 years ago

This is more a suggestion than an actual bug.

I've been using python-mip and cvxpy on the same MI problems, and even though their actual problem-solving performances are usually the same (or even slightly better for python-mip), the problem set-up (defining the objective and the constraints) is much too slow with python-mip.

For instance, let's consider the following two implementations:

A : np.ndarray of shape (2000, 3069)
d : np.ndarray of shape (3069,)

# python-mip
# the following block takes 11 seconds to run
model = mip.Model("disorders")
x = model.add_var_tensor(shape=(len(d),), name="x", var_type=mip.BINARY)
model.objective = mip.minimize(x.T @ d)
model += (A @ x == 1)

# this line takes 0.13s to run
status = model.optimize()

# cvxpy
# this blocks takes 0.14s
x = cp.Variable(shape=len(d), boolean=True)
obj = cvxpy.Minimize(d.T @ x)
constraints = [A @ x == 1]
prob = cvxpy.Problem(obj, constraints)
optimal = prob.solve()

It's quite a shame that python-mip is limited by such a non-central part of MIP-solving, that is, setting up the problem. I used CProfile to figure out what's taking so much time when building the problem, and it seems that the repeated instance checks when multiplying a variable with a scalar are the culprits.

Maybe you could exploit the structured nature of LinExprTensor to have only one typecheck when multiplying then with arrays... or something like this.

h-g-s commented 3 years ago

Hi Yuri, Could you post your profile results indicating which lines are the bottleneck ?

Em qui, 1 de out de 2020 às 08:37, hadware notifications@github.com escreveu:

This is more a suggestion than an actual bug.

I've been using python-mip and cvxpy on the same MI problems, and even though their actual problem-solving performances are usually the same (or even slightly better for python-mip), the problem set-up (defining the objective and the constraints) is much too slow with python-mip.

For instance, let's consider the following two implementations:

A : np.ndarray of shape (2000, 3069) d : np.ndarray of shape (3069,)

python-mip

the following block takes 11 seconds to run

model = mip.Model("optimal_disorder") x = model.add_var_tensor(shape=(len(disorders),), name="x", var_type=mip.BINARY) model.objective = mip.minimize(x.T @ disorders) model += (A @ x == 1)

this line takes 0.13s to run

status = model.optimize()

cvxpy

this blocks takes 0.14s

x = cp.Variable(shape=len(disorders), boolean=True) obj = cvxpy.Minimize(d.T @ x) constraints = [A @ x == 1] prob = cvxpy.Problem(obj, constraints) optimal = prob.solve() It's quite a shame that python-mip is limited by such a non-central part of MIP-solving, that is, setting up the problem. I used CProfile to figure out what's taking so much time when building the problem, and it seems that the repeated instance checks when multiplying a variable with a scalar are the culprits.

Maybe you could exploit the structured nature of LinExprTensor to have only one typecheck when multiplying then with arrays... or something like this.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/coin-or/python-mip/issues/146, or unsubscribe https://github.com/notifications/unsubscribe-auth/AB4VZOWRM4HTM3KE62DBRFLSISO2LANCNFSM4SAS5ZPA.

hadware commented 3 years ago

In my case, it's this one: https://github.com/coin-or/python-mip/blob/2f23cec23f886847738f7fe2bbe9064e2ed86919/mip/entities.py#L591 (since i'm multiplying a LinExprTensor with a numpy array), but I'm guessing it could be any of the instance checks from the Variable class.

jurasofish commented 3 years ago

What sort of sparsity is the A matrix? Would be great if you could provide a complete functioning example. I had a go at reproducing it but sparsity has a massive impact on it.

My experience is that CVXPY tends to run a massive "transformation" step before passing data into CBC, which can only(?) be worked around with a parameterised model. Very interested if you another know a way of mitigating this. Might be impacted by which solver is in use (CBC vs whatever CVXPY has defaulted to in your case).

hadware commented 3 years ago

The A matrix was relatively sparse (I'm sending the two matrices as attachment).

I kind of looked into the cvxpy internals to see how their way of doing things could be transposed to python-mip, but ended up losing my patience after having starting to get lost in an endless stack of decorated decorators. My understanding of MIP is very light (to say the least), thus I'm not entirely sure of what would be the best way to solve this while keeping python-mip's code clean (which it currently is).

My current understanding of python-mip is that doing operations with Variables builds a tree-like nested structure of the operations involving each variables. Each operation, which each and every variable, involves one or several costly instance checks at "run-time" to build that structure. In the case of vectorized operations, this structure should likely be built at "solve"-time, while doing some mininal checking at run-time based on the strong homogeneity of the vectors involved in the operations. Does that look right to you?

jurasofish commented 3 years ago

Thanks I'll have a play with that later.

Yeah I've found that instance checking is a bit of a bottleneck for me too. In my problems model construction time is the vast majority of the work so I'm really looking for ways to cut it down short of moving to pypy/julia. cylp looks good but its internal postfix notation makes me extremely uncomfortable, and it doesn't seem to be able to handle constraints constructed from linear expressions (and it has very spartan documentation)

performing operations on variables in python-mip results in a linear expression which is a constant float plus a dictionary map from a variable object to the variable object's coefficient in the linear expression (plus sense for use in constraints). Operations between two linear expressions result in merging the dictionaries.

The dictionary is, essentially, a DOK sparse vector of coefficients for the linear expression. I think the linear expressions could be sped up considerably by implementing them as sparse matrices:

Because it's all just array sum and multiplication there shouldn't need to be any type checking (perhaps make it optional).

Would be challenging to shoehorn it into python-mip at this point, and vectorising it down to the level of adding constraints in the underlying cbc/gurobi solvers might be tricky.

I'm thinking of building this as a proof of concept with scipy sparse matrices. If I can get that working nicely I might try and bring some of the learnings into python-mip.

jurasofish commented 3 years ago

I'm thinking of building this as a proof of concept with scipy sparse matrices. If I can get that working nicely I might try and bring some of the learnings into python-mip.

That didn't really work, too much overhead in the sparse matrices