Jump.jl: Ease the computation of IIS?

[dourouc05]

Many solvers present an option to compute an IIS in case of infeasibility (Gurobi, CPLEX, BARON at least, probably others that I don't know of). However, accessing this functionality from JuMP is a bit of a mess (like the example iis.jl).

I think such functionality should be directly available from JuMP (when the underlying solver supports it), with a function like computeIIS, but I have seen no previous work in this direction. It could have two kinds of outputs:

• either directly the IIS in the console (as in the example)
• or a new JuMP model with only the constraints and variables of the IIS (so that it could be exported as an LP file), without console output (other than what the solver does when computing the said IIS)

Which solution would be preferable? I think that the second is more general (you could get the first one by just printing the resulting model), although I am not sure having the whole model is really useful for the user (except probably that it would be consistent).

Of course, I propose myself as a potential implementer of the needed changes (in JuMP, MathProgBase, and a few solvers).

[mlubin]

I would think of an IIS as a mapping from variables and constraint references to true/false indicating if they participate in the IIS or not. I would be willing to consider making this mapping easier to access through JuMP/MPB.

Pretty printing and constructing submodels is a level of complexity that I would prefer not to be responsible for maintaining and could be included in a separate model diagnostics package.

[mlubin]

For the record this is a duplicate of https://github.com/JuliaOpt/JuMP.jl/issues/235.

[dourouc05]

In what form would you return this mapping? I can think of several ways:

• a inIIS function that might be called on each and every variable and constraint (or array) — for a variable, it could also return, as second and third arguments, the upper and lower bounds —, after something like computeIIS has been called on the model. This would be very close to how Gurobi works.
• computeIIS calls the solver and returns two vectors, one with the variables (tuple: variable, lower bound, upper bound), one with the constraints. This would be very close to how CPLEX works.

What elements should be considered when making a choice here? Does one affect more performance than the other (is this kind of question relevant for IISes, which are not really used that often)?

A bit of loud thinking about the implementation.

• A removed bound would be represented as $\pm\infty$.
• The data returned when the solver has found the IIS would be stored in new fields of Model (linconstrIis, quadconstrIis, sosconstrIis, socconstrIis, sdpconstrIis), arrays of Booleans.
• Nonlinear models would not be supported at first.

For MathProgBase, I would need define a few functions: one to start the computation of the IIS, one to get the UP/LB/IISness of a variable, IISness of a constraint (following the conventions of MathProgBase, which seem to be having functions that work on one item at a time for add… and get…). Would it be enough to add them like this?
https://github.com/JuliaOpt/MathProgBase.jl/blob/master/src/SolverInterface/SolverInterface.jl#L56

• iis!
  https://github.com/JuliaOpt/MathProgBase.jl/blob/master/src/SolverInterface/LinearQuadratic.jl#L10
• getvarIISUB
• getvarIISLB
• getconstrinIIS (particularly unsure about this name, though: isconstrinIIS?)

(I prefer to discuss this kind of details now before having to rewrite thing… Well, functions names are easy to change, but not really the overall architecture.)

About the "enhanced" functionalities, would it be possible to create a package such as JuMPExtensions with this kind of code? (Plus things like sensitivity analysis, as CPLEX does https://www.ibm.com/support/knowledgecenter/SSSA5P_12.7.0/ilog.odms.cplex.help/CPLEX/GettingStarted/topics/tutorials/InteractiveOptimizer/sensitivityAnalysis.html, possibly.)

[joehuchette]

I think we should focus with the MathProgBase interface portion first. I don't see a compelling reason to query only part of an IIS (correct me if wrong), so I think having only a single function that computes an IIS for a model would be appropriate. It could be called getIIS, and returns vectors of Bools that indicate whether a constraint, variable LB, variable UB, etc. are in the IIS.

E.g. MathProgBase.getIIS(m::Model) --> (constraintindices, varLBindices, varUBindices).

Does this sound reasonable?

[dourouc05]

That indeed sounds reasonable. The only problem I have with the interface is when this getIIS function is getting used outside linear constraints, but I don't think current solvers can do anything beyond (MI)LP concerning IIS. That could be made future-proof by returning a structure holding vectors of Booleans, to which we may add new fields for new kinds of constraints.

Another problem I have is with the name of the function, as it seems to imply that the results can be obtained quickly (like all existing get functions: it's just reading in memory something that's computed), as opposed to solve calls, which may take a while before returning.

[joehuchette]

I think it's alright to tailor the interface to LP in this case. We can also go with computeIIS if you prefer.

[mlubin]

Gurobi's IIS also includes SOS and quadratic constraints.

[dourouc05]

Hence, what would be the best? Returning a structure that may be extended in the future (or even a dictionary, indexed by the type of things in the IIS: variable bounds, linear/quadratic constraint, etc.; but that seems horrible to me)? That structure would already be concrete within MPB.

[mlubin]

A concrete type seems reasonable to me, e.g.,

type IISInfo
    varlb::Vector{Bool}
    varub::Vector{Bool}
    linconstr::Vector{Bool}
    quadconstr::Vector{Bool}
    sosconstr::Vector{Bool}
end

(don't take the names seriously, it's just a sketch)
This seems to call for an attributes concept in MPB, but we don't currently have that. If we're working at this low level, a vector of bools is probably closer to the solver representation than a list of indices. At the JuMP level we could choose to return a list of indices.

[dourouc05]

I agree that the vector of Booleans fits probably much better the level of abstraction of MPB; then, for JuMP, maybe both the vector of Booleans and the vector of indices make sense (they are just a find away from each other).

(I'm just not getting your point with attributes in MPB, but that's most likely due to my lack of acquaintance with that code.)

OK, so I'm planning to have drafts of this (first MPB, then solvers, then JuMP) on GitHub rather soonish so we may speak about actual code. (I think I should be targetting Julia 0.6 for this? That would just mean struct instead of type, if I followed things close enough.)

[mlubin]

We haven't dropped support for 0.5 at this point, so we would not be able to merge a PR that requires 0.6.

[dourouc05]

OK, no problem, I'll keep that in mind (but I've completely skipped version 0.5, migrating directly from 0.4 to 0.6).

[mlubin]

FWIW I expect we'll end up dropping 0.5 sooner rather than later.

[mlubin]

@dourouc05, following some very recent discussions, we're about to propose some very big changes to MathProgBase, so you may want to put your implementation on pause. More will come next week.

[dourouc05]

OK, no problem! (And I also see that I can use Julia 0.6 more freely!)

[dourouc05]

I guess I should wait for the branch https://github.com/JuliaOpt/MathProgBase.jl/tree/break_everything to be merged before considering working on this?

[mlubin]

Yep, that branch and corresponding updates to JuMP will make exposing the IIS much more straightforward (as a variable/constraint attribute).

[dourouc05]

OK, thank you! Do you have an ETA for this?

[mlubin]

Mid-July

[bbrunaud]

Interested in the discussion. I am using the following code for MILPs with Gurobi. I don't remeber who posted it, but it has been useful so far. Maybe it can be generalized.

function getIIS(m::JuMP.Model)
    grb_model = m.internalModel.inner
    num_constrs = Gurobi.num_constrs(grb_model)
    Gurobi.computeIIS(grb_model)
    iis_constrs = Gurobi.get_intattrarray(grb_model, "IISConstr",  1, num_constrs)
    m.linconstr[find(iis_constrs)]
end

[dourouc05]

Actually, the goal is to generalise such code!

[mlubin]

We now have a good solution for this issue with MOI, so let's actually put the pieces together so that we can consider the issue resolved. Gurobi needs to define an MOI attribute for IIS, and we need to make sure that it gets passed cleanly through JuMP.

[dourouc05]

To be sure of the steps to follow, here is what I understand now (I followed MOI's development, but only at a very high level), with the idea of a function that indicates, for each constraint, if it participates in the IIS:

• Define the needed attributes in MOI (around https://github.com/JuliaOpt/MathOptInterface.jl/blob/master/src/attributes.jl#L582)
• For each solver (say, Gurobi at first), implement supports and get for these attributes
• The first call to get must ask Gurobi to compute the IIS (I am thinking of adding a new attribute in https://github.com/JuliaOpt/Gurobi.jl/blob/2ba66a47a19050cca16f719eb1735a5620bb535a/src/grb_model.jl#L9 to indicate this fact. I have not yet understood how to reset this flag when the model is modified
• In JuMP, it would be nice to have a function to export the IIS as a new model (which the user could print or export as LP, for instance, using MathOptFormat)

@mlubin : any thoughts?

[dourouc05]

@mlubin up?

[odow]

The first step is probably to implement in in Gurobi.jl. Just an attribute in Gurobi.jl instead of MOI so you can go

model = JuMP.Model(with_optimizer(Gurobi.Optimizer))
@variable(model, x >= 1)
@constraint(model, c, 2x <= 0)
optimize!(model)
iis = MOI.get(model, Gurobi.IIS())
# typeof(iis) ???

Once we play around with Gurobi's IIS and sort out the API etc. We could add it to MOI.

[dourouc05]

Thanks @odow for your answer! (I may go first with CPLEX, as I don't have a Gurobi license on this computer, though, but I guess that's OK.)

[mlubin]

See https://github.com/blegat/JuMP-dev_2019_tutorial/blob/master/Defining%20new%20attributes.ipynb for how to define a new IIS attribute. (Steps 1 and 2 above.) After this is in place, we can figure out how to provide nice interfaces at the JuMP level.

[dourouc05]

While I'm working on this, I though of another layer of complexity: what to do for solvers which have no IIS functionality (like Cbc)? It would be nice to provide a pure Julia implementation of an IIS solver, that could be used to compute the IIS if it is not available. While providing this IIS solver is completely out of scope of what I'm doing right now, I think it would be best if the interface could fit this use case (i.e. MOI attributes not provided by the solver, but by an external layer on top of the MILP solver).

There could be other uses of wrappers to provide missing functionalities, but probably more edge cases than normal use. I thought of finding an optimal basic solution from an optimal solution to a LP (in case it is not solved by the simplex algorithm), but that's pretty much it…

[odow]

MOI attributes not provided by the solver, but by an external layer on top of the MILP solver

The attribute should be provided by the solver. You're half there on an external layer, but it would look something like:

model = Cbc.Optimizer()
# infeasible definition
MOI.get(model, MOI.IIS())  # errors

iis_model = IISPackage.Optimizer(model)
optimize!(model)
MOI.get(iis_model, MOI.IIS())

For now, we should just focus on solvers with a native IIS.

[dourouc05]

It perfectly makes sense, thanks!

Another thing on my mind: CPLEX has a conflict refiner, which is more generic than IIS. I guess the name of the attributes should thus not refer to IIS, to keep as much generality as possible?

[dourouc05]

For the record, the modifications for CPLEX are done: https://github.com/JuliaOpt/CPLEX.jl/pull/233. There is still one question left for now (how to report that constraints have not been proved to be in/out the conflict/IIS? --- there is no such concept for Gurobi, AFAIK).

Is this good enough to be ported to MOI?

[odow]

We should repeat the process for Gurobi to check that things are similar. There isn't a rush to port it to MOI.

[dourouc05]

No problem, except I have no Gurobi license on any computer I have access to right now. However, I may have access to an Xpress license at work, and it seems to have the feature (https://www.fico.com/fico-xpress-optimization/docs/latest/evalguide2/dhtml/eg2sec4_sec_eg2ssec42.html).

[dourouc05]

Process repeated for Xpress!

[dourouc05]

And also for Gurobi! (Once the first one is done, it's becoming easier, it seems!)

[dourouc05]

PRs are now merged for both CPLEX and Gurobi (Xpress still waiting). What would be the next step? Moving the new attributes to MOI?

[mlubin]

Yes, please write up a short proposal for the new attributes in the form of MOI issue. I'd be interested in the lessons learned from the implementations for Gurobi, CPLEX, and Xpress.

There is still one question left for now (how to report that constraints have not been proved to be in/out the conflict/IIS? --- there is no such concept for Gurobi, AFAIK).

Is this question resolved?

Thanks for your work on this!

[dourouc05]

This is kind of solved: either the IIS solver is done and the question makes no sense (everything is proved to be part of the IIS or not), or it is not and no distinction is made between things that are proved and those that have not been excluded yet. From a user perspective, I think this is OK, albeit it might be nice to have a specific way to query the exact status (but the attribute would then be a union of types: either a boolean or a missing value -- the latter indicating that the constraint has not been proved to be part of the IIS, which is closest to the meaning of missing, I think).

https://github.com/JuliaOpt/CPLEX.jl/pull/233#discussion_r278754808

There is still a question that was only partially decided: what about forcing the user to call compute_conflict each time an IIS is needed? The other option is hiding this function, and calling it under the hood when required; the current IIS would then be dropped each time the model is modified.

https://github.com/JuliaOpt/CPLEX.jl/pull/233#issuecomment-484262378

About the current limitations, there are mostly solver-specific things:

• for CPLEX, only LP (including integer) models are implemented. refineconflictext would need to be used for more than LPs: indicator constraints (are they implemented in CPLEX.jl?), quadratic constraints, special ordered set, piecewise-linear functions.
• for CPLEX, the notion of group is not handled (cannot tell a constraint to be part of the IIS: decisions are taken at the level of the group); this would require a binding to refineconflictext
• for Xpress, only one IIS is accessible
• for Xpress, no isolation is available (it seems to be a stronger property than "part of the IIS"): supporting it in MOI would mean than the attribute would have more than two/three possible values, so that would mean another enumeration to support them all (is there any nonneligible performance penalty for an enumeration vs. a boolean?)

For the non-solver-specific things (currently, all three implementation work for MILP models), requiring no change in MOI:

• for Gurobi, the attribute is only implemented for linear constraints, but implementing it for other kinds of constraints (SOS, quadratic) is just a matter of a few lines of code to add
• for Xpress, nothing is done for SOS1 and SOS2 constraints (I don't think model having quadratic constraints can be used for IIS computations)

[ghost]

How do I correctly apply MOI.get on constraints created via JuMP to identify those equations within the IIS? This is my code that throws an error:

model = Model(with_optimizer(Gurobi.Optimizer))
@variable(model, x >= 0)
@variable(model, y >= 0)

@constraint(model, con1, x + y <= 3)
@constraint(model, con2, x + y >= 5)
@objective(model, Min, x)
optimize!(model)

grb_model = model.moi_backend.optimizer.model
computeIIS(grb_model.inner)
MOI.get(grb_model, Gurobi.ConstraintConflictStatus(), con2.index)

(also see my post regarding this in the Julia board)

[blegat]

con2.index corresponds to thte index of the cache, which may not match the index of the optimizer, you should do

MOI.get(model, Gurobi.ConstraintConflictStatus(), con2)

[mlubin]

I'm surprised that direct_model isn't required for this. Has anyone confirmed that it works with Model(with_optimizer(Gurobi.Optimizer))?

[ghost]

I got my minimal example working like this now:

Gurobi.compute_conflict(model.moi_backend.optimizer.model)
MOI.get(model.moi_backend, Gurobi.ConstraintConflictStatus(), con1.index)

[blegat]

What's the error for MOI.get(model, Gurobi.ConstraintConflictStatus(), con2) ?

[ghost]

It works as well. I came up with the posted code, because your suggestion didn't work for me first, but I couldn't reproduce that error.

[dourouc05]

Follow-up on the topic. As the required bits have been accepted into MOI age should make it in release 0.9.14, I updated the support for the solvers:
https://github.com/JuliaOpt/Gurobi.jl/pull/326
https://github.com/JuliaOpt/CPLEX.jl/pull/299

I guess the next and final step would be to add something in JuMP to retrieve the IIS, I was thinking of a way to copy only the parts of the model that are in the IIS.

[dourouc05]

Now that everything needed is merged in MOI and that the PR for two solvers seem to have reached an acceptable state, I guess the time is right to discuss about what we could need in JuMP.

A simple thing would be to provide a function in JuMP to call MOI.compute_conflict! (just to have something easier to work with from the user point of view, exactly like optimize!), then the attributes defined in MOI are directly accessible. In other words, the user can directly determine whether a constraint is part of the conflict or not.

Then, the next thing is probably to allow copying only the relevant parts of the conflict into a new model. It would almost be a clone of copy_model, but only for the things that belong to the conflict (copy_conflict?). With this, you directly get a way to export only the conflict part of the model in an LP file (we could also copy Gurobi's behaviour: when exporting a model into a file with the .ilp extension, only consider the conflict).

I had a new look at other optimisation packages, but they don't seem to offer an interface for conflicts (Pyomo, GAMS, AIMMS, AMPL): for some solvers, it looks like you have access to the underlying functionality, but without any kind of unification. While it's a very good point for selling JuMP once it's done, there are not many sources for inspiration…

Do you have any other ideas, or comments about these ideas?

[odow]

What's the workflow for people using the IIS? I haven't really used it, so I don't have a strong opinion on what to do.

I imagine people would need to write out the IIS to an LP file for inspection, so having a way to copy_to seems like a good idea. That could happen at the MOI level though.

If there is no generic interface in other modeling packages, it might be a good idea to post on Discourse and say "Have you used the IIS feature of a solver? What is your workflow?/How would you design the IIS interface for JuMP?"

[dourouc05]

Whenever I needed that functionality, that was exactly my workflow. It worked okayish, but I couldn't find any better way to do it…

Here is the discussion: https://discourse.julialang.org/t/have-you-used-the-iis-conflict-refiner-feature-of-a-solver-how-would-you-like-to-access-it-from-jump/41862.

On the side, I'll also work on a copy_conflict_to function in MOI.

[dourouc05]

I've been digging in MOI.copy_to… and it's much more complicated than expected! Only pass_constraints would need to be modified to focus on constraints participating in the conflict, so I was thinking about adding one (keyword?) parameter to that function (so that it would be quite easy to use the new functionality without having to set all the existing optional positional arguments). That would mean adding this parameter to quite a few functions: default_copy_to (supports_default_copy_to?), automatic_copy_to, allocate_load (supports_allocate_load?), copy_to. copy_conflict_to could then also be overridden in other models, if need be.

What do you think of this? Or do you have other/better ideas?

[dourouc05]

I've just thought of another solution: add (more) parameters to copy_to and the others that take lambdas as arguments and are used to filter the elements to copy (e.g., one lambda for the variables, one for the constraints, one for the attributes, one for the NLP part). Before an element is copied, the relevant lambda is called: if it returns true, the element is copied (exactly like Base.filter). This would give the user a lot of flexibility, but with potential performance penalty (or is the Julia compiler able to generate a version of a function for the default value of a parameter, i.e. a lambda that always returns true?).

[dourouc05]

@odow @mlubin @joehuchette do you have opinions on how to implement this? As there has been no answer on Discourse, I guess that no other functionality is really required…

[mlubin]

A simple thing would be to provide a function in JuMP to call MOI.compute_conflict! (just to have something easier to work with from the user point of view, exactly like optimize!), then the attributes defined in MOI are directly accessible. In other words, the user can directly determine whether a constraint is part of the conflict or not.

Let's go ahead with this as a first step, it will be a big improvement over not having any JuMP-level API for IIS. As you observed, creating a sub-model with only the constraints that participate in the IIS is a bit complex, and it's not obvious what users want.

[mlubin]

2300 addresses the bulk of this issue, so I'm going to close this. New issues may be opened to discuss improvements to the IIS API. A big thanks to @dourouc05 for seeing this through over three years!

[dourouc05]

I was waiting for MOI 0.9.15 to be tagged (for https://github.com/jump-dev/MathOptInterface.jl/pull/1135), so that JuMP's interface is on the same level as Gurobi's API, for instance, but that's the last step I see for this :)!