Ask Your Question

Find the center of a maximum inscribed sphere of a convex polytope

asked 2020-04-23 11:57:31 +0100

Sébastien gravatar image

updated 2020-04-23 13:31:17 +0100

Is there a way in SageMath to find the center of the (or of some, if not unique) sphere of maximum radius which is inside of a convex polytope? I ask the question in general dimension, but I am currently interested in a solution in dimension 2.

edit retag flag offensive close merge delete

1 Answer

Sort by » oldest newest most voted

answered 2020-04-23 22:52:44 +0100

Sébastien gravatar image

It's ok, I found the solution using MILP:

def center_inscribed_sphere_polytope(polytope, solver=None):
    Return the center and radius of maximal inscribed sphere


    - ``polytope`` -- polytope


    a 2-tuple (center, radius)


        sage: P = polytopes.associahedron(['A',3])
        sage: center_inscribed_sphere_polytope(P)
        ([0.03553390593273766, 0.5355339059327378, 0.03553390593273766],
        sage: center_inscribed_sphere_polytope(P + vector((10,10,10)))
        ([10.035533905932738, 10.535533905932738, 10.035533905932738],


        sage: P = Polyhedron([(0,0), (1,0), (0,10)])
        sage: center_inscribed_sphere_polytope(P)
        ([0.47506218943955486, 0.47506218943955486], 0.47506218943955486)

    from math import sqrt
    from sage.numerical.mip import MixedIntegerLinearProgram
    p = MixedIntegerLinearProgram(solver=solver)
    x = p.new_variable()
    d = p.new_variable()
    for ineq in polytope.inequalities_list():
        constant = ineq[0]
        coeffs = ineq[1:]
        norm = sqrt(sum(a**2 for a in coeffs))
        S = p.sum(a/norm*x[i] for (i,a) in enumerate(coeffs))
        p.add_constraint(constant/norm + S >= d[0])
    radius = p.solve()
    soln = p.get_values(x)
    center = [soln[k] for k in sorted(x.keys())]
    return center, radius
edit flag offensive delete link more


I don't know if the above can be changed to find the center of the John ellipsoid instead?

Sébastien gravatar imageSébastien ( 2020-04-29 09:47:25 +0100 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower


Asked: 2020-04-23 11:57:31 +0100

Seen: 286 times

Last updated: Apr 23 '20