Quadratically constrained MIQP with julia and gurobi

Question

This is an attempt to answer the following question: https://matheducators.stackexchange.com/questions/11757/small-data-sets-with-integral-sample-standard-deviations

So the intent of the following code is to find examples of small datasets with integer standard deviation. That can be formulated as a quadratically constrained mixed integer quadratic program, so I try to use Gurobin from Julia. Following is my code:

using JuMP  
using Gurobi

m = Model(solver = GurobiSolver() )
@variable(m,  0<= x[1:20] <= 100,  Int)
@variable(m,  Gj,  Int)
@constraint(m,  Gj == sum(x[1:20])/20 )
@variable(m,  Var,  Int)
@constraint(m,  Var == sum( (x[1:20]-Gj).^2/19) )
@variable(m,  sd,  Int)
@constraint(m, sd * sd == Var)
### We need some restrictions to avoid all equal, < or zero,  solutions:
@constraint(m,  sd >= 5)
@objective(m, Min, sd)

print(m)

status = solve(m)

println("Objective value: ", getobjectivevalue(m) )

x = getvalue(x)

Running this results in:

ERROR: Gurobi.GurobiError(10021, "Quadratic equality constraints")
Stacktrace:
 [1] optimize(::Gurobi.Model) at /home/kjetil/.julia/v0.6/Gurobi/src/grb_solve.jl:7
 [2] optimize!(::Gurobi.GurobiMathProgModel) at /home/kjetil/.julia/v0.6/Gurobi/src/GurobiSolverInterface.jl:294
 [3] #solve#101(::Bool, ::Bool, ::Bool, ::Array{Any,1}, ::Function, ::JuMP.Model) at /home/kjetil/.julia/v0.6/JuMP/src/solvers.jl:173
 [4] solve(::JuMP.Model) at /home/kjetil/.julia/v0.6/JuMP/src/solvers.jl:148

Any ideas?

Answer 1

A math programming solver like Gurobi Optimizer cannot solve models with quadratic equality constraints. Here are the types of constraints that Gurobi Optimizer can solve. To solve your model using Gurobi Optimizer, you must transform your constraints into one of these forms, such as quadratic inequality constraints.

Answer 2

The major problem is that, in general, a quadratic equality is not convex, and most solvers only work for convex problems (plus integer constraints). A product of two binary variables is easy to linearise (it's the equivalent of a logical AND), that of one binary variable and one continuous variable is easy too; the rest is not so easy.

Since Gurobi 9, you can solve nonconvex bilinear problems, in particular those having quadratic equality constraints. You just have to add the right parameter. With Gurobi.jl, if m is your JuMP model, you can do this:

set_optimizer_attribute(m, "NonConvex", 2)