search

jump-dev / Gurobi.jl

A Julia interface to the Gurobi Optimizer
http://www.gurobi.com/
MIT License
219 stars 80 forks source link

Gurobi sometimes fails to return duals in a QCP #415

Closed odow closed 1 year ago

odow commented 3 years ago

@adow031 has seen an example where Gurobi solves a QCP, experiences numerical error and fails to return the duals, but DualStatus is FEASIBLE_POINT.

The error is 

Barrier performed 54 iterations in 0.08 seconds
Sub-optimal termination - objective 2.48377287e+05

Warning: failed to compute QCP dual solution due to inaccurate barrier solution
         Try decreasing BarQCPConvTol for more accuracy

(This issue is to remind me that the problem exists. I'll try to find a reproducible example that isn't "run this large SDDP problem". @adow031 will post a .mof.json file when he next triggers it.)

I guess it isn't sufficient to just consider the SolCount. We should probably check that we can query RC. https://github.com/jump-dev/Gurobi.jl/blob/32b15ed6af8c32ff8b85436d9e3fec947f6d559a/src/MOI_wrapper/MOI_wrapper.jl#L2833-L2836

odow commented 2 years ago

@haoxiangyang89 found the likely culprit in https://github.com/odow/SDDP.jl/issues/477: the dual isn't defined at the point of a SecondOrderCone. I'll see if I can reproduce.

@adow031 / @drandyphilpott you should add a constraint to avoid the case 0 >= 0 in a quadratic constraint.

adow031 commented 2 years ago

@haoxiangyang89 found the likely culprit in odow/SDDP.jl#477: the dual isn't defined at the point of a SecondOrderCone. I'll see if I can reproduce.

@adow031 / @DrAndyPhilpott you should add a constraint to avoid the case 0 >= 0 in a quadratic constraint.

We found that running our model with CPLEX didn't have any errors; so it may be a Gurobi specific issue. However, running SDDP.jl with a QP took a very long time, so we have gone back to using a piecewise-linear approximation.

MathiasPBerger commented 2 years ago

Hi guys,

I don't know if this is still relevant but I just ran into the same problem on a fairly simple QCQP example, please see here.

The example was run under Julia 1.6.5 and uses JuMP v0.23.2 and Gurobi v0.11.0. Termination status is 'OPTIMAL' and dual status is 'FEASIBLE POINT' but the 'Pi' attribute cannot be retrieved.

odow commented 2 years ago

Oh great! I've been struggling to reproduce this.

Can you try adding constraints such that

(p_max[g] - p[g] + sum(alpha[g, i]*mu[i] for i = 1:n_w)) >= 0.001

My hunch is that the problem occurs at the point of the second-order cone.

MathiasPBerger commented 2 years ago

Will do this and let you know.

MathiasPBerger commented 2 years ago

The error indeed disappears when this constraint is added for each generator (g index). The optimal objectives are almost equal in both cases (with and without the constraint that you proposed).

I'm a bit perplexed by the fact that the matrix I'm using to build the second-order cone constraint is positive definite and the vector that multiplies it has nonzero entries in the solution (including in the case where the error occurs), so I don't exactly see how we could end up with a 0 >= 0 situation.

odow commented 2 years ago

If you don't include the constraints, what is value.(max_prod) and value.(min_prod) of the optimal solution?

MathiasPBerger commented 2 years ago

This is what I get:

max_prod: [50.27933145432515, 1.0638298606949077, 1.0638298606977323] [68.25933392790412, 0.7978734083340456, 0.7978734083360246] [79.04733539688482, 0.638296730977071, 0.6382967309721246]

min_prod: [71.31205574336671, 1.063829860695034, 1.0638298606982575] [53.33205327073633, 0.7978734083327054, 0.7978734083346006] [42.54405178742139, 0.6382967309756538, 0.6382967309741281]

I solved the same problem with Mosek and no error occurred, which supports the idea that this is a Gurobi-related issue.

odow commented 2 years ago

This is what I get

Interesting. We obviously need to dig a bit deeper to understand what's going on. I don't really have time to do this, unfortunately. cc @simonbowly

I solved the same problem with Mosek and no error occurred, which supports the idea that this is Gurobi-related issue

The problem is that Gurobi doesn't actually solve a second-order cone, it solves t^2 >= x^2 + y^2.

MathiasPBerger commented 2 years ago

Hi guys,

In case somebody picks it up later, I updated the formulation a bit and turned the original SOC constraints into quadratic constraints of the form x^T A x <= y^2 (along with a couple of affine equality constraints), where y is a non-negative scalar variable and A is positive definite.

Gurobi manages to solve the updated formulation without any problem, including for parameter values that proved problematic for the original SOCP formulation. The objectives and primal solutions are the same for both formulations.

Best,

Mathias