coin-or / CyLP

A Python interface to CLP, CBC, and CGL to solve LPs and MIPs.
Other
181 stars 67 forks source link

Big M problem #122

Closed rsemenoff closed 3 years ago

rsemenoff commented 3 years ago

So if you take the example in the documentation (http://coin-or.github.io/CyLP/modules/CyCbcModel.html) and add the following constraints:

isin = model.addVariable('isin', 3, isInt=True) model+= isin*1.01e6-x>=0 model+=0<=isin<=1

, it should arrive at the same optimum, but it doesn't.

The second line is just constraining x<1e6, effectively, i.e. 1e6 is the big M. You can comment that one line out and see the result being different. Here is the entire listing of the sample with my additional constraints:

import numpy as np from cylp.cy import CyCbcModel, CyClpSimplex from cylp.py.modeling.CyLPModel import CyLPModel, CyLPArray model = CyLPModel()

x = model.addVariable('x', 3, isInt=True)##RWS added to show issue isin = model.addVariable('isin', 3, isInt=True)##RWS added to show issue

y = model.addVariable('y', 2) A = np.matrix([[1., 2., 0],[1., 0, 1.]]) B = np.matrix([[1., 0, 0], [0, 0, 1.]]) D = np.matrix([[1., 2.],[0, 1]]) a = CyLPArray([5, 2.5]) b = CyLPArray([4.2, 3]) x_u= CyLPArray([2., 3.5])

model+= isin*1.01e6-x>=0.0 ##RWS added to show issue model+=0<=isin<=1 ##RWS added to show issue

model += Ax <= a model += 2 <= B x + D * y <= b model += y >= 0 model += 1.1 <= x[1:3] <= x_u

c = CyLPArray([1., -2., 3.]) model.objective = c x + 2 y.sum() s = CyClpSimplex(model)

cbcModel = s.getCbcModel() cbcModel.solve()

print (cbcModel.status) sol_x = cbcModel.primalVariableSolution['x']

(abs(sol_x - np.array([0, 2, 2])) <= 10**-6).all() print(sol_x)

tkralphs commented 3 years ago

Can you try dumping the model(s) to MPS and solving in stand-alone Cbc? If you see the same behavior, then this is a bug in Cbc and should be reported there.

rsemenoff commented 3 years ago

Ted, Thanks for the suggestion - I get the same variable values reported from the command line and the CyLP,  but the command line says infeasible whereas CyLP says solution.  Could CyLP be generating bad models? I wonder this because I notice that if I try to solve a pure linear problem, I getvery different results from the Simplex model vs the CBC model.I've attached the .mps files and the solution files - I can't tell by looking at them if the model is represented as expressedin the python code. If you have any suggestions as to how to verify the .mps file, I will try it out. I was thinking about making a sample model in a human-readable format like AMPL or GAMS and run it through NEOS,as I see they have the latest CBC installed, I suppose either one of these languages should be easy to learn.

On Monday, May 17, 2021, 08:31:17 PM PDT, Ted Ralphs ***@***.***> wrote:  

Can you try dumping the model(s) to MPS and solving in stand-alone Cbc? If you see the same behavior, then this is a bug in Cbc and should be reported there.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub, or unsubscribe.

tkralphs commented 3 years ago

I don't know exactly what's going on internally, but now that I'm looking more carefully, I see that you're multiplying a vector by a constant. I think I've seen this cause problems before, perhaps because the coefficients are usually expected to be left-multiplied, not right . If I change to

1.01e6*isin-x>=0.0

then it behaves as expected. Technically, I guess right-multiplying is the correct thing to do if you think of the constant as a 1x1, but I guess CyLP doesn't treat constants as a 1x1. We should probably throw an error when trying to do what you're doing if it's not going to result in a correct model, but I don't really have time to dig into it at the moment. This will get you going anyway.

rsemenoff commented 3 years ago

Ah ha, yes that works now.  Thank you! It makes sense. If the constant value is getting broadcast to a CyLPArray (and doing dot product), its already documented that the CyLPArray comes first in the expression, so actually might not be an issue. On Tuesday, May 18, 2021, 08:35:35 PM PDT, Ted Ralphs @.***> wrote:

I don't know exactly what's going on internally, but now that I'm looking more carefully, I see that you're multiplying a vector by a constant. I think I've seen this cause problems before, perhaps because the coefficients are left-multiplied, not right . If I change to 1.01e6*isin-x>=0.0

then it behaves as expected. Technically, I guess that right-multiplying is the right thing to do if you think of the constant as a 1x1, but I guess CyLP treats constants as a special case. It should probably throw an error when trying to do what you're doing. I don't really have time to dig into it at the moment, but this will get you going anyway.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub, or unsubscribe.

tkralphs commented 3 years ago

I'm going to close this issue for now, feel free to reopen if needed.