Getting proper ParameterDuals

[sstroemer]

I'm once again (after the close of the related issue) trying out POI, and stumbled across some questions related to "duals" of parameters. Four test cases - based on the README - below, all using:

using JuMP, HiGHS
import ParametricOptInterface as POI

I tried digging into the code but didn't really find a possible culprit (or maybe I'm just doing something wrong?). I assume all those test cases to be linked to the same underlying issue.

Edit: There is a small typo in duals.jl, on the example of getting duals. That should be MOI.get(model, POI.ParameterDual(), p).

TC1

model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
@variable(model, x)
@variable(model, p in POI.Parameter(1.0))
@constraint(model, cons, x * p >= 3)
@objective(model, Min, 2x)
optimize!(model)
MOI.get(model, POI.ParameterDual(), p) # 0.0

I assume this could/should at least warn the user (that the returned value of 0.0 is not a "dual"), throw an error, or return something else? (in this special case (a single variable in the LHS) there could actually be a reasonable value (depending on the sign of p though...; see also cache_multiplicative_params, although I assume the comment there talks about p * q instead))

TC2

model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
@variable(model, x)
@variable(model, p in POI.Parameter(1.0))
@variable(model, q); fix(q, 1.0)
@constraint(model, cons1, x >= p)
@constraint(model, cons2, x >= q)
@objective(model, Min, 2x)
optimize!(model)

MOI.get(model, POI.ParameterDual(), p) # 2.0
reduced_cost(q) # 0.0

Did I overlook something or shouldn't those two return the same value? This by the way does not happen with @variable(model, x >= 1.0) instead of cons2. (I know that depends on whether we think of an infinitesimal change as increasing or decreasing p but...)

TC3

model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
@variable(model, x)
@variable(model, p in POI.Parameter(1.0))
@constraint(model, cons, x >= 3 * p)
@objective(model, Min, 2x + p)
optimize!(model)

MOI.get(model, POI.ParameterDual(), p) # 5.0

Result should be 2.0 * 3.0 + 1.0 = 7.0 where the 2.0 comes from the objective coefficient of x, the 3.0 from the coefficient of p in the RHS, and the 1.0 from p occurring in the objective.

Remark: Re-reading everything after I typed it, I noticed that TC3 and TC4 are probably the same (different signs in RHS vs. objective)... so please ignore the overhead.

TC4

model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
@variable(model, x)
@variable(model, p in POI.Parameter(1.0))
@constraint(model, cons, x >= 2)
@objective(model, Min, 2x + 2*p)
optimize!(model)

MOI.get(model, POI.ParameterDual(), p) # -2.0

I'm not sure what the actual sign convention is on ParameterDuals, but I assume this should (?) have the same sign as the baseline given below.

# baseline TC4
model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
@variable(model, x)
@variable(model, p in POI.Parameter(1.0))
@constraint(model, cons, x >= p)
@objective(model, Min, 2x)
optimize!(model)

MOI.get(model, POI.ParameterDual(), p) # 2.0

I assume the calculation also causes the following result, where the RHS and the objective contribution therefore "cancel out" eachother (instead of summing up to 4.0, or -4.0).

# ad TC4
model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
@variable(model, x)
@variable(model, p in POI.Parameter(1.0))
@constraint(model, cons, x >= p)
@objective(model, Min, 2x + 2*p)
optimize!(model)

MOI.get(model, POI.ParameterDual(), p) # 0.0

---

edit1: changed TC2 since I messed it up when cleaning it up.

[joaquimg]

TC1 was a bug that was solved by #138

TC2 is not a bug...
You can look at HiGHS solution:

julia> model.moi_backend.optimizer.model.optimizer.solution.rowdual
2-element Vector{Float64}:
  2.0
 -0.0

both constraints are identical. The problem is dual degenerate. Therefore, the dual solution from HiGHS (and from POI) is valid.

TC3 and TC4 have the same underlying issue, also a bug, in which the difference between min and max was not treated correctly. Also solved in #138