Quadratic constraint accuracy
AnsweredThe following 2D python model finds the greatest ball inscribed inside an ellipsoid; optimize retrieves the radius within 1.e8 accuracy however the optimal centre should be [0,0] but barrier algorithm yields
x 0.133747
y 1.09291e11
radius 0.707107
which is "far" from 0. Is there anyway to get a more accurate solution?
#!/usr/bin/env python3.7
# Copyright 2021, Gurobi Optimization, LLC
# This example formulates and solves the following simple model:
# maximize radius
# subject to x^2 + 2 y^2 <= 1
# x + radius <= 1
# x  radius >= 1
# y + radius <= sqrt(2)/2
# y  radius >= sqrt(2)/2
# radius nonnegative
import gurobipy as gp
from gurobipy import GRB
# Create a new model
m = gp.Model("inscribed")
# Create variables
x = m.addVar(name="x")
y = m.addVar(name="y")
radius = m.addVar(name="radius")
# Set objective: maximize radius
m.setObjective(radius, GRB.MAXIMIZE)
from scipy import sqrt
bd=sqrt(2)/2
# Add ellipsoid bounding box
m.addConstr(x + radius <= 1, "ellipsoid_UBx ")
m.addConstr(x  radius >= 1, "ellipsoid_LBx ")
m.addConstr(y + radius <= bd, "ellipsoid_UBy ")
m.addConstr(y  radius <= bd, "ellipsoid_UBy ")
# Add ellipsoid interior
m.addConstr(x**2 + 2*y**2 <= 1, "ellipsoid")
# First optimize() call will fail  need to set NonConvex to 2
try:
m.optimize()
except gp.GurobiError:
print("Optimize failed ")
for v in m.getVars():
print('%s %lg' % (v.varName, v.x))

Hi Dominique,
In your code, the line
m.addConstr(y  radius <= bd, "ellipsoid_UBy ")
should read
m.addConstr(y  radius >= bd, "ellipsoid_UBy ")
Moreover, the default lower bound for optimization variables is \(0\). Once you make your coordinate variables \(x,y\) free variables by setting
x = m.addVar(lb=GRB.INFINITY,name="x")
y = m.addVar(lb=GRB.INFINITY,name="y")you get the desired solution.
Best regards,
Jaromił
Please sign in to leave a comment.
Comments
1 comment