X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=src%2Fdunshire%2Fcones.py;h=61feab453bd0866b54c053d2871d722a63ec3730;hb=d1f749bdc31dcaef1dce62f303c1889f40ec0fe2;hp=b60811e548de76d95abe0933a245b5a52141b606;hpb=56ea961887d507114174af5f92b8c3c77b0b7a50;p=dunshire.git diff --git a/src/dunshire/cones.py b/src/dunshire/cones.py index b60811e..61feab4 100644 --- a/src/dunshire/cones.py +++ b/src/dunshire/cones.py @@ -1,10 +1,11 @@ """ Class definitions for all of the symmetric cones (and their superclass, -SymmetricCone) supported by CVXOPT. +:class:`SymmetricCone`) supported by CVXOPT. """ from cvxopt import matrix -from matrices import norm +from matrices import eigenvalues, norm +import options class SymmetricCone: """ @@ -12,41 +13,70 @@ class SymmetricCone: There are three types of symmetric cones supported by CVXOPT: - * The nonnegative orthant in the real n-space. + 1. The nonnegative orthant in the real n-space. + 2. The Lorentz "ice cream" cone, or the second-order cone. + 3. The cone of symmetric positive-semidefinite matrices. - * The Lorentz "ice cream" cone, or the second-order cone. + This class is intended to encompass them all. - * The cone of symmetric positive-semidefinite matrices. + When constructing a single symmetric cone (i.e. not a + :class:`CartesianProduct` of them), the only information that we + need is its dimension. We take that dimension as a parameter, and + store it for later. + + Parameters + ---------- + + dimension : int + The dimension of this cone. + + Raises + ------ + + ValueError + If you try to create a cone with dimension zero or less. - This class is intended to encompass them all. """ def __init__(self, dimension): """ A generic constructor for symmetric cones. - - When constructing a single symmetric cone (i.e. not a cartesian - product of them), the only information that we need is its - dimension. We take that dimension as a parameter, and store it - for later. """ if dimension <= 0: raise ValueError('cones must have dimension greater than zero') self._dimension = dimension + def __contains__(self, point): """ Return whether or not ``point`` belongs to this cone. - """ - raise NotImplementedError - def contains_strict(self, point): - """ - Return whether or not ``point`` belongs to the interior - of this cone. + Parameters + ---------- + + point : matrix + The point to test for membership in this cone. + + Raises + ------ + + NotImplementedError + Always, this method must be implemented in subclasses. + + Examples + -------- + + >>> K = SymmetricCone(5) + >>> matrix([1,2]) in K + Traceback (most recent call last): + ... + NotImplementedError + """ raise NotImplementedError + + def dimension(self): """ Return the dimension of this symmetric cone. @@ -56,15 +86,31 @@ class SymmetricCone: should not need to be overridden in subclasses. We prefer to do any special computation in ``__init__()`` and record the result in ``self._dimension``. + + Returns + ------- + + int + The stored dimension (from when this cone was constructed) + of this cone. + + Examples + -------- + + >>> K = SymmetricCone(5) + >>> K.dimension() + 5 + """ return self._dimension class NonnegativeOrthant(SymmetricCone): """ - The nonnegative orthant in ``n`` dimensions. + The nonnegative orthant in the given number of dimensions. - EXAMPLES: + Examples + -------- >>> K = NonnegativeOrthant(3) >>> print(K) @@ -78,78 +124,69 @@ class NonnegativeOrthant(SymmetricCone): tpl = 'Nonnegative orthant in the real {:d}-space' return tpl.format(self.dimension()) + def __contains__(self, point): """ Return whether or not ``point`` belongs to this cone. - INPUT: - - An instance of the ``cvxopt.base.matrix`` class having - dimensions ``(n,1)`` where ``n`` is the dimension of this cone. - - EXAMPLES: + Since this test is expected to work on points whose components + are floating point numbers, it doesn't make any sense to + distinguish between strict and non-strict containment -- the + test uses a tolerance parameter. - >>> K = NonnegativeOrthant(3) - >>> matrix([1,2,3]) in K - True - - >>> K = NonnegativeOrthant(3) - >>> matrix([1,-0.1,3]) in K - False + Parameters + ---------- - >>> K = NonnegativeOrthant(3) - >>> [1,2,3] in K - Traceback (most recent call last): - ... - TypeError: the given point is not a cvxopt.base.matrix + point : matrix + A :class:`cvxopt.base.matrix` having dimensions ``(n,1)`` + where ``n`` is the :meth:`dimension` of this cone. - >>> K = NonnegativeOrthant(3) - >>> matrix([1,2]) in K - Traceback (most recent call last): - ... - TypeError: the given point has the wrong dimensions + Returns + ------- - """ - 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') + bool - return all([x >= 0 for x in point]) + ``True`` if ``point`` belongs to this cone, ``False`` otherwise. + Raises + ------ - def contains_strict(self, point): - """ - Return whether or not ``point`` belongs to the interior of this - cone. + TypeError + If ``point`` is not a :class:`cvxopt.base.matrix`. - INPUT: + TypeError + If ``point`` has the wrong dimensions. - An instance of the ``cvxopt.base.matrix`` class having - dimensions ``(n,1)`` where ``n`` is the dimension of this cone. + Examples + -------- - EXAMPLES: + All of these coordinates are positive enough: >>> K = NonnegativeOrthant(3) - >>> K.contains_strict(matrix([1,2,3])) + >>> matrix([1,2,3]) in K True + The one negative coordinate pushes this point outside of ``K``: + >>> K = NonnegativeOrthant(3) - >>> K.contains_strict(matrix([1,0,1])) + >>> matrix([1,-0.1,3]) in K False + A boundary point is considered inside of ``K``: >>> K = NonnegativeOrthant(3) - >>> K.contains_strict(matrix([1,-0.1,3])) - False + >>> matrix([1,0,3]) in K + True + + Junk arguments don't work: >>> K = NonnegativeOrthant(3) - >>> K.contains_strict([1,2,3]) + >>> [1,2,3] in K Traceback (most recent call last): ... TypeError: the given point is not a cvxopt.base.matrix >>> K = NonnegativeOrthant(3) - >>> K.contains_strict(matrix([1,2])) + >>> matrix([1,2]) in K Traceback (most recent call last): ... TypeError: the given point has the wrong dimensions @@ -160,15 +197,16 @@ class NonnegativeOrthant(SymmetricCone): if not point.size == (self.dimension(), 1): raise TypeError('the given point has the wrong dimensions') - return all([x > 0 for x in point]) + return all([x > -options.ABS_TOL for x in point]) class IceCream(SymmetricCone): """ - The nonnegative orthant in ``n`` dimensions. + The Lorentz "ice cream" cone in the given number of dimensions. - EXAMPLES: + Examples + -------- >>> K = IceCream(3) >>> print(K) @@ -187,25 +225,57 @@ class IceCream(SymmetricCone): """ Return whether or not ``point`` belongs to this cone. - INPUT: + Since this test is expected to work on points whose components + are floating point numbers, it doesn't make any sense to + distinguish between strict and non-strict containment -- the + test uses a tolerance parameter. + + Parameters + ---------- + + point : matrix + A :class:`cvxopt.base.matrix` having dimensions ``(n,1)`` + where ``n`` is the :meth:`dimension` of this cone. + + Returns + ------- + + bool + + ``True`` if ``point`` belongs to this cone, ``False`` otherwise. - An instance of the ``cvxopt.base.matrix`` class having - dimensions ``(n,1)`` where ``n`` is the dimension of this cone. + Raises + ------ - EXAMPLES: + TypeError + If ``point`` is not a :class:`cvxopt.base.matrix`. + + TypeError + If ``point`` has the wrong dimensions. + + Examples + -------- + + This point lies well within the ice cream cone: >>> K = IceCream(3) >>> matrix([1,0.5,0.5]) in K True + This one lies on its boundary: + >>> K = IceCream(3) >>> matrix([1,0,1]) in K True + This point lies entirely outside of the ice cream cone: + >>> K = IceCream(3) >>> matrix([1,1,1]) in K False + Junk arguments don't work: + >>> K = IceCream(3) >>> [1,2,3] in K Traceback (most recent call last): @@ -228,44 +298,140 @@ class IceCream(SymmetricCone): if self.dimension() == 1: # In one dimension, the ice cream cone is the nonnegative # orthant. - return height >= 0 + return height > -options.ABS_TOL else: radius = point[1:] - return height >= norm(radius) + return norm(radius) < (height + options.ABS_TOL) - def contains_strict(self, point): + +class SymmetricPSD(SymmetricCone): + r""" + The cone of real symmetric positive-semidefinite matrices. + + This cone has a dimension ``n`` associated with it, but we let ``n`` + refer to the dimension of the domain of our matrices and not the + dimension of the (much larger) space in which the matrices + themselves live. In other words, our ``n`` is the ``n`` that appears + in the usual notation :math:`S^{n}` for symmetric matrices. + + As a result, the cone ``SymmetricPSD(n)`` lives in a space of dimension + :math:`\left(n^{2} + n\right)/2)`. + + Examples + -------- + + >>> K = SymmetricPSD(3) + >>> print(K) + Cone of symmetric positive-semidefinite matrices on the real 3-space + >>> K.dimension() + 3 + + """ + def __str__(self): """ - Return whether or not ``point`` belongs to the interior - of this cone. + Output a human-readable description of myself. + """ + tpl = 'Cone of symmetric positive-semidefinite matrices ' \ + 'on the real {:d}-space' + return tpl.format(self.dimension()) + - INPUT: + def __contains__(self, point): + """ + Return whether or not ``point`` belongs to this cone. - An instance of the ``cvxopt.base.matrix`` class having - dimensions ``(n,1)`` where ``n`` is the dimension of this cone. + Since this test is expected to work on points whose components + are floating point numbers, it doesn't make any sense to + distinguish between strict and non-strict containment -- the + test uses a tolerance parameter. - EXAMPLES: + Parameters + ---------- - >>> K = IceCream(3) - >>> K.contains_strict(matrix([1,0.5,0.5])) + point : matrix + A :class:`cvxopt.base.matrix` having dimensions ``(n,n)`` + where ``n`` is the :meth:`dimension` of this cone. + + Returns + ------- + + bool + + ``True`` if ``point`` belongs to this cone, ``False`` otherwise. + + Raises + ------ + + TypeError + If ``point`` is not a :class:`cvxopt.base.matrix`. + + TypeError + If ``point`` has the wrong dimensions. + + Examples + -------- + + These all lie in the interior of the Symmetric PSD cone: + + >>> K = SymmetricPSD(2) + >>> matrix([[1,0],[0,1]]) in K True - >>> K = IceCream(3) - >>> K.contains_strict(matrix([1,0,1])) - False + >>> K = SymmetricPSD(3) + >>> matrix([[2,-1,0],[-1,2,-1],[0,-1,2]]) in K + True - >>> K = IceCream(3) - >>> K.contains_strict(matrix([1,1,1])) - False + >>> K = SymmetricPSD(5) + >>> A = matrix([[5,4,3,2,1], + ... [4,5,4,3,2], + ... [3,4,5,4,3], + ... [2,3,4,5,4], + ... [1,2,3,4,5]]) + >>> A in K + True - >>> K = IceCream(3) - >>> K.contains_strict([1,2,3]) + Boundary points lie in the cone as well: + + >>> K = SymmetricPSD(2) + >>> matrix([[0,0],[0,0]]) in K + True + + >>> K = SymmetricPSD(5) + >>> A = matrix([[1,0,0,0,0], + ... [0,1,0,0,0], + ... [0,0,0,0,0], + ... [0,0,0,1,0], + ... [0,0,0,0,1]]) + >>> A in K + True + + However, this matrix has a negative eigenvalue: + + >>> K = SymmetricPSD(2) + >>> A = matrix([[ 1, -2], + ... [-2, 1]]) + >>> A in K + False + + An asymmetric cone with positive eigenvalues is not in the cone: + + >>> K = SymmetricPSD(2) + >>> A = matrix([[10, 2], + ... [4, 8]]) + >>> A in K + False + + Junk arguments don't work: + + >>> K = SymmetricPSD(2) + >>> [[1,2],[2,3]] in K Traceback (most recent call last): ... TypeError: the given point is not a cvxopt.base.matrix - >>> K = IceCream(3) - >>> K.contains_strict(matrix([1,2])) + >>> K = SymmetricPSD(3) + >>> matrix([[1,2],[3,4]]) in K Traceback (most recent call last): ... TypeError: the given point has the wrong dimensions @@ -273,37 +439,15 @@ class IceCream(SymmetricCone): """ if not isinstance(point, matrix): raise TypeError('the given point is not a cvxopt.base.matrix') - if not point.size == (self.dimension(), 1): + if not point.size == (self.dimension(), self.dimension()): raise TypeError('the given point has the wrong dimensions') + if not point.typecode == 'd': + point = matrix(point, (self.dimension(), self.dimension()), 'd') + if not norm(point - point.trans()) < options.ABS_TOL: + # It's not symmetric. + return False + return all([e > -options.ABS_TOL for e in eigenvalues(point)]) - height = point[0] - if self.dimension() == 1: - # In one dimension, the ice cream cone is the nonnegative - # orthant. - return height > 0 - else: - radius = point[1:] - return height > norm(radius) - - -class SymmetricPSD(SymmetricCone): - """ - The nonnegative orthant in ``n`` dimensions. - - EXAMPLES: - - >>> K = SymmetricPSD(3) - >>> print(K) - Cone of symmetric positive-semidefinite matrices on the real 3-space - - """ - def __str__(self): - """ - Output a human-readable description of myself. - """ - tpl = 'Cone of symmetric positive-semidefinite matrices ' \ - 'on the real {:d}-space' - return tpl.format(self.dimension()) class CartesianProduct(SymmetricCone): @@ -311,7 +455,8 @@ class CartesianProduct(SymmetricCone): A cartesian product of symmetric cones, which is itself a symmetric cone. - EXAMPLES: + Examples + -------- >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(2)) >>> print(K) @@ -325,6 +470,7 @@ class CartesianProduct(SymmetricCone): super().__init__(my_dimension) self._factors = factors + def __str__(self): """ Output a human-readable description of myself. @@ -335,121 +481,77 @@ class CartesianProduct(SymmetricCone): format_args += list(self.factors()) return tpl.format(*format_args) + def __contains__(self, point): """ Return whether or not ``point`` belongs to this cone. - INPUT: + The ``point`` is expected to be a tuple of points which will be + tested for membership in this cone's factors. If each point in + the tuple belongs to its corresponding factor, then the whole + point belongs to this cone. Otherwise, it doesn't. - An instance of the ``cvxopt.base.matrix`` class having - dimensions ``(n,1)`` where ``n`` is the dimension of this cone. + Parameters + ---------- - EXAMPLES: + point : tuple of matrix + A tuple of :class:`cvxopt.base.matrix` corresponding to the + :meth:`factors` of this cartesian product. - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([1,2,3,1,0.5,0.5]) in K - True + Returns + ------- - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([0,0,0,1,0,1]) in K - True + bool - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([1,1,1,1,1,1]) in K - False + ``True`` if ``point`` belongs to this cone, ``False`` otherwise. - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([1,-1,1,1,0,1]) in K - False + Raises + ------ - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> [1,2,3,4,5,6] in K - Traceback (most recent call last): - ... - TypeError: the given point is not a cvxopt.base.matrix + TypeError + If ``point`` is not a tuple of :class:`cvxopt.base.matrix`. - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([1,2]) in K - Traceback (most recent call last): - ... - TypeError: the given point has the wrong dimensions + TypeError + If any element of ``point`` has the wrong dimensions. - """ - 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') - - for factor in self.factors(): - # Split off the components of ``point`` corresponding to - # ``factor``. - factor_part = point[0:factor.dimension()] - if not factor_part in factor: - return False - point = point[factor.dimension():] - - return True - - - def contains_strict(self, point): - """ - Return whether or not ``point`` belongs to the interior - of this cone. - - INPUT: - - An instance of the ``cvxopt.base.matrix`` class having - dimensions ``(n,1)`` where ``n`` is the dimension of this cone. + Examples + -------- - EXAMPLES: + The result depends on how containment is defined for our factors: >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> K.contains_strict(matrix([1,2,3,1,0.5,0.5])) + >>> (matrix([1,2,3]), matrix([1,0.5,0.5])) in K True >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> K.contains_strict(matrix([1,2,3,1,0,1])) - False + >>> (matrix([0,0,0]), matrix([1,0,1])) in K + True >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> K.contains_strict(matrix([0,1,1,1,0.5,0.5])) + >>> (matrix([1,1,1]), matrix([1,1,1])) in K False >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> K.contains_strict(matrix([1,1,1,1,1,1])) + >>> (matrix([1,-1,1]), matrix([1,0,1])) in K False - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> K.contains_strict(matrix([1,-1,1,1,0,1])) - False + Junk arguments don't work: >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> K.contains_strict([1,2,3,4,5,6]) + >>> [[1,2,3],[4,5,6]] in K Traceback (most recent call last): ... TypeError: the given point is not a cvxopt.base.matrix >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> K.contains_strict(matrix([1,2])) + >>> (matrix([1,2]), matrix([3,4,5,6])) in K Traceback (most recent call last): ... TypeError: the given point has the wrong dimensions """ - 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') - - for factor in self.factors(): - # Split off the components of ``point`` corresponding to - # ``factor``. - factor_part = point[0:factor.dimension()] - if not factor.contains_strict(factor_part): - return False - point = point[factor.dimension():] + return all([p in f for (p, f) in zip(point, self.factors())]) - return True def factors(self): @@ -457,7 +559,14 @@ class CartesianProduct(SymmetricCone): Return a tuple containing the factors (in order) of this cartesian product. - EXAMPLES: + Returns + ------- + + tuple of :class:`SymmetricCone`. + The factors of this cartesian product. + + Examples + -------- >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(2)) >>> len(K.factors()) @@ -466,6 +575,7 @@ class CartesianProduct(SymmetricCone): """ return self._factors + def cvxopt_dims(self): """ Return a dictionary of dimensions corresponding to the factors @@ -474,7 +584,14 @@ class CartesianProduct(SymmetricCone): http://cvxopt.org/userguide/coneprog.html#linear-cone-programs - EXAMPLES: + Returns + ------- + + dict + A dimension dictionary suitable to feed to CVXOPT. + + Examples + -------- >>> K = CartesianProduct(NonnegativeOrthant(3), ... IceCream(2),