Add support to query the reduced cost of a variable

[mtanneau]

Querying the reduced cost of a variable is very common in linear programming (and the corresponding function is in all linear solvers' APIs), however there is no straightforward way of doing this through MOI.
Currently, to access the reduced cost of a variable x, one has to manually query the dual value of every possible bound constraint on x (if any), and aggregate the results.
That's very cumbersome, given that this feature is supported by all linear solvers.

Proposal:
Add a variable attribute to query the reduced cost of a variable, e.g. ReducedCost, or (more verbose) VariableReducedCost. For now, this would only apply to linear programming solvers (see list below).

The syntax would look like:

MOI.get(model, MOI.ReducedCost(), x)

Solvers wrappers to be updated:

• Clp.jl
• CPLEX.jl
• GLPK.jl
• Gurobi.jl
• MosekTools.jl
• SCIP.jl
• Xpress.jl

If it makes sense, support could be extended to non-linear solvers as well.

[mlubin]

This was discussed early in the design of MOI. I am not in favor of a ReducedCost attribute because it's redundant with ConstraintDual on the bound constraints. Another attribute means additional code on all solver interfaces and another set of tests to maintain. Furthermore solvers like Mosek, ECOS, SCS, OSQP, and Ipopt naturally return separate duals for the upper and lower bound constraints.

I would be in favor of having a helper function in MOI.Utilities that does the "query the dual value of every possible bound constraint on x (if any), and aggregate the results" because I agree that it's painful for to have to roll one's own version of this function.

[blegat]

The utility in MOI.Utilities could be get_fallback(::ModelLike, ::MOI.ReducedCost, ::MOI.VariableRef) allowing solvers to choose between the fallback or calling the solver API. That would require us to add it to MOI.Test though.

[mtanneau]

Having a universal fallback + a way to call the solver's API sounds like a good compromise. I also understand that adding an attribute and tests means more code maintenance.

One more thing I'd like to mention here (more for future readers of this issue).
As Miles pointed out, reduced costs only make sense for LP or when using a simplex-like algorithm (whenever there's some complementary slackness condition actually). I am not aware of the reduced cost concept being used in non-linear programming.
Therefore, if one were to use reduced costs, I would expect they do so in an LP context, for which they can use Clp, GLPK, or any commercial LP solver. If one really needs to use, say, SCS or ECOS and access reduced costs, then, IMO, it should be implemented at the solver level first, rather than MOI.
Expecting this to be supported by any and all solvers sounds (to me) like putting an unnecessary toll on MOI.

[mtanneau]

If the MOIU fallback + possibility to call the solver API directly (if one wishes to) is acceptable to all, then I can volunteer a PR for this. Most likely some time this summer.
It's probably best to wait for 0.9 to be released first?

[mlubin]

A ReducedCost attribute with a get_fallback implementation is reasonable. That's a more uniform way to handle the issue than having a one-off function in MOIU that won't be transparently accessible from JuMP.

It's probably best to wait for 0.9 to be released first?

The release wouldn't block the get_fallback work, only the implementation in the solvers, so up to you when to start.

[mtanneau]

@chriscoey and I discussed it a bit further and, rather than ReducedCost (which is only defined for LP + simplex), the name VariableDual might make more sense. (would it?)

• For LP, VariableDual would refer to the traditional reduced cost of said variable.
• For conic problems, VariableDual of a variable x ∈ K would be the dual variable s ∈ K* of that conic constraint. This does raise the question of whether VariableDual for a single scalar variable would make sense when x is a vector. @chriscoey @juan-pablo-vielma what are your thoughts on this?
• I don't know how much sense it makes for general NLP ; I guess the default fallback would be to query the dual variable associated to the bound(s) constraint(s) on x.

The release wouldn't block the get_fallback work, only the implementation in the solvers, so up to you when to start.

Time is not a pressing matter here (at least not for me). I'd rather think this through than rush into things.

[blegat]

This is related to constrained variables: #710.
For a constrained variable, it's dual is the dual of the constrained index returned by add_constrained_variables. Otherwise, it is unclear if the dual is part of several SingleVariable and/or VectorOfVariables constraints

[mtanneau]

This is related to constrained variables: #710.
For a constrained variable, it's dual is the dual of the constrained index returned by add_constrained_variables. Otherwise, it is unclear if the dual is part of several SingleVariable and/or VectorOfVariables constraints

I'm not sure I fully understand what you meant. I could not find add_constrained_variables in the doc (is that in a PR that hasn't been merged yet?).

Math-wise, do we agree that, for a conic (vector) variable x ∈ K, the corresponding VariableDual would be the dual (vector) s ∈ K* w.r.t that conic constraint? I have not seen any specific name for s in Conic optimization literature, which is why I was asking for some conic people's opinion.

Code-wise, I agree, if VariableDual is defined for a conic variable, then there should be a non-ambiguous way to access that conic constraint. Because vectors are involved, that's likely to be a bit trickier.

For example, if I have a conic variable (x1, x2, x3) ∈ SOC, and denote (s1, s2, s3) ∈ SOC* the corresponding dual vector, then

MOI.get(model, MOI.VariableDual(), x1)

would return s1. Obviously, s1 does not make sense by itself, so one would have to extract all components of s before doing anything meaningful.

[blegat]

is that in a PR that hasn't been merged yet?

Yes, see #759

Math-wise, do we agree that, for a conic (vector) variable x ∈ K, the corresponding VariableDual would be the dual (vector) s ∈ K* w.r.t that conic constraint?

Yes but currently (i.e. without add_constrained_variables), there are no conic variables. Variables are created free and some constraints are added on them. With add_constrained_variables, the variables are conic and their conic constraint is returned.

[mtanneau]

Then it's probably best to wait until #759 is merged. Happy to keep the discussion open though.

[mlubin]

I'm unsure I understand the motivation for defining another way to access the duals on x ∈ K constraints. You can already access duals on conic variables by using ConstraintDual on the VectorOfVariables-in-Cone constraint. #759 returns both a collection of variables and the corresponding constraint index that lets you access duals.

This also seems quite far away from the original issue of "Why can't I access reduced costs easily?".

[henriquebecker91]

May I contribute to this in some way? I am using the following excerpt in my code:

function reduced_cost(var) :: Float64
  rc = 0.0 
  has_upper_bound(var) && (rc += shadow_price(UpperBoundRef(var)))
  has_lower_bound(var) && (rc += shadow_price(LowerBoundRef(var)))
  is_fixed(var) && (rc += shadow_price(FixRef(var)))
  !has_upper_bound(var) && !has_lower_bound(var) && !is_fixed(var) &&
    @warn "CAUTION: reduce_cost was called over a variable with no bounds"
  rc  
end

But I would prefer to have a 'official' method to call instead. While I need to use reduced_costs on my research I have not a deep understanding of simplex and the maths related to it. So I would have a greater peace of mind if I could call a method already named reduced_cost (or taking an attribute with this name) and that was examined by people that understand more about the topic.

[odow]

Agreed that rather than at the MOI level, this might make sense at the JuMP level. We have shadow_price, for example.

[henriquebecker91]

What you would suggest to change in my method to qualify it as a contribution? With some direction I can make a PR.

[odow]

Probably just the following plus docs and tests. There's a question whether we want to use shadow_price or dual. I would probably vote for shadow_price since this function is likely going to be used by people first encountering JuMP.

function reduced_cost(x::VariableRef)::Float64
    if !has_duals(owner_model(x))
        error("... something descriptive ...")
    end
    if is_fixed(x)
        return shadow_price(FixRef(x))
    end
    rc = 0.0
    if has_upper_bound(x)
        rc += shadow_price(UpperBoundRef(x)))
    end
    if has_lower_bound(x)
        rc += shadow_price(LowerBoundRef(x))
    end
    return rc  
end