src/dunshire/cones.py

   1 """
   2 Class definitions for all of the symmetric cones (and their superclass,
   3 :class:`SymmetricCone`) supported by CVXOPT.
   4 """
   5
   6 from cvxopt import matrix
   7 from matrices import eigenvalues, norm
   8 import options
   9
  10 class SymmetricCone:
  11     """
  12     An instance of a symmetric (self-dual and homogeneous) cone.
  13
  14     There are three types of symmetric cones supported by CVXOPT:
  15
  16       1. The nonnegative orthant in the real n-space.
  17       2. The Lorentz "ice cream" cone, or the second-order cone.
  18       3. The cone of symmetric positive-semidefinite matrices.
  19
  20     This class is intended to encompass them all.
  21
  22     When constructing a single symmetric cone (i.e. not a
  23     :class:`CartesianProduct` of them), the only information that we
  24     need is its dimension. We take that dimension as a parameter, and
  25     store it for later.
  26
  27     Parameters
  28     ----------
  29
  30     dimension : int
  31         The dimension of this cone.
  32
  33     Raises
  34     ------
  35
  36     ValueError
  37         If you try to create a cone with dimension zero or less.
  38
  39     """
  40     def __init__(self, dimension):
  41         """
  42         A generic constructor for symmetric cones.
  43         """
  44         if dimension <= 0:
  45             raise ValueError('cones must have dimension greater than zero')
  46
  47         self._dimension = dimension
  48
  49
  50     def __contains__(self, point):
  51         """
  52         Return whether or not ``point`` belongs to this cone.
  53
  54         Parameters
  55         ----------
  56
  57         point : matrix
  58             The point to test for membership in this cone.
  59
  60         Raises
  61         ------
  62
  63         NotImplementedError
  64             Always, this method must be implemented in subclasses.
  65
  66         Examples
  67         --------
  68
  69             >>> K = SymmetricCone(5)
  70             >>> matrix([1,2]) in K
  71             Traceback (most recent call last):
  72             ...
  73             NotImplementedError
  74
  75         """
  76         raise NotImplementedError
  77
  78
  79
  80     def dimension(self):
  81         """
  82         Return the dimension of this symmetric cone.
  83
  84         The dimension of this symmetric cone is recorded during its
  85         creation. This method simply returns the recorded value, and
  86         should not need to be overridden in subclasses. We prefer to do
  87         any special computation in ``__init__()`` and record the result
  88         in ``self._dimension``.
  89
  90         Returns
  91         -------
  92
  93         int
  94             The stored dimension (from when this cone was constructed)
  95             of this cone.
  96
  97         Examples
  98         --------
  99
 100             >>> K = SymmetricCone(5)
 101             >>> K.dimension()
 102             5
 103
 104         """
 105         return self._dimension
 106
 107
 108 class NonnegativeOrthant(SymmetricCone):
 109     """
 110     The nonnegative orthant in the given number of dimensions.
 111
 112     Examples
 113     --------
 114
 115         >>> K = NonnegativeOrthant(3)
 116         >>> print(K)
 117         Nonnegative orthant in the real 3-space
 118
 119     """
 120     def __str__(self):
 121         """
 122         Output a human-readable description of myself.
 123         """
 124         tpl = 'Nonnegative orthant in the real {:d}-space'
 125         return tpl.format(self.dimension())
 126
 127
 128     def __contains__(self, point):
 129         """
 130         Return whether or not ``point`` belongs to this cone.
 131
 132         Since this test is expected to work on points whose components
 133         are floating point numbers, it doesn't make any sense to
 134         distinguish between strict and non-strict containment -- the
 135         test uses a tolerance parameter.
 136
 137         This test will err on the side of caution, and return ``True``
 138         only when the ``point`` lies inside this cone *past* the
 139         tolerance guarantee. That means if you ask whether or not a
 140         boundary point lies in this cone, you will get ``False`` as your
 141         answer.
 142
 143         Parameters
 144         ----------
 145
 146         point : matrix
 147             A :class:`cvxopt.base.matrix` having dimensions ``(n,1)``
 148             where ``n`` is the :meth:`dimension` of this cone.
 149
 150         Returns
 151         -------
 152
 153         bool
 154
 155            ``True`` if ``point`` belongs to this cone, ``False`` otherwise.
 156
 157         Raises
 158         ------
 159
 160         TypeError
 161             If ``point`` is not a :class:`cvxopt.base.matrix`.
 162
 163         TypeError
 164             If ``point`` has the wrong dimensions.
 165
 166         Examples
 167         --------
 168
 169         All of these coordinates are positive enough:
 170
 171             >>> K = NonnegativeOrthant(3)
 172             >>> matrix([1,2,3]) in K
 173             True
 174
 175         The one negative coordinate pushes this point outside of ``K``:
 176
 177             >>> K = NonnegativeOrthant(3)
 178             >>> matrix([1,-0.1,3]) in K
 179             False
 180
 181         However, not even a boundary point is considered inside of ``K``:
 182             >>> K = NonnegativeOrthant(3)
 183             >>> matrix([1,0,3]) in K
 184             False
 185
 186         Junk arguments don't work:
 187
 188             >>> K = NonnegativeOrthant(3)
 189             >>> [1,2,3] in K
 190             Traceback (most recent call last):
 191             ...
 192             TypeError: the given point is not a cvxopt.base.matrix
 193
 194             >>> K = NonnegativeOrthant(3)
 195             >>> matrix([1,2]) in K
 196             Traceback (most recent call last):
 197             ...
 198             TypeError: the given point has the wrong dimensions
 199
 200         """
 201         if not isinstance(point, matrix):
 202             raise TypeError('the given point is not a cvxopt.base.matrix')
 203         if not point.size == (self.dimension(), 1):
 204             raise TypeError('the given point has the wrong dimensions')
 205
 206         return all([x > options.ABS_TOL for x in point])
 207
 208
 209
 210 class IceCream(SymmetricCone):
 211     """
 212     The Lorentz "ice cream" cone in the given number of dimensions.
 213
 214     Examples
 215     --------
 216
 217         >>> K = IceCream(3)
 218         >>> print(K)
 219         Lorentz "ice cream" cone in the real 3-space
 220
 221     """
 222     def __str__(self):
 223         """
 224         Output a human-readable description of myself.
 225         """
 226         tpl = 'Lorentz "ice cream" cone in the real {:d}-space'
 227         return tpl.format(self.dimension())
 228
 229
 230     def __contains__(self, point):
 231         """
 232         Return whether or not ``point`` belongs to this cone.
 233
 234         Since this test is expected to work on points whose components
 235         are floating point numbers, it doesn't make any sense to
 236         distinguish between strict and non-strict containment -- the
 237         test uses a tolerance parameter.
 238
 239         This test will err on the side of caution, and return ``True``
 240         only when the ``point`` lies inside this cone *past* the
 241         tolerance guarantee. That means if you ask whether or not a
 242         boundary point lies in this cone, you will get ``False`` as your
 243         answer.
 244
 245         Parameters
 246         ----------
 247
 248         point : matrix
 249             A :class:`cvxopt.base.matrix` having dimensions ``(n,1)``
 250             where ``n`` is the :meth:`dimension` of this cone.
 251
 252         Returns
 253         -------
 254
 255         bool
 256
 257            ``True`` if ``point`` belongs to this cone, ``False`` otherwise.
 258
 259         Raises
 260         ------
 261
 262         TypeError
 263             If ``point`` is not a :class:`cvxopt.base.matrix`.
 264
 265         TypeError
 266             If ``point`` has the wrong dimensions.
 267
 268         Examples
 269         --------
 270
 271         This point lies well within the ice cream cone:
 272
 273             >>> K = IceCream(3)
 274             >>> matrix([1,0.5,0.5]) in K
 275             True
 276
 277         But this one lies on its boundary and runs foul of the tolerance:
 278
 279             >>> K = IceCream(3)
 280             >>> matrix([1,0,1]) in K
 281             False
 282
 283         This point lies entirely outside of the ice cream cone:
 284
 285             >>> K = IceCream(3)
 286             >>> matrix([1,1,1]) in K
 287             False
 288
 289         Junk arguments don't work:
 290
 291             >>> K = IceCream(3)
 292             >>> [1,2,3] in K
 293             Traceback (most recent call last):
 294             ...
 295             TypeError: the given point is not a cvxopt.base.matrix
 296
 297             >>> K = IceCream(3)
 298             >>> matrix([1,2]) in K
 299             Traceback (most recent call last):
 300             ...
 301             TypeError: the given point has the wrong dimensions
 302
 303         """
 304         if not isinstance(point, matrix):
 305             raise TypeError('the given point is not a cvxopt.base.matrix')
 306         if not point.size == (self.dimension(), 1):
 307             raise TypeError('the given point has the wrong dimensions')
 308
 309         height = point[0]
 310         if self.dimension() == 1:
 311             # In one dimension, the ice cream cone is the nonnegative
 312             # orthant.
 313             return height > options.ABS_TOL
 314         else:
 315             radius = point[1:]
 316             return norm(radius) < (height - options.ABS_TOL)
 317
 318
 319
 320 class SymmetricPSD(SymmetricCone):
 321     r"""
 322     The cone of real symmetric positive-semidefinite matrices.
 323
 324     This cone has a dimension ``n`` associated with it, but we let ``n``
 325     refer to the dimension of the domain of our matrices and not the
 326     dimension of the (much larger) space in which the matrices
 327     themselves live. In other words, our ``n`` is the ``n`` that appears
 328     in the usual notation :math:`S^{n}` for symmetric matrices.
 329
 330     As a result, the cone ``SymmetricPSD(n)`` lives in a space of dimension
 331     :math:`\left(n^{2} + n\right)/2)`.
 332
 333     Examples
 334     --------
 335
 336         >>> K = SymmetricPSD(3)
 337         >>> print(K)
 338         Cone of symmetric positive-semidefinite matrices on the real 3-space
 339         >>> K.dimension()
 340         3
 341
 342     """
 343     def __str__(self):
 344         """
 345         Output a human-readable description of myself.
 346         """
 347         tpl = 'Cone of symmetric positive-semidefinite matrices ' \
 348               'on the real {:d}-space'
 349         return tpl.format(self.dimension())
 350
 351
 352     def __contains__(self, point):
 353         """
 354         Return whether or not ``point`` belongs to this cone.
 355
 356         Since this test is expected to work on points whose components
 357         are floating point numbers, it doesn't make any sense to
 358         distinguish between strict and non-strict containment -- the
 359         test uses a tolerance parameter.
 360
 361         This test will err on the side of caution, and return ``True``
 362         only when the ``point`` lies inside this cone *past* the
 363         tolerance guarantee. That means if you ask whether or not a
 364         boundary point lies in this cone, you will get ``False`` as your
 365         answer.
 366
 367         Parameters
 368         ----------
 369
 370         point : matrix
 371             A :class:`cvxopt.base.matrix` having dimensions ``(n,n)``
 372             where ``n`` is the :meth:`dimension` of this cone.
 373
 374         Returns
 375         -------
 376
 377         bool
 378
 379            ``True`` if ``point`` belongs to this cone, ``False`` otherwise.
 380
 381         Raises
 382         ------
 383
 384         TypeError
 385             If ``point`` is not a :class:`cvxopt.base.matrix`.
 386
 387         TypeError
 388             If ``point`` has the wrong dimensions.
 389
 390         Examples
 391         --------
 392
 393         These all lie in the interior of the Symmetric PSD cone:
 394
 395             >>> K = SymmetricPSD(2)
 396             >>> matrix([[1,0],[0,1]]) in K
 397             True
 398
 399             >>> K = SymmetricPSD(3)
 400             >>> matrix([[2,-1,0],[-1,2,-1],[0,-1,2]]) in K
 401             True
 402
 403             >>> K = SymmetricPSD(5)
 404             >>> A = matrix([[5,4,3,2,1],
 405             ...            [4,5,4,3,2],
 406             ...            [3,4,5,4,3],
 407             ...            [2,3,4,5,4],
 408             ...            [1,2,3,4,5]])
 409             >>> A in K
 410             True
 411
 412         But any matrix with a zero eigenvalue will lie on the boundary
 413         of the Symmetric PSD cone and run foul of our tolerance:
 414
 415             >>> K = SymmetricPSD(2)
 416             >>> matrix([[0,0],[0,0]]) in K
 417             False
 418
 419             >>> K = SymmetricPSD(5)
 420             >>> A = matrix([[1,0,0,0,0],
 421             ...            [0,1,0,0,0],
 422             ...            [0,0,0,0,0],
 423             ...            [0,0,0,1,0],
 424             ...            [0,0,0,0,1]])
 425             >>> A in K
 426             False
 427
 428         Junk arguments don't work:
 429
 430             >>> K = SymmetricPSD(2)
 431             >>> [[1,2],[2,3]] in K
 432             Traceback (most recent call last):
 433             ...
 434             TypeError: the given point is not a cvxopt.base.matrix
 435
 436             >>> K = SymmetricPSD(3)
 437             >>> matrix([[1,2],[3,4]]) in K
 438             Traceback (most recent call last):
 439             ...
 440             TypeError: the given point has the wrong dimensions
 441
 442         """
 443         if not isinstance(point, matrix):
 444             raise TypeError('the given point is not a cvxopt.base.matrix')
 445         if not point.size == (self.dimension(), self.dimension()):
 446             raise TypeError('the given point has the wrong dimensions')
 447         if not point.typecode == 'd':
 448             point = matrix(point, (self.dimension(), self.dimension()), 'd')
 449         return all([e > options.ABS_TOL for e in eigenvalues(point)])
 450
 451
 452
 453 class CartesianProduct(SymmetricCone):
 454     """
 455     A cartesian product of symmetric cones, which is itself a symmetric
 456     cone.
 457
 458     Examples
 459     --------
 460
 461         >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(2))
 462         >>> print(K)
 463         Cartesian product of dimension 5 with 2 factors:
 464           * Nonnegative orthant in the real 3-space
 465           * Lorentz "ice cream" cone in the real 2-space
 466
 467     """
 468     def __init__(self, *factors):
 469         my_dimension = sum([f.dimension() for f in factors])
 470         super().__init__(my_dimension)
 471         self._factors = factors
 472
 473
 474     def __str__(self):
 475         """
 476         Output a human-readable description of myself.
 477         """
 478         tpl = 'Cartesian product of dimension {:d} with {:d} factors:'
 479         tpl += '\n  * {!s}' * len(self.factors())
 480         format_args = [self.dimension(), len(self.factors())]
 481         format_args += list(self.factors())
 482         return tpl.format(*format_args)
 483
 484
 485     def __contains__(self, point):
 486         """
 487         Return whether or not ``point`` belongs to this cone.
 488
 489         The ``point`` is expected to be a tuple of points which will be
 490         tested for membership in this cone's factors. If each point in
 491         the tuple belongs to its corresponding factor, then the whole
 492         point belongs to this cone. Otherwise, it doesn't.
 493
 494         Parameters
 495         ----------
 496
 497         point : tuple of matrix
 498             A tuple of :class:`cvxopt.base.matrix` corresponding to the
 499             :meth:`factors` of this cartesian product.
 500
 501         Returns
 502         -------
 503
 504         bool
 505
 506            ``True`` if ``point`` belongs to this cone, ``False`` otherwise.
 507
 508         Raises
 509         ------
 510
 511         TypeError
 512             If ``point`` is not a tuple of :class:`cvxopt.base.matrix`.
 513
 514         TypeError
 515             If any element of ``point`` has the wrong dimensions.
 516
 517         Examples
 518         --------
 519
 520         The result depends on how containment is defined for our factors:
 521
 522             >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3))
 523             >>> (matrix([1,2,3]), matrix([1,0.5,0.5])) in K
 524             True
 525
 526             >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3))
 527             >>> (matrix([0,0,0]), matrix([1,0,1])) in K
 528             False
 529
 530             >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3))
 531             >>> (matrix([1,1,1]), matrix([1,1,1])) in K
 532             False
 533
 534             >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3))
 535             >>> (matrix([1,-1,1]), matrix([1,0,1])) in K
 536             False
 537
 538         Junk arguments don't work:
 539
 540             >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3))
 541             >>> [[1,2,3],[4,5,6]] in K
 542             Traceback (most recent call last):
 543             ...
 544             TypeError: the given point is not a cvxopt.base.matrix
 545
 546             >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(3))
 547             >>> (matrix([1,2]), matrix([3,4,5,6])) in K
 548             Traceback (most recent call last):
 549             ...
 550             TypeError: the given point has the wrong dimensions
 551
 552         """
 553         return all([p in f for (p,f) in zip(point, self.factors())])
 554
 555
 556
 557     def factors(self):
 558         """
 559         Return a tuple containing the factors (in order) of this
 560         cartesian product.
 561
 562         Returns
 563         -------
 564
 565         tuple of :class:`SymmetricCone`.
 566             The factors of this cartesian product.
 567
 568         Examples
 569         --------
 570
 571             >>> K = CartesianProduct(NonnegativeOrthant(3), IceCream(2))
 572             >>> len(K.factors())
 573             2
 574
 575         """
 576         return self._factors
 577
 578
 579     def cvxopt_dims(self):
 580         """
 581         Return a dictionary of dimensions corresponding to the factors
 582         of this cartesian product. The format of this dictionary is
 583         described in the CVXOPT user's guide:
 584
 585           http://cvxopt.org/userguide/coneprog.html#linear-cone-programs
 586
 587         Returns
 588         -------
 589
 590         dict
 591             A dimension dictionary suitable to feed to CVXOPT.
 592
 593         Examples
 594         --------
 595
 596             >>> K = CartesianProduct(NonnegativeOrthant(3),
 597             ...                      IceCream(2),
 598             ...                      IceCream(3))
 599             >>> d = K.cvxopt_dims()
 600             >>> (d['l'], d['q'], d['s'])
 601             (3, [2, 3], [])
 602
 603         """
 604         dims = {'l':0, 'q':[], 's':[]}
 605         dims['l'] += sum([K.dimension()
 606                           for K in self.factors()
 607                           if isinstance(K, NonnegativeOrthant)])
 608         dims['q'] = [K.dimension()
 609                      for K in self.factors()
 610                      if isinstance(K, IceCream)]
 611         dims['s'] = [K.dimension()
 612                      for K in self.factors()
 613                      if isinstance(K, SymmetricPSD)]
 614         return dims
 615
 616