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
514 stars 89 forks source link

Changing `Model.sense` makes (feasible) relaxation infeasible #285

Open Markus28 opened 2 years ago

Markus28 commented 2 years ago

I am trying to determine whether the feasibility region of an MIP-relaxation is bounded or not. To do that, I pick a random objective and both minimize and maximize it, checking whether one of the problems is unbounded. Unfortunately, python-mip tells me during minimization that an optimal solution has been found, during maximization it tells me that the problem is infeasible. I'm not sure why that would be the case. Is this a bug or is this somehow expected behavior?

To reproduce

I can reproduce this problem with instances from MIPLIB, e.g. dano3mip:

# `model` is an MIP from MIPLIB, with the original objective etc.
status_min = model.optimize(relax=True) # OptimizationStatus.OPTIMAL
model.sense = mip.MAXIMIZE
status_max = model.optimize(relax=True) # OptimizationStatus.INFEASIBLE

I am happy to provide a complete minimal working example if that should be necessary.

System info

I am running:

sebheger commented 2 years ago

Hey @Markus28,

sounds strange. I know that there is/was some issue with simply changing the sense of the model.

Could you please use the latest version 1.14.0 and rerun your application? Is it happening with CBC or gurobi or maybe both? If CBC, pls also post the logs of both runs. What happens if you change the sense before the first .optimize() call?

Depending on the results and answers I will make some research.

Cheers, Sebastian

Markus28 commented 2 years ago

Thanks for the quick response

Here are the CBC logs for danoint:

The initial minimization (feasible):

Welcome to the CBC MILP Solver 
Version: Trunk
Build Date: Oct 24 2021 

Starting solution of the Linear programming relaxation problem using Dual Simplex

Coin0506I Presolve 600 (-64) rows, 457 (-64) columns and 3440 (208) elements
Clp0014I Perturbing problem by 0.001% of 1.1861368 - largest nonzero change 0.00069889106 ( 0.058921623%) - largest zero change 0.00068473773
Clp0000I Optimal - objective value 62.63728
Coin0511I After Postsolve, objective 62.63728, infeasibilities - dual 0 (0), primal 0 (0)
Clp0032I Optimal objective 62.63728042 - 611 iterations time 0.022, Presolve 0.00

and the subsequent maximization (infeasible):

Clp0014I Perturbing problem by 1e-07% of 1.1861368 - largest nonzero change 0 ( 0%) - largest zero change 0.00099428649
Clp0002I Dual infeasible - objective value 84.66802

I'm not sure whether I am interpreting the last log correctly, but knowing that the dual is infeasible isn't sufficient to conclude that the primal is infeasible, right? I don't know how the solvers work under the hood but to my naive mind, it looks like CBC only determines that the primal is either unbounded or infeasible (because dual is infeasible) and python-mip might not communicate this ambiguity correctly?

I also tried forcing CBC to use the primal simplex method by setting model.lp_method=mip.constants.LP_Method.PRIMAL but the logs stayed the same (in particular it still said that the Dual Simplex method was being used).

rschwarz commented 2 years ago

To be clear: the expectation is that in one case you get OPTIMAL, and in the other case UNBOUNDED?

Markus28 commented 2 years ago

I don't know whether the maximization problem is actually unbounded, that's what I want to find out (it might also be optimal, which occurs for some MIPLIB instances). But maximization being infeasible while minimization is optimal makes no sense.

rschwarz commented 2 years ago

OK, clear. I asked because it could be that maybe CBC does not differentiate between UNBOUNDED and INFEASIBLE (but I think it does).

Markus28 commented 2 years ago

I have tested it with Gurobi and maximization should indeed be unbounded.

odow commented 2 years ago

I know the answer to this, because it has bitten me before. This is the problem: https://github.com/coin-or/python-mip/blob/08b3ad46259406b45fa1845dfa89ea1557a4fdd0/mip/cbc.py#L1251-L1252 For LPs isProvenInfeasible can mean dual infeasible.

Here's what we do in Cbc.jl: https://github.com/jump-dev/Cbc.jl/blob/3db1f16751a18d4f200b15bfd6116925fbdfddb0/src/MOI_wrapper/MOI_wrapper.jl#L644-L650

tkralphs commented 2 years ago

For LPs isProvenInfeasible can mean dual infeasible.

Then this is just a straight-up bug in Cbc. I would assume it knows when it has an LP that is dual infeasible and should properly report unboundedness. Even in the case of an MILP, if the LP relaxation is dual infeasible, this means the problem is unbounded. This should be reported over in Cbc. I guess it should be an easy fix.

odow commented 2 years ago

Even in the case of an MILP, if the LP relaxation is dual infeasible, this means the problem is unbounded

Does it?

min y
0.1 <= x <= 0.9
x in Integer

is MILP infeasible, but the LP relaxation is unbounded.

tkralphs commented 2 years ago

Sorry, yes, you're right (of course). In this case, the recession cone associated with the convex hull has a ray with negative cost (assuming minimization), so if the convex hull itself is non-empty, the problem is unbounded, otherwise infeasible.... which is exactly the reason for that "unbounded or infeasible" return status. Still, this is a bug in Cbc if it says straight infeasible.

odow commented 2 years ago

Actually I wonder if the problem is just that I should have checked ContinuousUnbounded first?

    elseif Cbc_isContinuousUnbounded(model) != 0
        return MOI.INFEASIBLE_OR_UNBOUNDED

I think it's okay to return true for isProvenInfeasible, because it isn't isProvenPrimalInfeasible. This could arguably be a documentation issue.

tkralphs commented 2 years ago

Ah yes, that should be correct, now that you say it. It's not ideal that you need to check two different statuses.

odow commented 2 years ago

Maybe not:

julia> model = Cbc_newModel()
Ptr{Nothing} @0x00007f962d106c00

julia> Cbc_addCol(model, "x", -Inf, Inf, 1.0, 0, 0, C_NULL, C_NULL)

julia> Cbc_solve(model)
Dual infeasible - objective value 0
DualInfeasible objective 0 - 0 iterations time 0.002
1

julia> Cbc_isContinuousUnbounded(model)
0

julia> Cbc_isProvenInfeasible(model)
1

julia> Cbc_deleteModel(model)

julia> unsafe_string(Cbc_getVersion())
"2.10.5"