test/randomgen.py

   1 """
   2 Random thing generators used in the rest of the test suite.
   3 """
   4 from random import randint, uniform
   5
   6 from math import sqrt
   7 from cvxopt import matrix
   8 from dunshire.cones import NonnegativeOrthant, IceCream
   9 from dunshire.games import SymmetricLinearGame
  10 from dunshire.matrices import (append_col, append_row, identity)
  11
  12 MAX_COND = 125
  13 """
  14 The maximum condition number of a randomly-generated game. When the
  15 condition number of the games gets too high, we start to see
  16 :class:`PoorScalingException` being thrown. There's no science to
  17 choosing the upper bound -- it got lowered until those exceptions
  18 stopped popping up. It's at ``125`` because ``129`` doesn't work.
  19 """
  20
  21 RANDOM_MAX = 10
  22 """
  23 When generating random real numbers or integers, this is used as the
  24 largest allowed magnitude. It keeps our condition numbers down and other
  25 properties within reason.
  26 """
  27
  28 def random_scalar():
  29     """
  30     Generate a random scalar.
  31
  32     Returns
  33     -------
  34
  35     float
  36         A random real number between negative and positive
  37         :const:`RANDOM_MAX`, inclusive.
  38
  39     Examples
  40     --------
  41
  42         >>> abs(random_scalar()) <= RANDOM_MAX
  43         True
  44
  45     """
  46     return uniform(-RANDOM_MAX, RANDOM_MAX)
  47
  48
  49 def random_nn_scalar():
  50     """
  51     Generate a random nonnegative scalar.
  52
  53     Returns
  54     -------
  55
  56     float
  57         A random nonnegative real number between zero and
  58         :const:`RANDOM_MAX`, inclusive.
  59
  60     Examples
  61     --------
  62
  63         >>> 0 <= random_nn_scalar() <= RANDOM_MAX
  64         True
  65
  66     """
  67     return abs(random_scalar())
  68
  69
  70 def random_natural():
  71     """
  72     Generate a random natural number.
  73
  74     Returns
  75     -------
  76
  77     int
  78         A random natural number between ``1`` and :const:`RANDOM_MAX`,
  79         inclusive.
  80
  81     Examples
  82     --------
  83
  84         >>> 1 <= random_natural() <= RANDOM_MAX
  85         True
  86
  87     """
  88     return randint(1, RANDOM_MAX)
  89
  90
  91 def random_matrix(row_count, column_count=None):
  92     """
  93     Generate a random matrix.
  94
  95     Parameters
  96     ----------
  97
  98     row_count : int
  99         The number of rows you want in the returned matrix.
 100
 101     column_count: int
 102         The number of columns you want in the returned matrix (default:
 103         the same as ``row_count``).
 104
 105     Returns
 106     -------
 107
 108     matrix
 109         A new matrix whose entries are random floats chosen uniformly
 110         between negative and positive :const:`RANDOM_MAX`.
 111
 112     Examples
 113     --------
 114
 115         >>> A = random_matrix(3)
 116         >>> A.size
 117         (3, 3)
 118
 119         >>> A = random_matrix(3,2)
 120         >>> A.size
 121         (3, 2)
 122
 123     """
 124     if column_count is None:
 125         column_count = row_count
 126
 127     entries = [random_scalar() for _ in range(row_count*column_count)]
 128     return matrix(entries, (row_count, column_count))
 129
 130
 131 def random_nonnegative_matrix(row_count, column_count=None):
 132     """
 133     Generate a random matrix with nonnegative entries.
 134
 135     Parameters
 136     ----------
 137
 138     row_count : int
 139         The number of rows you want in the returned matrix.
 140
 141     column_count : int
 142         The number of columns you want in the returned matrix (default:
 143         the same as ``row_count``).
 144
 145     Returns
 146     -------
 147
 148     matrix
 149         A new matrix whose entries are chosen by :func:`random_nn_scalar`.
 150
 151     Examples
 152     --------
 153
 154         >>> A = random_nonnegative_matrix(3)
 155         >>> A.size
 156         (3, 3)
 157         >>> all([entry >= 0 for entry in A])
 158         True
 159
 160         >>> A = random_nonnegative_matrix(3,2)
 161         >>> A.size
 162         (3, 2)
 163         >>> all([entry >= 0 for entry in A])
 164         True
 165
 166     """
 167     if column_count is None:
 168         column_count = row_count
 169
 170     entries = [random_nn_scalar() for _ in range(row_count*column_count)]
 171     return matrix(entries, (row_count, column_count))
 172
 173
 174 def random_diagonal_matrix(dims):
 175     """
 176     Generate a random square matrix with zero off-diagonal entries.
 177
 178     These matrices are Lyapunov-like on the nonnegative orthant, as is
 179     fairly easy to see.
 180
 181     Parameters
 182     ----------
 183
 184     dims : int
 185         The number of rows/columns you want in the returned matrix.
 186
 187     Returns
 188     -------
 189
 190     matrix
 191         A new matrix whose diagonal entries are random floats chosen
 192         using :func:`random_scalar` and whose off-diagonal entries are
 193         zero.
 194
 195     Examples
 196     --------
 197
 198         >>> A = random_diagonal_matrix(3)
 199         >>> A.size
 200         (3, 3)
 201         >>> A[0,1] == A[0,2] == A[1,0] == A[2,0] == A[1,2] == A[2,1] == 0
 202         True
 203
 204     """
 205     return matrix([[random_scalar()*int(i == j)
 206                     for i in range(dims)]
 207                    for j in range(dims)])
 208
 209
 210 def random_skew_symmetric_matrix(dims):
 211     """
 212     Generate a random skew-symmetrix matrix.
 213
 214     Parameters
 215     ----------
 216
 217     dims : int
 218         The number of rows/columns you want in the returned matrix.
 219
 220     Returns
 221     -------
 222
 223     matrix
 224         A new skew-matrix whose strictly above-diagonal entries are
 225         random floats chosen with :func:`random_scalar`.
 226
 227     Examples
 228     --------
 229
 230         >>> A = random_skew_symmetric_matrix(3)
 231         >>> A.size
 232         (3, 3)
 233
 234         >>> from dunshire.options import ABS_TOL
 235         >>> from dunshire.matrices import norm
 236         >>> A = random_skew_symmetric_matrix(random_natural())
 237         >>> norm(A + A.trans()) < ABS_TOL
 238         True
 239
 240     """
 241     strict_ut = [[random_scalar()*int(i < j)
 242                   for i in range(dims)]
 243                  for j in range(dims)]
 244
 245     strict_ut = matrix(strict_ut, (dims, dims))
 246     return strict_ut - strict_ut.trans()
 247
 248
 249 def random_lyapunov_like_icecream(dims):
 250     r"""
 251     Generate a random matrix Lyapunov-like on the ice-cream cone.
 252
 253     The form of these matrices is cited in Gowda and Tao
 254     [GowdaTao]_. The scalar ``a`` and the vector ``b`` (using their
 255     notation) are easy to generate. The submatrix ``D`` is a little
 256     trickier, but it can be found noticing that :math:`C + C^{T} = 0`
 257     for a skew-symmetric matrix :math:`C` implying that :math:`C + C^{T}
 258     + \left(2a\right)I = \left(2a\right)I`. Thus we can stick an
 259     :math:`aI` with each of :math:`C,C^{T}` and let those be our
 260     :math:`D,D^{T}`.
 261
 262     Parameters
 263     ----------
 264
 265     dims : int
 266         The dimension of the ice-cream cone (not of the matrix you want!)
 267         on which the returned matrix should be Lyapunov-like.
 268
 269     Returns
 270     -------
 271
 272     matrix
 273         A new matrix, Lyapunov-like on the ice-cream cone in ``dims``
 274         dimensions, whose free entries are random floats chosen uniformly
 275         between negative and positive :const:`RANDOM_MAX`.
 276
 277     References
 278     ----------
 279
 280     .. [GowdaTao] M. S. Gowda and J. Tao. On the bilinearity rank of a
 281        proper cone and Lyapunov-like transformations. Mathematical
 282        Programming, 147:155-170, 2014.
 283
 284     Examples
 285     --------
 286
 287         >>> L = random_lyapunov_like_icecream(3)
 288         >>> L.size
 289         (3, 3)
 290
 291         >>> from dunshire.options import ABS_TOL
 292         >>> from dunshire.matrices import inner_product
 293         >>> x = matrix([1,1,0])
 294         >>> s = matrix([1,-1,0])
 295         >>> abs(inner_product(L*x, s)) < ABS_TOL
 296         True
 297
 298     """
 299     a = matrix([random_scalar()], (1, 1))
 300     b = matrix([random_scalar() for _ in range(dims-1)], (dims-1, 1))
 301     D = random_skew_symmetric_matrix(dims-1) + a*identity(dims-1)
 302     row1 = append_col(a, b.trans())
 303     row2 = append_col(b, D)
 304     return append_row(row1, row2)
 305
 306
 307 def random_orthant_game():
 308     """
 309     Generate a random game over the nonnegative orthant.
 310
 311     We generate each of ``L``, ``K``, ``e1``, and ``e2`` randomly within
 312     the constraints of the nonnegative orthant, and then construct a
 313     game from them. The process is repeated until we generate a game with
 314     a condition number under :const:`MAX_COND`.
 315
 316     Returns
 317     -------
 318
 319     SymmetricLinearGame
 320         A random game over some nonnegative orthant.
 321
 322     Examples
 323     --------
 324
 325         >>> random_orthant_game()
 326         <dunshire.games.SymmetricLinearGame object at 0x...>
 327
 328     """
 329     ambient_dim = random_natural() + 1
 330     K = NonnegativeOrthant(ambient_dim)
 331     e1 = [0.1 + random_nn_scalar() for _ in range(K.dimension())]
 332     e2 = [0.1 + random_nn_scalar() for _ in range(K.dimension())]
 333     L = random_matrix(K.dimension())
 334     G = SymmetricLinearGame(L, K, e1, e2)
 335
 336     if G.condition() <= MAX_COND:
 337         return G
 338     else:
 339         return random_orthant_game()
 340
 341
 342 def random_icecream_game():
 343     """
 344     Generate a random game over the ice-cream cone.
 345
 346     We generate each of ``L``, ``K``, ``e1``, and ``e2`` randomly within
 347     the constraints of the ice-cream cone, and then construct a game
 348     from them. The process is repeated until we generate a game with a
 349     condition number under :const:`MAX_COND`.
 350
 351     Returns
 352     -------
 353
 354     SymmetricLinearGame
 355         A random game over some ice-cream cone.
 356
 357     Examples
 358     --------
 359
 360         >>> random_icecream_game()
 361         <dunshire.games.SymmetricLinearGame object at 0x...>
 362
 363     """
 364     # Use a minimum dimension of two to avoid divide-by-zero in
 365     # the fudge factor we make up later.
 366     ambient_dim = random_natural() + 1
 367     K = IceCream(ambient_dim)
 368     e1 = [1] # Set the "height" of e1 to one
 369     e2 = [1] # And the same for e2
 370
 371     # If we choose the rest of the components of e1,e2 randomly
 372     # between 0 and 1, then the largest the squared norm of the
 373     # non-height part of e1,e2 could be is the 1*(dim(K) - 1). We
 374     # need to make it less than one (the height of the cone) so
 375     # that the whole thing is in the cone. The norm of the
 376     # non-height part is sqrt(dim(K) - 1), and we can divide by
 377     # twice that.
 378     fudge_factor = 1.0 / (2.0*sqrt(K.dimension() - 1.0))
 379     e1 += [fudge_factor*uniform(0, 1) for _ in range(K.dimension() - 1)]
 380     e2 += [fudge_factor*uniform(0, 1) for _ in range(K.dimension() - 1)]
 381     L = random_matrix(K.dimension())
 382     G = SymmetricLinearGame(L, K, e1, e2)
 383
 384     if G.condition() <= MAX_COND:
 385         return G
 386     else:
 387         return random_icecream_game()
 388
 389
 390 def random_game():
 391     """
 392     Return a random game.
 393
 394     One of the functions,
 395
 396       1. :func:`random_orthant_game`
 397       2. :func:`random_icecream_game`
 398
 399     is chosen at random and used to generate a random game.
 400
 401     Returns
 402     -------
 403
 404     SymmetricLinearGame
 405         A random game.
 406
 407     Examples
 408     --------
 409
 410         >>> random_game()
 411         <dunshire.games.SymmetricLinearGame object at 0x...>
 412
 413     """
 414     cone_type = randint(0,1)
 415     if cone_type == 0:
 416         return random_orthant_game()
 417     elif cone_type == 1:
 418         return random_icecream_game()
 419
 420 def random_ll_orthant_game():
 421     """
 422     Return a random Lyapunov game over some nonnegative orthant.
 423
 424     We first construct a :func:`random_orthant_game` and then modify it
 425     to have a :func:`random_diagonal_matrix` as its operator. Such
 426     things are Lyapunov-like on the nonnegative orthant. That process is
 427     repeated until the condition number of the resulting game is within
 428     :const:`MAX_COND`.
 429
 430     Returns
 431     -------
 432
 433     SymmetricLinearGame
 434
 435         A random game over some nonnegative orthant whose
 436         :meth:`dunshire.games.SymmetricLinearGame.payoff` method is
 437         based on a Lyapunov-like
 438         :meth:`dunshire.games.SymmetricLinearGame.L` operator.
 439
 440     Examples
 441     --------
 442
 443         >>> random_ll_orthant_game()
 444         <dunshire.games.SymmetricLinearGame object at 0x...>
 445
 446     """
 447     G = random_orthant_game()
 448     L = random_diagonal_matrix(G.dimension())
 449
 450     # Replace the totally-random ``L`` with random Lyapunov-like one.
 451     G = SymmetricLinearGame(L, G.K(), G.e1(), G.e2())
 452
 453     while G.condition() > MAX_COND:
 454         # Try again until the condition number is satisfactory.
 455         G = random_orthant_game()
 456         L = random_diagonal_matrix(G.dimension())
 457         G = SymmetricLinearGame(L, G.K(), G.e1(), G.e2())
 458
 459     return G
 460
 461
 462 def random_ll_icecream_game():
 463     """
 464     Return a random Lyapunov game over some ice-cream cone.
 465
 466     We first construct a :func:`random_icecream_game` and then modify it
 467     to have a :func:`random_lyapunov_like_icecream` operator. That
 468     process is repeated until the condition number of the resulting game
 469     is within :const:`MAX_COND`.
 470
 471     Returns
 472     -------
 473
 474     SymmetricLinearGame
 475         A random game over some ice-cream cone whose
 476         :meth:`dunshire.games.SymmetricLinearGame.payoff` method
 477         is based on a Lyapunov-like
 478         :meth:`dunshire.games.SymmetricLinearGame.L` operator.
 479
 480     Examples
 481     --------
 482
 483         >>> random_ll_icecream_game()
 484         <dunshire.games.SymmetricLinearGame object at 0x...>
 485
 486     """
 487     G = random_icecream_game()
 488     L = random_lyapunov_like_icecream(G.dimension())
 489
 490     # Replace the totally-random ``L`` with random Lyapunov-like one.
 491     G = SymmetricLinearGame(L, G.K(), G.e1(), G.e2())
 492
 493     while G.condition() > MAX_COND:
 494         # Try again until the condition number is satisfactory.
 495         G = random_icecream_game()
 496         L = random_lyapunov_like_icecream(G.dimension())
 497         G = SymmetricLinearGame(L, G.K(), G.e1(), G.e2())
 498
 499     return G
 500
 501
 502 def random_positive_orthant_game():
 503     """
 504     Return a random game over the nonnegative orthant with a positive
 505     operator.
 506
 507     We first construct a :func:`random_orthant_game` and then modify it
 508     to have a :func:`random_nonnegative_matrix` as its operator. That
 509     process is repeated until the condition number of the resulting game
 510     is within :const:`MAX_COND`.
 511
 512     Returns
 513     -------
 514
 515     SymmetricLinearGame
 516         A random game over some nonnegative orthant whose
 517         :meth:`dunshire.games.SymmetricLinearGame.payoff` method
 518         is based on a positive
 519         :meth:`dunshire.games.SymmetricLinearGame.L` operator.
 520
 521     Examples
 522     --------
 523
 524         >>> random_positive_orthant_game()
 525         <dunshire.games.SymmetricLinearGame object at 0x...>
 526
 527     """
 528
 529     G = random_orthant_game()
 530     L = random_nonnegative_matrix(G.dimension())
 531
 532     # Replace the totally-random ``L`` with the random nonnegative one.
 533     G = SymmetricLinearGame(L, G.K(), G.e1(), G.e2())
 534
 535     while G.condition() > MAX_COND:
 536         # Try again until the condition number is satisfactory.
 537         G = random_orthant_game()
 538         L = random_nonnegative_matrix(G.dimension())
 539         G = SymmetricLinearGame(L, G.K(), G.e1(), G.e2())
 540
 541     return G
 542
 543
 544 def random_nn_scaling(G):
 545     """
 546     Scale the given game by a random nonnegative amount.
 547
 548     We re-attempt the scaling with a new random number until the
 549     resulting scaled game has an acceptable condition number.
 550
 551     Parameters
 552     ----------
 553
 554     G : SymmetricLinearGame
 555         The game that you would like to scale.
 556
 557     Returns
 558     -------
 559     (float, SymmetricLinearGame)
 560         A pair containing the both the scaling factor and the new scaled game.
 561
 562     Examples
 563     --------
 564
 565         >>> from dunshire.matrices import norm
 566         >>> from dunshire.options import ABS_TOL
 567         >>> G = random_orthant_game()
 568         >>> (alpha, H) = random_nn_scaling(G)
 569         >>> alpha >= 0
 570         True
 571         >>> G.K() == H.K()
 572         True
 573         >>> norm(G.e1() - H.e1()) < ABS_TOL
 574         True
 575         >>> norm(G.e2() - H.e2()) < ABS_TOL
 576         True
 577
 578     """
 579     alpha = random_nn_scalar()
 580     H = SymmetricLinearGame(alpha*G.L().trans(), G.K(), G.e1(), G.e2())
 581
 582     while H.condition() > MAX_COND:
 583         # Loop until the condition number of H doesn't suck.
 584         alpha = random_nn_scalar()
 585         H = SymmetricLinearGame(alpha*G.L().trans(), G.K(), G.e1(), G.e2())
 586
 587     return (alpha, H)
 588
 589
 590 def random_translation(G):
 591     """
 592     Translate the given game by a random amount.
 593
 594     We re-attempt the translation with new random scalars until the
 595     resulting translated game has an acceptable condition number.
 596
 597     Parameters
 598     ----------
 599
 600     G : SymmetricLinearGame
 601         The game that you would like to translate.
 602
 603     Returns
 604     -------
 605     (float, SymmetricLinearGame)
 606         A pair containing the both the translation distance and the new
 607         scaled game.
 608
 609     Examples
 610     --------
 611
 612         >>> from dunshire.matrices import norm
 613         >>> from dunshire.options import ABS_TOL
 614         >>> G = random_orthant_game()
 615         >>> (alpha, H) = random_translation(G)
 616         >>> G.K() == H.K()
 617         True
 618         >>> norm(G.e1() - H.e1()) < ABS_TOL
 619         True
 620         >>> norm(G.e2() - H.e2()) < ABS_TOL
 621         True
 622
 623     """
 624     alpha = random_scalar()
 625     tensor_prod = G.e1() * G.e2().trans()
 626     M = G.L() + alpha*tensor_prod
 627
 628     H = SymmetricLinearGame(M.trans(), G.K(), G.e1(), G.e2())
 629     while H.condition() > MAX_COND:
 630         # Loop until the condition number of H doesn't suck.
 631         alpha = random_scalar()
 632         M = G.L() + alpha*tensor_prod
 633         H = SymmetricLinearGame(M.trans(), G.K(), G.e1(), G.e2())
 634
 635     return (alpha, H)