Gurobi Bug v6.5.2

main.py

from gurobipy import *

# This function mimics the JuMP functionality whereby it always
# sets the objective coefficients to the JuMP representation of
# what the coefficients should be...
def JuMPsolve(m, objective_coefficients):
    # Set the objective coefficients
    m.setAttr("Obj", m.getVars(), objective_coefficients)
    # Update the model just to be safe
    m.update()
    # Optimize
    m.optimize()

def exampleproblem(c):
    # Intialise model
    m = Model("Gurobi Bug")
    
    # A variable in [0, 1]
    X1 = m.addVar(lb=0, ub=1)
    # A variable in [0, inf)
    X2 = m.addVar(lb=0, ub=GRB.INFINITY)
    # The special variable in (-inf, u] needed to trigger the bug
    X3 = m.addVar(lb=-GRB.INFINITY, ub=1, obj=1)

    # Choose our objective coefficients to always use
    objective_coefficients = [0., 0., 1.]

    # Lets maximise (although minimisation also works)
    m.modelSense = GRB.MAXIMIZE

    # Add all this to the Gurobi Model
    m.update()

    # A normal constraint
    C1 = m.addConstr(  X1 + X2      >= 0)
    
    # We need a few more constraints. You can just duplicate C1
    C2 = m.addConstr(  X1 + X2      >= 0)
    # although other constraints work like
    # C2 = m.addConstr(  X1 - X2      <= 1)

    # These are the two problem constraints
    C3 = m.addConstr(2*X1 - X2      >= 0)

    # Now say X1 >= X3
    # We intend to add an almost colinear constraint later
    C4 = m.addConstr(  X1      - X3 >= 0)

    # If you include this, it returns optimal instead of unbounded
    # The -10 is totally arbitrary. It should just be a valid bound
    # m.addConstr(X3 >= -10)

    # Solve the base problem
    JuMPsolve(m, objective_coefficients)

    # Add a new constraint
    C5 = m.addConstr(c*X1 + X2 >= X3)
    
    # Resolve the problem (and introduce the reformulation at the end)
    JuMPsolve(m, objective_coefficients)

    # If you include this, it returns optimal
    # The rhs is totally arbitrary
    m.addConstr(X3 <= 20)

    # Trigger the bad solve
    JuMPsolve(m, objective_coefficients)

    return

# Any numbers between
#exampleproblem(-1.000001)
# and
exampleproblem(-1.0000019)
# will trigger this.
# It's to do with the reformulation of almost co-linear constraints