mjo/cone/tests.py

   1 """
   2 Additional tests for the mjo.cone.cone module. These are extra
   3 properties that we'd like to check, but which are overkill for inclusion
   4 into Sage.
   5 """
   6
   7 # Sage doesn't load ~/.sage/init.sage during testing (sage -t), so we
   8 # have to explicitly mangle our sitedir here so that "mjo.cone"
   9 # resolves.
  10 from os.path import abspath
  11 from site import addsitedir
  12 addsitedir(abspath('../../'))
  13
  14 from sage.all import *
  15
  16 # The double-import is needed to get the underscore methods.
  17 from mjo.cone.cone import *
  18
  19
  20 def _restrict_to_subspace(K, W):
  21     r"""
  22     Restrict ``K`` (up to linear isomorphism) to a vector subspace.
  23
  24     This operation not only restricts the cone to a subspace of its
  25     ambient space, but also represents the rays of the cone in a new
  26     (smaller) lattice corresponding to the subspace. The resulting
  27     cone will be linearly isomorphic (but not equal) to the
  28     desired restriction, since it has likely undergone a change of
  29     basis.
  30
  31     To explain the difficulty, consider the cone ``K =
  32     Cone([(1,1,1)])`` having a single ray. The span of ``K`` is a
  33     one-dimensional subspace containing ``K``, yet we have no way to
  34     perform operations like "dual of" in the subspace. To represent
  35     ``K`` in the space ``K.span()``, we must perform a change of basis
  36     and write its sole ray as ``(1,0,0)``. Now the restricted
  37     ``Cone([(1,)])`` is linearly isomorphic (but of course not equal) to
  38     ``K`` interpreted as living in ``K.span()``.
  39
  40     INPUT:
  41
  42     - ``K`` -- The cone to restrict.
  43
  44     - ``W`` -- The subspace into which ``K`` will be restricted.
  45
  46     OUTPUT:
  47
  48     A new cone in a sublattice corresponding to ``W``.
  49
  50     REFERENCES:
  51
  52     M. Orlitzky. The Lyapunov rank of an improper cone.
  53     http://www.optimization-online.org/DB_HTML/2015/10/5135.html
  54
  55     EXAMPLES:
  56
  57     Restricting a solid cone to its own span returns a cone linearly
  58     isomorphic to the original::
  59
  60         sage: K = Cone([(1,2,3),(-1,1,0),(9,0,-2)])
  61         sage: K.is_solid()
  62         True
  63         sage: _restrict_to_subspace(K, K.span()).rays()
  64         N(-1,  1,  0),
  65         N( 1,  0,  0),
  66         N( 9, -6, -1)
  67         in 3-d lattice N
  68
  69     A single ray restricted to its own span has the same
  70     representation regardless of the ambient space::
  71
  72         sage: K = Cone([(1,0)])
  73         sage: K_S = _restrict_to_subspace(K, K.span()).rays()
  74         sage: K_S
  75         N(1)
  76         in 1-d lattice N
  77         sage: K = Cone([(1,1,1)])
  78         sage: K_S = _restrict_to_subspace(K, K.span()).rays()
  79         sage: K_S
  80         N(1)
  81         in 1-d lattice N
  82
  83     Restricting to a trivial space gives the trivial cone::
  84
  85         sage: K = Cone([(8,3,-1,0),(9,2,2,0),(-4,6,7,0)])
  86         sage: trivial_space = K.lattice().vector_space().span([])
  87         sage: _restrict_to_subspace(K, trivial_space)
  88         0-d cone in 0-d lattice N
  89
  90     TESTS:
  91
  92     Restricting a cone to its own span results in a solid cone::
  93
  94         sage: set_random_seed()
  95         sage: K = random_cone(max_ambient_dim = 8)
  96         sage: K_S = _restrict_to_subspace(K, K.span())
  97         sage: K_S.is_solid()
  98         True
  99
 100     Restricting a cone to its span should not affect the number of
 101     rays in the cone::
 102
 103         sage: set_random_seed()
 104         sage: K = random_cone(max_ambient_dim = 8)
 105         sage: K_S = _restrict_to_subspace(K, K.span())
 106         sage: K.nrays() == K_S.nrays()
 107         True
 108
 109     Restricting a cone to its span should not affect its dimension::
 110
 111         sage: set_random_seed()
 112         sage: K = random_cone(max_ambient_dim = 8)
 113         sage: K_S = _restrict_to_subspace(K, K.span())
 114         sage: K.dim() == K_S.dim()
 115         True
 116
 117     Restricting a cone to its span should not affects its lineality::
 118
 119         sage: set_random_seed()
 120         sage: K = random_cone(max_ambient_dim = 8)
 121         sage: K_S = _restrict_to_subspace(K, K.span())
 122         sage: K.lineality() == K_S.lineality()
 123         True
 124
 125     Restricting a cone to its span should not affect the number of
 126     facets it has::
 127
 128         sage: set_random_seed()
 129         sage: K = random_cone(max_ambient_dim = 8)
 130         sage: K_S = _restrict_to_subspace(K, K.span())
 131         sage: len(K.facets()) == len(K_S.facets())
 132         True
 133
 134     Restricting a solid cone to its span is a linear isomorphism
 135     and should not affect the dimension of its ambient space::
 136
 137         sage: set_random_seed()
 138         sage: K = random_cone(max_ambient_dim = 8, solid = True)
 139         sage: K_S = _restrict_to_subspace(K, K.span())
 140         sage: K.lattice_dim() == K_S.lattice_dim()
 141         True
 142
 143     Restricting a solid cone to its span is a linear isomorphism
 144     that establishes a one-to-one correspondence of discrete
 145     complementarity sets::
 146
 147         sage: set_random_seed()
 148         sage: K = random_cone(max_ambient_dim = 8, solid = True)
 149         sage: K_S = _restrict_to_subspace(K, K.span())
 150         sage: dcs1 = K.discrete_complementarity_set()
 151         sage: dcs2 = K_S.discrete_complementarity_set()
 152         sage: len(dcs1) == len(dcs2)
 153         True
 154
 155     Restricting a solid cone to its span is a linear isomorphism
 156     under which Lyapunov rank (the length of a Lyapunov-like basis)
 157     is invariant::
 158
 159         sage: set_random_seed()
 160         sage: K = random_cone(max_ambient_dim = 8, solid = True)
 161         sage: K_S = _restrict_to_subspace(K, K.span())
 162         sage: LL1 = K.lyapunov_like_basis()
 163         sage: LL2 = K_S.lyapunov_like_basis()
 164         sage: len(LL1) == len(LL2)
 165         True
 166
 167     If we restrict a cone to a subspace of its span, the resulting
 168     cone should have the same dimension as the subspace::
 169
 170         sage: set_random_seed()
 171         sage: K = random_cone(max_ambient_dim = 8)
 172         sage: W_basis = random_sublist(K.rays(), 0.5)
 173         sage: W = K.lattice().vector_space().span(W_basis)
 174         sage: K_W = _restrict_to_subspace(K,W)
 175         sage: K_W.lattice_dim() == W.dimension()
 176         True
 177
 178     Through a series of restrictions, any closed convex cone can be
 179     reduced to a cartesian product with a proper factor [Orlitzky]_::
 180
 181         sage: set_random_seed()
 182         sage: K = random_cone(max_ambient_dim = 8)
 183         sage: K_S = _restrict_to_subspace(K, K.span())
 184         sage: P = K_S.dual().span()
 185         sage: K_SP = _restrict_to_subspace(K_S, P)
 186         sage: K_SP.is_proper()
 187         True
 188     """
 189     # We want to intersect this cone with ``W``. We can do that via
 190     # cone intersection, so we first turn the space ``W`` into a cone.
 191     W_rays = W.basis() + [ -b for b in W.basis() ]
 192     W_cone = Cone(W_rays, lattice=K.lattice())
 193     K = K.intersection(W_cone)
 194
 195     # Now every generator of ``K`` should belong to ``W``.
 196     K_W_rays = [ W.coordinate_vector(r) for r in K.rays() ]
 197
 198     L = ToricLattice(W.dimension())
 199     return Cone(K_W_rays, lattice=L)
 200
 201
 202
 203 #
 204 # Tests for _restrict_to_subspace.
 205 #
 206 def _look_isomorphic(K1, K2):
 207     r"""
 208     Test whether or not ``K1`` and ``K2`` look linearly isomorphic.
 209
 210     This is a hack to get around the fact that it's difficult to tell
 211     when two cones are linearly isomorphic. Instead, we check a list of
 212     properties that should be preserved under linear isomorphism.
 213
 214     OUTPUT:
 215
 216     ``True`` if ``K1`` and ``K2`` look isomorphic, or ``False``
 217     if we can prove that they are not isomorphic.
 218
 219     EXAMPLES:
 220
 221     Any proper cone with three generators in `\mathbb{R}^{3}` is
 222     isomorphic to the nonnegative orthant::
 223
 224         sage: K1 = Cone([(1,0,0), (0,1,0), (0,0,1)])
 225         sage: K2 = Cone([(1,2,3), (3, 18, 4), (66, 51, 0)])
 226         sage: _look_isomorphic(K1, K2)
 227         True
 228
 229     Negating a cone gives you an isomorphic cone::
 230
 231         sage: K = Cone([(0,2,-5), (-6, 2, 4), (0, 51, 0)])
 232         sage: _look_isomorphic(K, -K)
 233         True
 234
 235     TESTS:
 236
 237     Any cone is isomorphic to itself::
 238
 239         sage: K = random_cone(max_ambient_dim = 8)
 240         sage: _look_isomorphic(K, K)
 241         True
 242
 243     After applying an invertible matrix to the rows of a cone, the
 244     result should is isomorphic to the cone we started with::
 245
 246         sage: K1 = random_cone(max_ambient_dim = 8)
 247         sage: A = random_matrix(QQ, K1.lattice_dim(), algorithm='unimodular')
 248         sage: K2 = Cone( [ A*r for r in K1.rays() ], lattice=K1.lattice())
 249         sage: _look_isomorphic(K1, K2)
 250         True
 251
 252     """
 253     if K1.lattice_dim() != K2.lattice_dim():
 254         return False
 255
 256     if K1.nrays() != K2.nrays():
 257         return False
 258
 259     if K1.dim() != K2.dim():
 260         return False
 261
 262     if K1.lineality() != K2.lineality():
 263         return False
 264
 265     if K1.is_solid() != K2.is_solid():
 266         return False
 267
 268     if K1.is_strictly_convex() != K2.is_strictly_convex():
 269         return False
 270
 271     if len(K1.lyapunov_like_basis()) != len(K2.lyapunov_like_basis()):
 272         return False
 273
 274     C_of_K1 = K1.discrete_complementarity_set()
 275     C_of_K2 = K2.discrete_complementarity_set()
 276     if len(C_of_K1) != len(C_of_K2):
 277         return False
 278
 279     if len(K1.facets()) != len(K2.facets()):
 280         return False
 281
 282     return True
 283
 284
 285 """
 286 Apply _restrict_to_subspace according to our paper (to obtain our main
 287 result). Test all four parameter combinations::
 288
 289     sage: set_random_seed()
 290     sage: K = random_cone(max_ambient_dim = 8,
 291     ....:                 strictly_convex=False,
 292     ....:                 solid=False)
 293     sage: K_S = _restrict_to_subspace(K, K.span())
 294     sage: K_S2 = K.solid_restriction()
 295     sage: _look_isomorphic(K_S, K_S2)
 296     True
 297     sage: K_SP = _restrict_to_subspace(K_S.dual(), K_S.dual().span()).dual()
 298     sage: K_SP2 = K_S.strict_quotient()
 299     sage: K_SP.is_proper()
 300     True
 301     sage: K_SP2.is_proper()
 302     True
 303     sage: _look_isomorphic(K_SP, K_SP2)
 304     True
 305     sage: K_SP = _restrict_to_subspace(K_S, K_S.dual().span())
 306     sage: K_SP.is_proper()
 307     True
 308     sage: _look_isomorphic(K_SP, K_SP2)
 309     True
 310
 311 ::
 312
 313     sage: set_random_seed()
 314     sage: K = random_cone(max_ambient_dim = 8,
 315     ....:                 strictly_convex=False,
 316     ....:                 solid=True)
 317     sage: K_S = _restrict_to_subspace(K, K.span())
 318     sage: K_S2 = K.solid_restriction()
 319     sage: _look_isomorphic(K_S, K_S2)
 320     True
 321     sage: K_SP = _restrict_to_subspace(K_S.dual(), K_S.dual().span()).dual()
 322     sage: K_SP2 = K_S.strict_quotient()
 323     sage: K_SP.is_proper()
 324     True
 325     sage: K_SP2.is_proper()
 326     True
 327     sage: _look_isomorphic(K_SP, K_SP2)
 328     True
 329     sage: K_SP = _restrict_to_subspace(K_S, K_S.dual().span())
 330     sage: K_SP.is_proper()
 331     True
 332     sage: _look_isomorphic(K_SP, K_SP2)
 333     True
 334
 335 ::
 336
 337     sage: set_random_seed()
 338     sage: K = random_cone(max_ambient_dim = 8,
 339     ....:                 strictly_convex=True,
 340     ....:                 solid=False)
 341     sage: K_S = _restrict_to_subspace(K, K.span())
 342     sage: K_S2 = K.solid_restriction()
 343     sage: _look_isomorphic(K_S, K_S2)
 344     True
 345     sage: K_SP = _restrict_to_subspace(K_S.dual(), K_S.dual().span()).dual()
 346     sage: K_SP2 = K_S.strict_quotient()
 347     sage: K_SP.is_proper()
 348     True
 349     sage: K_SP2.is_proper()
 350     True
 351     sage: _look_isomorphic(K_SP, K_SP2)
 352     True
 353     sage: K_SP = _restrict_to_subspace(K_S, K_S.dual().span())
 354     sage: K_SP.is_proper()
 355     True
 356     sage: _look_isomorphic(K_SP, K_SP2)
 357     True
 358
 359 ::
 360
 361     sage: set_random_seed()
 362     sage: K = random_cone(max_ambient_dim = 8,
 363     ....:                 strictly_convex=True,
 364     ....:                 solid=True)
 365     sage: K_S = _restrict_to_subspace(K, K.span())
 366     sage: K_S2 = K.solid_restriction()
 367     sage: _look_isomorphic(K_S, K_S2)
 368     True
 369     sage: K_SP = _restrict_to_subspace(K_S.dual(), K_S.dual().span()).dual()
 370     sage: K_SP2 = K_S.strict_quotient()
 371     sage: K_SP.is_proper()
 372     True
 373     sage: K_SP2.is_proper()
 374     True
 375     sage: _look_isomorphic(K_SP, K_SP2)
 376     True
 377     sage: K_SP = _restrict_to_subspace(K_S, K_S.dual().span())
 378     sage: K_SP.is_proper()
 379     True
 380     sage: _look_isomorphic(K_SP, K_SP2)
 381     True
 382
 383 Test the proposition in our paper concerning the duals and
 384 restrictions. Generate a random cone, then create a subcone of
 385 it. The operation of dual-taking should then commute with rho. Test
 386 all parameter combinations::
 387
 388     sage: set_random_seed()
 389     sage: J = random_cone(max_ambient_dim = 8,
 390     ....:                 solid=False,
 391     ....:                 strictly_convex=False)
 392     sage: K = Cone(random_sublist(J.rays(), 0.5), lattice=J.lattice())
 393     sage: K_W_star = _restrict_to_subspace(K, J.span()).dual()
 394     sage: K_star_W = _restrict_to_subspace(K.dual(), J.span())
 395     sage: _look_isomorphic(K_W_star, K_star_W)
 396     True
 397
 398 ::
 399
 400     sage: set_random_seed()
 401     sage: J = random_cone(max_ambient_dim = 8,
 402     ....:                 solid=True,
 403     ....:                 strictly_convex=False)
 404     sage: K = Cone(random_sublist(J.rays(), 0.5), lattice=J.lattice())
 405     sage: K_W_star = _restrict_to_subspace(K, J.span()).dual()
 406     sage: K_star_W = _restrict_to_subspace(K.dual(), J.span())
 407     sage: _look_isomorphic(K_W_star, K_star_W)
 408     True
 409
 410 ::
 411
 412     sage: set_random_seed()
 413     sage: J = random_cone(max_ambient_dim = 8,
 414     ....:                 solid=False,
 415     ....:                 strictly_convex=True)
 416     sage: K = Cone(random_sublist(J.rays(), 0.5), lattice=J.lattice())
 417     sage: K_W_star = _restrict_to_subspace(K, J.span()).dual()
 418     sage: K_star_W = _restrict_to_subspace(K.dual(), J.span())
 419     sage: _look_isomorphic(K_W_star, K_star_W)
 420     True
 421
 422 ::
 423
 424     sage: set_random_seed()
 425     sage: J = random_cone(max_ambient_dim = 8,
 426     ....:                 solid=True,
 427     ....:                 strictly_convex=True)
 428     sage: K = Cone(random_sublist(J.rays(), 0.5), lattice=J.lattice())
 429     sage: K_W_star = _restrict_to_subspace(K, J.span()).dual()
 430     sage: K_star_W = _restrict_to_subspace(K.dual(), J.span())
 431     sage: _look_isomorphic(K_W_star, K_star_W)
 432     True
 433
 434 Ensure that ``__restrict_to_subspace(K, K.span())`` and
 435 ``K.solid_restriction()`` are actually equivalent::
 436
 437     sage: set_random_seed()
 438     sage: K = random_cone(max_ambient_dim=8)
 439     sage: K1 = _restrict_to_subspace(K, K.span())
 440     sage: K2 = K.solid_restriction()
 441     sage: _look_isomorphic(K1,K2)
 442     True
 443
 444 Ensure that ``K.__restrict_to_subspace(K,K.dual().span())`` and
 445 ``strict_quotient`` are actually equivalent::
 446
 447     sage: set_random_seed()
 448     sage: K = random_cone(max_ambient_dim=6)
 449     sage: K1 = _restrict_to_subspace(K, K.dual().span())
 450     sage: K2 = K.strict_quotient()
 451     sage: _look_isomorphic(K1,K2)
 452     True
 453 """
 454
 455
 456 #
 457 # Lyapunov rank tests
 458 #
 459 """
 460
 461 The Lyapunov rank is invariant under a linear isomorphism. Check all
 462 combinations of parameters::
 463
 464     sage: K1 = random_cone(max_ambient_dim=8,
 465     ....:                  strictly_convex=True,
 466     ....:                  solid=True)
 467     sage: A = random_matrix(QQ, K1.lattice_dim(), algorithm='unimodular')
 468     sage: K2 = Cone( [ A*r for r in K1.rays() ], lattice=K1.lattice())
 469     sage: K1.lyapunov_rank() == K2.lyapunov_rank()
 470     True
 471
 472 ::
 473
 474     sage: K1 = random_cone(max_ambient_dim=8,
 475     ....:                  strictly_convex=True,
 476     ....:                  solid=False)
 477     sage: A = random_matrix(QQ, K1.lattice_dim(), algorithm='unimodular')
 478     sage: K2 = Cone( [ A*r for r in K1.rays() ], lattice=K1.lattice())
 479     sage: K1.lyapunov_rank() == K2.lyapunov_rank()
 480     True
 481
 482 ::
 483
 484     sage: K1 = random_cone(max_ambient_dim=8,
 485     ....:                  strictly_convex=False,
 486     ....:                  solid=True)
 487     sage: A = random_matrix(QQ, K1.lattice_dim(), algorithm='unimodular')
 488     sage: K2 = Cone( [ A*r for r in K1.rays() ], lattice=K1.lattice())
 489     sage: K1.lyapunov_rank() == K2.lyapunov_rank()
 490     True
 491
 492 ::
 493
 494     sage: K1 = random_cone(max_ambient_dim=8,
 495     ....:                  strictly_convex=False,
 496     ....:                  solid=False)
 497     sage: A = random_matrix(QQ, K1.lattice_dim(), algorithm='unimodular')
 498     sage: K2 = Cone( [ A*r for r in K1.rays() ], lattice=K1.lattice())
 499     sage: K1.lyapunov_rank() == K2.lyapunov_rank()
 500     True
 501
 502 The Lyapunov rank of a dual cone should be the same as the original
 503 cone. Check all combinations of parameters::
 504
 505     sage: set_random_seed()
 506     sage: K = random_cone(max_ambient_dim=8,
 507     ....:                 strictly_convex=False,
 508     ....:                 solid=False)
 509     sage: K.lyapunov_rank() == K.dual().lyapunov_rank()
 510     True
 511
 512 ::
 513
 514     sage: set_random_seed()
 515     sage: K = random_cone(max_ambient_dim=8,
 516     ....:                 strictly_convex=False,
 517     ....:                 solid=True)
 518     sage: K.lyapunov_rank() == K.dual().lyapunov_rank()
 519     True
 520
 521 ::
 522
 523     sage: set_random_seed()
 524     sage: K = random_cone(max_ambient_dim=8,
 525     ....:                 strictly_convex=True,
 526     ....:                 solid=False)
 527     sage: K.lyapunov_rank() == K.dual().lyapunov_rank()
 528     True
 529
 530 ::
 531
 532     sage: set_random_seed()
 533     sage: K = random_cone(max_ambient_dim=8,
 534     ....:                 strictly_convex=True,
 535     ....:                 solid=True)
 536     sage: K.lyapunov_rank() == K.dual().lyapunov_rank()
 537     True
 538
 539 The Lyapunov rank of a cone ``K`` is the dimension of
 540 ``K.lyapunov_like_basis()``. Check all combinations of parameters::
 541
 542     sage: set_random_seed()
 543     sage: K = random_cone(max_ambient_dim=8,
 544     ....:                 strictly_convex=True,
 545     ....:                 solid=True)
 546     sage: K.lyapunov_rank() == len(K.lyapunov_like_basis())
 547     True
 548
 549 ::
 550
 551     sage: set_random_seed()
 552     sage: K = random_cone(max_ambient_dim=8,
 553     ....:                 strictly_convex=True,
 554     ....:                 solid=False)
 555     sage: K.lyapunov_rank() == len(K.lyapunov_like_basis())
 556     True
 557
 558 ::
 559
 560     sage: set_random_seed()
 561     sage: K = random_cone(max_ambient_dim=8,
 562     ....:                 strictly_convex=False,
 563     ....:                 solid=True)
 564     sage: K.lyapunov_rank() == len(K.lyapunov_like_basis())
 565     True
 566
 567 ::
 568
 569     sage: set_random_seed()
 570     sage: K = random_cone(max_ambient_dim=8,
 571     ....:                 strictly_convex=False,
 572     ....:                 solid=False)
 573     sage: K.lyapunov_rank() == len(K.lyapunov_like_basis())
 574     True
 575
 576 """