Update optimization model by changing variables to coefficients

Specific Domains Optimization (Mathematical)
1

Hi there, :slight_smile:
I am trying to extract shadow prices out of a MIQP model. In order to do this I want to “fix” the mixed integer variables after a first solving and then extract the duals.

Here’s a minimal working example:

using JuMP
using Gurobi

model = Model(Gurobi.Optimizer)

@variable(model, x, Bin)
@variable(model, y)

@constraint(model, x * y >= 1)
@constraint(model, c1, x + y <= 6)
@objective(model, Max, y)

optimize!(model)

solution = Dict(v => value(v) for v in all_variables(model))
for v in all_variables(model)
    is_binary(v) || continue
    unset_binary(v)
    fix(v, solution[v]; force=true)
end
new_model = copy(model)
set_optimizer(new_model, Gurobi.Optimizer)
set_optimizer_attribute(new_model, "QCPDual", 1)
optimize!(new_model)

JuMP.shadow_price(new_model[:c1])

However, for larger and more difficult models, QCPDual is not guaranteed to compute duals according to this post.

I was thinking if there is an option to get fully rid of the binary/integer variables and replace those by the result values in order to end up with a full LP. Is there a way to e.g. iterate through my constraints by using e.g. all_constraints(model; include_variable_in_set_constraints=true) replace the variables and afterwards delete the variables from the model?
In the end I would hope to get then a model in this case which looks like this:

max(y)
s.t.
1 * y >= 1
1 + y <= 6

So far I was not able to figure this out…

2

This is an interesting question. We don’t have a good way to solve this, unfortunately, and I can’t even recommend a simple work-around.

The fixed variable approach should be sufficient. Do you have an example where the QCP dual failed?

Also note that you can fix the variables with:

undo = fix_discrete_variables(model);

See Computing the duals of a mixed-integer program · JuMP

2 Likes
3

I don’t think this problem is well posed.

@max Dual variable (shadow prices) exists in convex conic programs, not in QP.
Because QP may be convex or nonconvex. Therefore, if you want a dual variable, you should first formulate your problem as a convex quadratic conic program, not only a QP.

If your intention is to get an LP with adjustable parameters, then just natively write a function in which you build this LP, and the function arguments being those parameters.

4

The problem is well posed, and Gurobi can compute duals for QPs (and QCPs, if you set QCPDual to 1).

Here’s how I would write the code in JuMP:

julia> using JuMP, Gurobi

julia> begin
           model = Model(Gurobi.Optimizer)
           set_silent(model)
           @variable(model, x, Bin)
           @variable(model, y)
           @constraint(model, x * y >= 1)
           @constraint(model, c1, x + y <= 6)
           @objective(model, Max, y)
           optimize!(model)
           undo = fix_discrete_variables(model)
           set_attribute(model, "QCPDual", 1)
           optimize!(model)
           @assert is_solved_and_feasible(model)
           shadow_price(c1)
       end
1.0
1 Like
5

I don’t think in a realistic operations research problem we will make use of the dual variable of a general QP (or QCQP).

We tend to write dual programs for convex conic programs (including convex SOCP). (Do we also write dual programs for QP? seems exists but I don’t encounter it often).

The Gurobi solver is usable but some times may have problems like Gurobi Error 10005: Unable to retrieve attribute 'Pi' - #7 by WalterMadelim

I have a bit of experience about the Gurobi solver. Sometimes I encounter strange behaviors. Perhaps tuning those parameters could help make it more stable, like Presolve, DualReductions

6

Do we also write dual programs for QP

See Duality · MathOptInterface

7

Thanks for your quick answers! Actually, I cannot provide an easy example where I encounter the failed computation of the Duals by QCPDual…

However, I thought it could be an “easy” way of getting rid of the QP formulations (at least in my usecase). I only have QP formulations in my model where I multiply binaries with linear continuous variables. That’s why I thought I could remove those formulations.

Maybe I need to consider reformulating my Model

8

You can reformulate x * y to:

using JuMP
function bilinear(model, x, y)
    @assert is_binary(x)
    @assert has_lower_bound(y)
    @assert has_upper_bound(y)
    y_L, y_U = lower_bound(y), upper_bound(y)
    xy = @variable(model)
    @constraint(model, xy <= y_U * x)
    @constraint(model, xy >= y_U * x)
    @constraint(model, xy <= y - y_L * (1 - x))
    @constraint(model, xy >= y - y_U * (1 - x))
    return xy
end

model = Model()
@variable(model, x, Bin)
@variable(model, -100 <= y <= 100)
xy = bilinear(model, x, y)
@constraint(model, xy >= 1)
1 Like