X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=dunshire%2Fcones.py;h=50c57e747952cd207005d82ef25b6817e29bef61;hb=779389fb623926afc411587fdc43f69bf63808a6;hp=7fc8bf94993e4df86b2f57c2e6b6d6aefa65d8c4;hpb=bdb596b84a06d0c97e39d42586a51fc36ba44186;p=dunshire.git diff --git a/dunshire/cones.py b/dunshire/cones.py index 7fc8bf9..50c57e7 100644 --- a/dunshire/cones.py +++ b/dunshire/cones.py @@ -77,6 +77,44 @@ class SymmetricCone: raise NotImplementedError + def ball_radius(self, point): + """ + Return the radius of a ball around ``point`` in this cone. + + Since a radius cannot be negative, the ``point`` must be + contained in this cone; otherwise, an error is raised. + + Parameters + ---------- + + point : matrix + A point contained in this cone. + + Returns + ------- + + float + A radius ``r`` such that the ball of radius ``r`` centered at + ``point`` is contained entirely within this cone. + + Raises + ------ + + NotImplementedError + Always, this method must be implemented in subclasses. + + Examples + -------- + + >>> K = SymmetricCone(5) + >>> K.ball_radius(matrix([1,1,1,1,1])) + Traceback (most recent call last): + ... + NotImplementedError + + """ + raise NotImplementedError + def dimension(self): """ @@ -201,6 +239,64 @@ class NonnegativeOrthant(SymmetricCone): return all([x > -options.ABS_TOL for x in point]) + def ball_radius(self, point): + """ + Return the radius of a ball around ``point`` in this cone. + + Since a radius cannot be negative, the ``point`` must be + contained in this cone; otherwise, an error is raised. + + The minimum distance from ``point`` to the complement of this + cone will always occur at its projection onto that set. It is + not hard to see that the projection is directly along one of the + coordinates, and so the minimum distance from ``point`` to the + boundary of this cone is the smallest coordinate of + ``point``. Therefore any radius less than that will work; we + divide it in half somewhat arbitrarily. + + Parameters + ---------- + + point : matrix + A point contained in this cone. + + Returns + ------- + + float + A radius ``r`` such that the ball of radius ``r`` centered at + ``point`` is contained entirely within this cone. + + Raises + ------ + + TypeError + If ``point`` is not a :class:`cvxopt.base.matrix`. + + TypeError + If ``point`` has the wrong dimensions. + + ValueError + if ``point`` is not contained in this cone. + + Examples + -------- + + >>> K = NonnegativeOrthant(5) + >>> K.ball_radius(matrix([1,2,3,4,5])) + 0.5 + + """ + if not isinstance(point, matrix): + raise TypeError('the given point is not a cvxopt.base.matrix') + if not point.size == (self.dimension(), 1): + raise TypeError('the given point has the wrong dimensions') + if not point in self: + raise ValueError('the given point does not lie in the cone') + + return min(list(point)) / 2.0 + + class IceCream(SymmetricCone): """ @@ -305,6 +401,81 @@ class IceCream(SymmetricCone): return norm(radius) < (height + options.ABS_TOL) + def ball_radius(self, point): + r""" + Return the radius of a ball around ``point`` in this cone. + + Since a radius cannot be negative, the ``point`` must be + contained in this cone; otherwise, an error is raised. + + The minimum distance from ``point`` to the complement of this + cone will always occur at its projection onto that set. It is + not hard to see that the projection is at a "down and out" angle + of :math:`\pi/4` towards the outside of the cone. If one draws a + right triangle involving the ``point`` and that projection, it + becomes clear that the distance between ``point`` and its + projection is a factor of :math:`1/\sqrt(2)` times the "horizontal" + distance from ``point`` to boundary of this cone. For simplicity + we take :math:`1/2` instead. + + Parameters + ---------- + + point : matrix + A point contained in this cone. + + Returns + ------- + + float + A radius ``r`` such that the ball of radius ``r`` centered at + ``point`` is contained entirely within this cone. + + Raises + ------ + + TypeError + If ``point`` is not a :class:`cvxopt.base.matrix`. + + TypeError + If ``point`` has the wrong dimensions. + + ValueError + if ``point`` is not contained in this cone. + + Examples + -------- + + The height of ``x`` below is one (its first coordinate), and so + the radius of the circle obtained from a height-one cross + section is also one. Note that the last two coordinates of ``x`` + are half of the way to the boundary of the cone, and in the + direction of a 30-60-90 triangle. If one follows those + coordinates, they hit at :math:`\left(1, \frac{\sqrt(3)}{2}, + \frac{1}{2}\right)` having unit norm. Thus the "horizontal" + distance to the boundary of the cone is :math:`1 - \left\lVert x + \right\rVert`, which simplifies to :math:`1/2`. And rather than + involve a square root, we divide by two for a final safe radius + of :math:`1/4`. + + >>> from math import sqrt + >>> K = IceCream(3) + >>> x = matrix([1, sqrt(3)/4.0, 1/4.0]) + >>> K.ball_radius(x) + 0.25 + + """ + if not isinstance(point, matrix): + raise TypeError('the given point is not a cvxopt.base.matrix') + if not point.size == (self.dimension(), 1): + raise TypeError('the given point has the wrong dimensions') + if not point in self: + raise ValueError('the given point does not lie in the cone') + + height = point[0] + radius = norm(point[1:]) + return (height - radius) / 2.0 + class SymmetricPSD(SymmetricCone): r""" @@ -579,11 +750,10 @@ class CartesianProduct(SymmetricCone): def cvxopt_dims(self): """ - Return a dictionary of dimensions corresponding to the factors - of this cartesian product. The format of this dictionary is - described in the CVXOPT user's guide: - - http://cvxopt.org/userguide/coneprog.html#linear-cone-programs + Return a dictionary of dimensions corresponding to the + factors of this cartesian product. The format of this dictionary + is described in the `CVXOPT user's guide + `_. Returns -------