Solver claims solution is optimal, but violates constraints

[matt-telstra]

tldr

I'm trying to add booleans together, hoping that the result is a boolean OR. I'm expecting a boolean OR, XOR, or an error, or INFEASIBLE. But what happens is that the solver violates other constraints to satisfy the addition.

Minimum working example

import cvxpy
from cvxopt.modeling import op
from cvxopt.modeling import variable, op, max, sum
import cylp

A = cvxpy.Bool(name='A')
B = cvxpy.Bool(name='B')
C = cvxpy.Bool(name='C')

constr = [] # constraints list

constr.append(A == 1)
constr.append(B == 1)

# trying to implement A or B
# not sure if addition works as boolean or
constr.append(C == A + B)

# I also tried this
#constr.append(1.0*C == 1.0*A + 1.0*B)

# I don't care what the objective is,
#all states are constrained to a single value
obj = cvxpy.Maximize(A)

problem = cvxpy.Problem(obj, constr)

ret = problem.solve(solver=cvxpy.CBC)

print('problem.status == %s' % problem.status)

assert(problem.status in [cvxpy.OPTIMAL, cvxpy.OPTIMAL_INACCURATE])

print('A:%f | B:%f | C:%f' % tuple([x.value for x in [A,B,C]]))

Expected output

• A == 1 and B == 1, because that's what the constraints explicitly state.
• If the solver does not know how to sum two true booleans (state C), it should throw an error, or return INFEASIBLE.

Observed output

• A == 0 and B == 0, even though I specified exactly that through the constraints
• Solver claims solution is 'optimal'

Comments

My third constraint seems to be causing issues. I understand that summing booleans is a bit odd, I was just wondering/hoping that it would work.
The issue here is that when the solver encounters this odd third constraint, it violates the first two constraints and claims the solution is optimal.

When I run this, I do see some FutureWarnings. They sound like they are not of consequence, but could they be the cause of the problem?

The output when I run this code is:

matt@machine:~/muckaround/opt$ python simple2.py
/usr/local/lib/python2.7/dist-packages/cvxpy/problems/solvers/cbc_intf.py:143: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  x = model.addVariable('x', n)
/usr/local/lib/python2.7/dist-packages/cylp/py/modeling/CyLPModel.py:193: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  if (other == None):
/usr/local/lib/python2.7/dist-packages/cylp/py/modeling/CyLPModel.py:516: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  if self.lower == None:
/usr/local/lib/python2.7/dist-packages/cvxpy/problems/solvers/cbc_intf.py:153: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  model += A[0:dims[s.EQ_DIM], :] * x == b[0:dims[s.EQ_DIM]]
problem.status == optimal
A:0.000000 | B:0.000000 | C:0.000000

[tkralphs]

Thanks for the detailed report. Just to be sure that this is a problem with cylp and not with the interface between cvxpy and cylp, can you verify that this happens when building the model directly in cylp? Otherwise, I guess the problem should be reported to cvxpy for investigation first. Once it can be verified that the translation of the model by cvxpy is really correct and the problem is with cylp, we can try to look into it further.

[matt-telstra]

How can I have booleans in cylp?

The documentation does not mention bool anywhere.

I tried just constraining integers to 0<=x<=1. But the code takes forever to run.

(I haven't written code directly for cylp before, so it's possibly wrong)

import numpy as np
from cylp.cy import CyClpSimplex
from cylp.py.modeling.CyLPModel import CyLPArray

s = CyClpSimplex()

# Add variables
x = s.addVariable('x', 4)
y = s.addVariable('y', 4)
z = s.addVariable('z', 4)

s.setInteger(x)
s.setInteger(y)
s.setInteger(z)

# Add constraints
s += x <= 1
s += x >= 0
s += y <= 1
s += y >= 0
s += z <= 1
s += z >= 0
s += 0 == x + y - z

# Set the objective function
s.objective = sum(z)

# Solve using primal Simplex
s.primal()
print(s.primalVariableSolution['x','y','z'])

[matt-telstra]

For the record, I found a way of implementing boolean ORs which doesn't involve summing booleans.

Basically for C = A or B
C <= A + B and 2C >= A + B

Here's the updated MWE.

Although the issue that I raised this ticket for is still valid.

import cvxpy
from cvxopt.modeling import op
from cvxopt.modeling import variable, op, max, sum
import cylp

num_states = 4

A = cvxpy.Bool(num_states,name='A')
B = cvxpy.Bool(num_states,name='B')
C = cvxpy.Bool(num_states,name='C')

constr = [] # constraints list

for i in range(num_states):
    constr.append(A[i] == i%2)
    constr.append(B[i] == (i%4 >= 2))

# trying to implement A or B
# not sure if addition works as boolean or
constr.append(1.0*C <= 1.0*A + 1.0*B)
constr.append(2.0*C >= 1.0*A + 1.0*B)

# I also tried this
#constr.append(1.0*C == 1.0*A + 1.0*B)

# I don't care what the objective is,
#all states are constrained to a single value
obj = cvxpy.Maximize(sum(C))

problem = cvxpy.Problem(obj, constr)

ret = problem.solve(solver=cvxpy.CBC)

print('problem.status == %s' % problem.status)

assert(problem.status in [cvxpy.OPTIMAL, cvxpy.OPTIMAL_INACCURATE])


for i in range(num_states):
    print('A:%f | B:%f | C:%f' % tuple([x[i].value for x in [A,B,C]]))

[tkralphs]

OK, I didn't catch before how you expected a bool variable to behave. As I understand (and this is the usual case in mathematical optimization), a bool variable is just a regular real-valued variable that is constrained to take values 0 or 1, not a logical variable, as you expected. With constraints A == 1 and B == 1, the expression A+B evaluates to 2. This is what the documentation says: http://www.cvxpy.org/en/latest/tutorial/advanced/#mixed-integer-programs.

Your original model should have been infeasible. Thus, there is a bug, but it's not clear to me whether it's a bug in cvxpy or cylp. I would report this over there first.

[matt-telstra]

If I use the default solver, then it returns infeasible_inaccurate.

How can I write this directly in cylp? How do I write bools into cylp? How does cvxpy convert it's booleans into cylp?

[tkralphs]

I'm not sure what you mean by the "default solver," but I infeasible is the proper return. To create a bool variable (in the sense that both cvxpy and CyLP mean that term), you would just create an integer variables with bounds of 0 and 1. Since cvxpy understands bool in the same way as CyLP, I assume it just does exactly the same thing, but this is a question and a discussion that should really be happening over in cvxpy because the issue doesn't seem to be with CyLP from what I gather.

[matt-telstra]

When I say 'default solver', I mean using

ret = problem.solve()

instead of

ret = problem.solve(solver=cvxpy.CBC)

With the former, it returns infeasible. With the latter, it returns a solution that violates constraints.
Doesn't that suggest that this is a cylp issue?

[GiorgioBalestrieri]

@matt-telstra not sure if you ever sorted this out, I believe @tkralphs was pointing to the fact that there are two interfaces involved (and therefore potentially two sources of errors): between cvxpy and cylp, and between cylp and CBC.

[bill-skinner]

I am also encountering this issue. CBC is returning a solution it claims is optimal, but in reality violates a constraint. I have constraints which sum binary variables.