X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=src%2Fdunshire%2Fcones.py;h=61feab453bd0866b54c053d2871d722a63ec3730;hb=23a5893ede3a22653128d3e7a66fb3f6b80616b8;hp=2609b16eea729651c7d005c226d5f47ab4cd616c;hpb=002b5370da24f083d2088c3482cf076615a13563;p=dunshire.git diff --git a/src/dunshire/cones.py b/src/dunshire/cones.py index 2609b16..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 eigenvalues, norm +import options class SymmetricCone: """ @@ -17,20 +18,28 @@ class SymmetricCone: 3. The cone of symmetric positive-semidefinite matrices. 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. + 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. - INPUT: + Parameters + ---------- - - ``dimension`` -- the dimension of this cone. + dimension : int + The dimension of this cone. + Raises + ------ + + ValueError + If you try to create a cone with dimension zero or less. + + """ + def __init__(self, dimension): + """ + A generic constructor for symmetric cones. """ if dimension <= 0: raise ValueError('cones must have dimension greater than zero') @@ -42,32 +51,32 @@ class SymmetricCone: """ Return whether or not ``point`` belongs to this cone. - EXAMPLES: + Parameters + ---------- - >>> K = SymmetricCone(5) - >>> matrix([1,2]) in K - Traceback (most recent call last): - ... - NotImplementedError + point : matrix + The point to test for membership in this cone. - """ - raise NotImplementedError + Raises + ------ - def contains_strict(self, point): - """ - Return whether or not ``point`` belongs to the interior - of this cone. + NotImplementedError + Always, this method must be implemented in subclasses. - EXAMPLES: + Examples + -------- >>> K = SymmetricCone(5) - >>> K.contains_strict(matrix([1,2])) + >>> matrix([1,2]) in K Traceback (most recent call last): ... NotImplementedError + """ raise NotImplementedError + + def dimension(self): """ Return the dimension of this symmetric cone. @@ -78,7 +87,15 @@ class SymmetricCone: any special computation in ``__init__()`` and record the result in ``self._dimension``. - EXAMPLES: + Returns + ------- + + int + The stored dimension (from when this cone was constructed) + of this cone. + + Examples + -------- >>> K = SymmetricCone(5) >>> K.dimension() @@ -90,9 +107,10 @@ class SymmetricCone: 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) @@ -106,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. + 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 = NonnegativeOrthant(3) - >>> matrix([1,2,3]) in K - True - - >>> K = NonnegativeOrthant(3) - >>> matrix([1,-0.1,3]) in K - False - - >>> 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 @@ -188,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) @@ -215,85 +225,65 @@ class IceCream(SymmetricCone): """ 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. + 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) - >>> matrix([1,0.5,0.5]) in K - True - - >>> K = IceCream(3) - >>> matrix([1,0,1]) in K - True + point : matrix + A :class:`cvxopt.base.matrix` having dimensions ``(n,1)`` + where ``n`` is the :meth:`dimension` of this cone. - >>> K = IceCream(3) - >>> matrix([1,1,1]) in K - False + Returns + ------- - >>> K = IceCream(3) - >>> [1,2,3] in K - Traceback (most recent call last): - ... - TypeError: the given point is not a cvxopt.base.matrix + bool - >>> K = IceCream(3) - >>> matrix([1,2]) in K - Traceback (most recent call last): - ... - TypeError: the given point has the wrong dimensions + ``True`` if ``point`` belongs to this cone, ``False`` otherwise. - """ - 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') + Raises + ------ - 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) + TypeError + If ``point`` is not a :class:`cvxopt.base.matrix`. + TypeError + If ``point`` has the wrong dimensions. - def contains_strict(self, point): - """ - Return whether or not ``point`` belongs to the interior - of this cone. + Examples + -------- - INPUT: + This point lies well within the ice cream cone: - An instance of the ``cvxopt.base.matrix`` class having - dimensions ``(n,1)`` where ``n`` is the dimension of this cone. + >>> K = IceCream(3) + >>> matrix([1,0.5,0.5]) in K + True - EXAMPLES: + This one lies on its boundary: >>> K = IceCream(3) - >>> K.contains_strict(matrix([1,0.5,0.5])) + >>> matrix([1,0,1]) in K True - >>> K = IceCream(3) - >>> K.contains_strict(matrix([1,0,1])) - False + This point lies entirely outside of the ice cream cone: >>> K = IceCream(3) - >>> K.contains_strict(matrix([1,1,1])) + >>> matrix([1,1,1]) in K False + Junk arguments don't work: + >>> K = IceCream(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 = IceCream(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 @@ -308,26 +298,28 @@ 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) + 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 `S^{n}` for symmetric matrices. + in the usual notation :math:`S^{n}` for symmetric matrices. As a result, the cone ``SymmetricPSD(n)`` lives in a space of dimension - ``(n**2 + n)/2)``. + :math:`\left(n^{2} + n\right)/2)`. - EXAMPLES: + Examples + -------- >>> K = SymmetricPSD(3) >>> print(K) @@ -349,116 +341,97 @@ class SymmetricPSD(SymmetricCone): """ Return whether or not ``point`` belongs to this cone. - INPUT: - - An instance of the ``cvxopt.base.matrix`` class having - dimensions ``(n,n)`` 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 = SymmetricPSD(2) - >>> matrix([[1,0],[0,1]]) in K - True - - >>> K = SymmetricPSD(2) - >>> matrix([[0,0],[0,0]]) in K - True - - >>> K = SymmetricPSD(3) - >>> matrix([[2,-1,0],[-1,2,-1],[0,-1,2]]) in K - True - - >>> 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 + point : matrix + A :class:`cvxopt.base.matrix` having dimensions ``(n,n)`` + where ``n`` is the :meth:`dimension` of this cone. - >>> 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 - - >>> K = SymmetricPSD(2) - >>> [[1,2],[2,3]] in K - Traceback (most recent call last): - ... - TypeError: the given point is not a cvxopt.base.matrix + Returns + ------- - >>> K = SymmetricPSD(3) - >>> matrix([[1,2],[3,4]]) in K - Traceback (most recent call last): - ... - TypeError: the given point has the wrong dimensions + bool - """ - if not isinstance(point, matrix): - raise TypeError('the given point is not a cvxopt.base.matrix') - 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') - return all([e >= 0 for e in eigenvalues(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,n)`` where ``n`` is the dimension of this cone. - Its type code must be 'd'. + Examples + -------- - EXAMPLES: + These all lie in the interior of the Symmetric PSD cone: >>> K = SymmetricPSD(2) - >>> K.contains_strict(matrix([[1,0],[0,1]])) + >>> matrix([[1,0],[0,1]]) in K True - >>> K = SymmetricPSD(2) - >>> K.contains_strict(matrix([[0,0],[0,0]])) - False - >>> K = SymmetricPSD(3) >>> matrix([[2,-1,0],[-1,2,-1],[0,-1,2]]) in K True >>> 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]]) + ... [4,5,4,3,2], + ... [3,4,5,4,3], + ... [2,3,4,5,4], + ... [1,2,3,4,5]]) >>> A in K True + 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]]) - >>> K.contains_strict(A) - False + ... [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) - >>> K.contains_strict([[1,2],[2,3]]) + >>> [[1,2],[2,3]] in K Traceback (most recent call last): ... TypeError: the given point is not a cvxopt.base.matrix >>> K = SymmetricPSD(3) - >>> K.contains_strict(matrix([[1,2],[3,4]])) + >>> matrix([[1,2],[3,4]]) in K Traceback (most recent call last): ... TypeError: the given point has the wrong dimensions @@ -470,7 +443,11 @@ class SymmetricPSD(SymmetricCone): raise TypeError('the given point has the wrong dimensions') if not point.typecode == 'd': point = matrix(point, (self.dimension(), self.dimension()), 'd') - return all([e > 0 for e in eigenvalues(point)]) + 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)]) + class CartesianProduct(SymmetricCone): @@ -478,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) @@ -492,6 +470,7 @@ class CartesianProduct(SymmetricCone): super().__init__(my_dimension) self._factors = factors + def __str__(self): """ Output a human-readable description of myself. @@ -502,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: - - An instance of the ``cvxopt.base.matrix`` class having - dimensions ``(n,1)`` where ``n`` is the dimension of this cone. - - EXAMPLES: - - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([1,2,3,1,0.5,0.5]) in K - True - - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([0,0,0,1,0,1]) in K - True - - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([1,1,1,1,1,1]) in K - False + 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. - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([1,-1,1,1,0,1]) in K - False - - >>> 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 + Parameters + ---------- - >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3)) - >>> matrix([1,2]) in K - Traceback (most recent call last): - ... - TypeError: the given point has the wrong dimensions + point : tuple of matrix + A tuple of :class:`cvxopt.base.matrix` corresponding to the + :meth:`factors` of this cartesian product. - """ - 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') + Returns + ------- - 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():] + bool - return True + ``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 tuple of :class:`cvxopt.base.matrix`. - INPUT: + TypeError + If any element of ``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: + 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') + return all([p in f for (p, f) in zip(point, self.factors())]) - 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 True def factors(self): @@ -624,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()) @@ -633,6 +575,7 @@ class CartesianProduct(SymmetricCone): """ return self._factors + def cvxopt_dims(self): """ Return a dictionary of dimensions corresponding to the factors @@ -641,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),