dunshire/games.py

   1 """
   2 Symmetric linear games and their solutions.
   3
   4 This module contains the main :class:`SymmetricLinearGame` class that
   5 knows how to solve a linear game.
   6 """
   7
   8 from cvxopt import matrix, printing, solvers
   9 from .cones import CartesianProduct
  10 from .errors import GameUnsolvableException, PoorScalingException
  11 from .matrices import append_col, append_row, identity
  12 from . import options
  13
  14 printing.options['dformat'] = options.FLOAT_FORMAT
  15 solvers.options['show_progress'] = options.VERBOSE
  16
  17
  18 class Solution:
  19     """
  20     A representation of the solution of a linear game. It should contain
  21     the value of the game, and both players' strategies.
  22
  23     Examples
  24     --------
  25
  26         >>> print(Solution(10, matrix([1,2]), matrix([3,4])))
  27         Game value: 10.0000000
  28         Player 1 optimal:
  29           [ 1]
  30           [ 2]
  31         Player 2 optimal:
  32           [ 3]
  33           [ 4]
  34
  35     """
  36     def __init__(self, game_value, p1_optimal, p2_optimal):
  37         """
  38         Create a new Solution object from a game value and two optimal
  39         strategies for the players.
  40         """
  41         self._game_value = game_value
  42         self._player1_optimal = p1_optimal
  43         self._player2_optimal = p2_optimal
  44
  45     def __str__(self):
  46         """
  47         Return a string describing the solution of a linear game.
  48
  49         The three data that are described are,
  50
  51           * The value of the game.
  52           * The optimal strategy of player one.
  53           * The optimal strategy of player two.
  54
  55         The two optimal strategy vectors are indented by two spaces.
  56         """
  57         tpl = 'Game value: {:.7f}\n' \
  58               'Player 1 optimal:{:s}\n' \
  59               'Player 2 optimal:{:s}'
  60
  61         p1_str = '\n{!s}'.format(self.player1_optimal())
  62         p1_str = '\n  '.join(p1_str.splitlines())
  63         p2_str = '\n{!s}'.format(self.player2_optimal())
  64         p2_str = '\n  '.join(p2_str.splitlines())
  65
  66         return tpl.format(self.game_value(), p1_str, p2_str)
  67
  68
  69     def game_value(self):
  70         """
  71         Return the game value for this solution.
  72
  73         Examples
  74         --------
  75
  76             >>> s = Solution(10, matrix([1,2]), matrix([3,4]))
  77             >>> s.game_value()
  78             10
  79
  80         """
  81         return self._game_value
  82
  83
  84     def player1_optimal(self):
  85         """
  86         Return player one's optimal strategy in this solution.
  87
  88         Examples
  89         --------
  90
  91             >>> s = Solution(10, matrix([1,2]), matrix([3,4]))
  92             >>> print(s.player1_optimal())
  93             [ 1]
  94             [ 2]
  95             <BLANKLINE>
  96
  97         """
  98         return self._player1_optimal
  99
 100
 101     def player2_optimal(self):
 102         """
 103         Return player two's optimal strategy in this solution.
 104
 105         Examples
 106         --------
 107
 108             >>> s = Solution(10, matrix([1,2]), matrix([3,4]))
 109             >>> print(s.player2_optimal())
 110             [ 3]
 111             [ 4]
 112             <BLANKLINE>
 113
 114         """
 115         return self._player2_optimal
 116
 117
 118 class SymmetricLinearGame:
 119     r"""
 120     A representation of a symmetric linear game.
 121
 122     The data for a symmetric linear game are,
 123
 124       * A "payoff" operator ``L``.
 125       * A symmetric cone ``K``.
 126       * Two points ``e1`` and ``e2`` in the interior of ``K``.
 127
 128     The ambient space is assumed to be the span of ``K``.
 129
 130     With those data understood, the game is played as follows. Players
 131     one and two choose points :math:`x` and :math:`y` respectively, from
 132     their respective strategy sets,
 133
 134     .. math::
 135         \begin{aligned}
 136             \Delta_{1}
 137             &=
 138             \left\{
 139               x \in K \ \middle|\ \left\langle x, e_{2} \right\rangle = 1
 140             \right\}\\
 141             \Delta_{2}
 142             &=
 143             \left\{
 144               y \in K \ \middle|\ \left\langle y, e_{1} \right\rangle = 1
 145             \right\}.
 146         \end{aligned}
 147
 148     Afterwards, a "payout" is computed as :math:`\left\langle
 149     L\left(x\right), y \right\rangle` and is paid to player one out of
 150     player two's pocket. The game is therefore zero sum, and we suppose
 151     that player one would like to guarantee himself the largest minimum
 152     payout possible. That is, player one wishes to,
 153
 154     .. math::
 155         \begin{aligned}
 156             \text{maximize }
 157               &\underset{y \in \Delta_{2}}{\min}\left(
 158                 \left\langle L\left(x\right), y \right\rangle
 159               \right)\\
 160             \text{subject to } & x \in \Delta_{1}.
 161         \end{aligned}
 162
 163     Player two has the simultaneous goal to,
 164
 165     .. math::
 166         \begin{aligned}
 167             \text{minimize }
 168               &\underset{x \in \Delta_{1}}{\max}\left(
 169                 \left\langle L\left(x\right), y \right\rangle
 170               \right)\\
 171             \text{subject to } & y \in \Delta_{2}.
 172         \end{aligned}
 173
 174     These goals obviously conflict (the game is zero sum), but an
 175     existence theorem guarantees at least one optimal min-max solution
 176     from which neither player would like to deviate. This class is
 177     able to find such a solution.
 178
 179     Parameters
 180     ----------
 181
 182     L : list of list of float
 183         A matrix represented as a list of ROWS. This representation
 184         agrees with (for example) SageMath and NumPy, but not with CVXOPT
 185         (whose matrix constructor accepts a list of columns).
 186
 187     K : :class:`SymmetricCone`
 188         The symmetric cone instance over which the game is played.
 189
 190     e1 : iterable float
 191         The interior point of ``K`` belonging to player one; it
 192         can be of any iterable type having the correct length.
 193
 194     e2 : iterable float
 195         The interior point of ``K`` belonging to player two; it
 196         can be of any enumerable type having the correct length.
 197
 198     Raises
 199     ------
 200
 201     ValueError
 202         If either ``e1`` or ``e2`` lie outside of the cone ``K``.
 203
 204     Examples
 205     --------
 206
 207         >>> from dunshire import *
 208         >>> K = NonnegativeOrthant(3)
 209         >>> L = [[1,-5,-15],[-1,2,-3],[-12,-15,1]]
 210         >>> e1 = [1,1,1]
 211         >>> e2 = [1,2,3]
 212         >>> SLG = SymmetricLinearGame(L, K, e1, e2)
 213         >>> print(SLG)
 214         The linear game (L, K, e1, e2) where
 215           L = [  1  -5 -15]
 216               [ -1   2  -3]
 217               [-12 -15   1],
 218           K = Nonnegative orthant in the real 3-space,
 219           e1 = [ 1]
 220                [ 1]
 221                [ 1],
 222           e2 = [ 1]
 223                [ 2]
 224                [ 3].
 225
 226     Lists can (and probably should) be used for every argument::
 227
 228         >>> from dunshire import *
 229         >>> K = NonnegativeOrthant(2)
 230         >>> L = [[1,0],[0,1]]
 231         >>> e1 = [1,1]
 232         >>> e2 = [1,1]
 233         >>> G = SymmetricLinearGame(L, K, e1, e2)
 234         >>> print(G)
 235         The linear game (L, K, e1, e2) where
 236           L = [ 1  0]
 237               [ 0  1],
 238           K = Nonnegative orthant in the real 2-space,
 239           e1 = [ 1]
 240                [ 1],
 241           e2 = [ 1]
 242                [ 1].
 243
 244     The points ``e1`` and ``e2`` can also be passed as some other
 245     enumerable type (of the correct length) without much harm, since
 246     there is no row/column ambiguity::
 247
 248         >>> import cvxopt
 249         >>> import numpy
 250         >>> from dunshire import *
 251         >>> K = NonnegativeOrthant(2)
 252         >>> L = [[1,0],[0,1]]
 253         >>> e1 = cvxopt.matrix([1,1])
 254         >>> e2 = numpy.matrix([1,1])
 255         >>> G = SymmetricLinearGame(L, K, e1, e2)
 256         >>> print(G)
 257         The linear game (L, K, e1, e2) where
 258           L = [ 1  0]
 259               [ 0  1],
 260           K = Nonnegative orthant in the real 2-space,
 261           e1 = [ 1]
 262                [ 1],
 263           e2 = [ 1]
 264                [ 1].
 265
 266     However, ``L`` will always be intepreted as a list of rows, even
 267     if it is passed as a :class:`cvxopt.base.matrix` which is
 268     otherwise indexed by columns::
 269
 270         >>> import cvxopt
 271         >>> from dunshire import *
 272         >>> K = NonnegativeOrthant(2)
 273         >>> L = [[1,2],[3,4]]
 274         >>> e1 = [1,1]
 275         >>> e2 = e1
 276         >>> G = SymmetricLinearGame(L, K, e1, e2)
 277         >>> print(G)
 278         The linear game (L, K, e1, e2) where
 279           L = [ 1  2]
 280               [ 3  4],
 281           K = Nonnegative orthant in the real 2-space,
 282           e1 = [ 1]
 283                [ 1],
 284           e2 = [ 1]
 285                [ 1].
 286         >>> L = cvxopt.matrix(L)
 287         >>> print(L)
 288         [ 1  3]
 289         [ 2  4]
 290         <BLANKLINE>
 291         >>> G = SymmetricLinearGame(L, K, e1, e2)
 292         >>> print(G)
 293         The linear game (L, K, e1, e2) where
 294           L = [ 1  2]
 295               [ 3  4],
 296           K = Nonnegative orthant in the real 2-space,
 297           e1 = [ 1]
 298                [ 1],
 299           e2 = [ 1]
 300                [ 1].
 301
 302     """
 303     def __init__(self, L, K, e1, e2):
 304         """
 305         Create a new SymmetricLinearGame object.
 306         """
 307         self._K = K
 308         self._e1 = matrix(e1, (K.dimension(), 1))
 309         self._e2 = matrix(e2, (K.dimension(), 1))
 310
 311         # Our input ``L`` is indexed by rows but CVXOPT matrices are
 312         # indexed by columns, so we need to transpose the input before
 313         # feeding it to CVXOPT.
 314         self._L = matrix(L, (K.dimension(), K.dimension())).trans()
 315
 316         if not self._e1 in K:
 317             raise ValueError('the point e1 must lie in the interior of K')
 318
 319         if not self._e2 in K:
 320             raise ValueError('the point e2 must lie in the interior of K')
 321
 322     def __str__(self):
 323         """
 324         Return a string representation of this game.
 325         """
 326         tpl = 'The linear game (L, K, e1, e2) where\n' \
 327               '  L = {:s},\n' \
 328               '  K = {!s},\n' \
 329               '  e1 = {:s},\n' \
 330               '  e2 = {:s}.'
 331         indented_L = '\n      '.join(str(self._L).splitlines())
 332         indented_e1 = '\n       '.join(str(self._e1).splitlines())
 333         indented_e2 = '\n       '.join(str(self._e2).splitlines())
 334         return tpl.format(indented_L, str(self._K), indented_e1, indented_e2)
 335
 336
 337     def solution(self):
 338         """
 339         Solve this linear game and return a :class:`Solution`.
 340
 341         Returns
 342         -------
 343
 344         :class:`Solution`
 345             A :class:`Solution` object describing the game's value and
 346             the optimal strategies of both players.
 347
 348         Raises
 349         ------
 350         GameUnsolvableException
 351             If the game could not be solved (if an optimal solution to its
 352             associated cone program was not found).
 353
 354         PoorScalingException
 355             If the game could not be solved because CVXOPT crashed while
 356             trying to take the square root of a negative number.
 357
 358         Examples
 359         --------
 360
 361         This example is computed in Gowda and Ravindran in the section
 362         "The value of a Z-transformation"::
 363
 364             >>> from dunshire import *
 365             >>> K = NonnegativeOrthant(3)
 366             >>> L = [[1,-5,-15],[-1,2,-3],[-12,-15,1]]
 367             >>> e1 = [1,1,1]
 368             >>> e2 = [1,1,1]
 369             >>> SLG = SymmetricLinearGame(L, K, e1, e2)
 370             >>> print(SLG.solution())
 371             Game value: -6.1724138
 372             Player 1 optimal:
 373               [ 0.5517241]
 374               [-0.0000000]
 375               [ 0.4482759]
 376             Player 2 optimal:
 377               [0.4482759]
 378               [0.0000000]
 379               [0.5517241]
 380
 381         The value of the following game can be computed using the fact
 382         that the identity is invertible::
 383
 384             >>> from dunshire import *
 385             >>> K = NonnegativeOrthant(3)
 386             >>> L = [[1,0,0],[0,1,0],[0,0,1]]
 387             >>> e1 = [1,2,3]
 388             >>> e2 = [4,5,6]
 389             >>> SLG = SymmetricLinearGame(L, K, e1, e2)
 390             >>> print(SLG.solution())
 391             Game value: 0.0312500
 392             Player 1 optimal:
 393               [0.0312500]
 394               [0.0625000]
 395               [0.0937500]
 396             Player 2 optimal:
 397               [0.1250000]
 398               [0.1562500]
 399               [0.1875000]
 400
 401         """
 402         # The cone "C" that appears in the statement of the CVXOPT
 403         # conelp program.
 404         C = CartesianProduct(self._K, self._K)
 405
 406         # The column vector "b" that appears on the right-hand side of
 407         # Ax = b in the statement of the CVXOPT conelp program.
 408         b = matrix([1], tc='d')
 409
 410         # A column of zeros that fits K.
 411         zero = matrix(0, (self._K.dimension(), 1), tc='d')
 412
 413         # The column vector "h" that appears on the right-hand side of
 414         # Gx + s = h in the statement of the CVXOPT conelp program.
 415         h = matrix([zero, zero])
 416
 417         # The column vector "c" that appears in the objective function
 418         # value <c,x> in the statement of the CVXOPT conelp program.
 419         c = matrix([-1, zero])
 420
 421         # The matrix "G" that appears on the left-hand side of Gx + s = h
 422         # in the statement of the CVXOPT conelp program.
 423         G = append_row(append_col(zero, -identity(self._K.dimension())),
 424                        append_col(self._e1, -self._L))
 425
 426         # The matrix "A" that appears on the right-hand side of Ax = b
 427         # in the statement of the CVXOPT conelp program.
 428         A = matrix([0, self._e2], (1, self._K.dimension() + 1), 'd')
 429
 430         # Actually solve the thing and obtain a dictionary describing
 431         # what happened.
 432         try:
 433             soln_dict = solvers.conelp(c, G, h, C.cvxopt_dims(), A, b)
 434         except ValueError as e:
 435             if str(e) == 'math domain error':
 436                 # Oops, CVXOPT tried to take the square root of a
 437                 # negative number. Report some details about the game
 438                 # rather than just the underlying CVXOPT crash.
 439                 raise PoorScalingException(self)
 440             else:
 441                 raise e
 442
 443         # The optimal strategies are named ``p`` and ``q`` in the
 444         # background documentation, and we need to extract them from
 445         # the CVXOPT ``x`` and ``z`` variables. The objective values
 446         # :math:`nu` and :math:`omega` can also be found in the CVXOPT
 447         # ``x`` and ``y`` variables; however, they're stored
 448         # conveniently as separate entries in the solution dictionary.
 449         p1_value = -soln_dict['primal objective']
 450         p2_value = -soln_dict['dual objective']
 451         p1_optimal = soln_dict['x'][1:]
 452         p2_optimal = soln_dict['z'][self._K.dimension():]
 453
 454         # The "status" field contains "optimal" if everything went
 455         # according to plan. Other possible values are "primal
 456         # infeasible", "dual infeasible", "unknown", all of which mean
 457         # we didn't get a solution. The "infeasible" ones are the
 458         # worst, since they indicate that CVXOPT is convinced the
 459         # problem is infeasible (and that cannot happen).
 460         if soln_dict['status'] in ['primal infeasible', 'dual infeasible']:
 461             raise GameUnsolvableException(self, soln_dict)
 462         elif soln_dict['status'] == 'unknown':
 463             # When we get a status of "unknown", we may still be able
 464             # to salvage a solution out of the returned
 465             # dictionary. Often this is the result of numerical
 466             # difficulty and we can simply check that the primal/dual
 467             # objectives match (within a tolerance) and that the
 468             # primal/dual optimal solutions are within the cone (to a
 469             # tolerance as well).
 470             if abs(p1_value - p2_value) > options.ABS_TOL:
 471                 raise GameUnsolvableException(self, soln_dict)
 472             if (p1_optimal not in self._K) or (p2_optimal not in self._K):
 473                 raise GameUnsolvableException(self, soln_dict)
 474
 475         return Solution(p1_value, p1_optimal, p2_optimal)
 476
 477
 478     def dual(self):
 479         r"""
 480         Return the dual game to this game.
 481
 482         If :math:`G = \left(L,K,e_{1},e_{2}\right)` is a linear game,
 483         then its dual is :math:`G^{*} =
 484         \left(L^{*},K^{*},e_{2},e_{1}\right)`. However, since this cone
 485         is symmetric, :math:`K^{*} = K`.
 486
 487         Examples
 488         --------
 489
 490             >>> from dunshire import *
 491             >>> K = NonnegativeOrthant(3)
 492             >>> L = [[1,-5,-15],[-1,2,-3],[-12,-15,1]]
 493             >>> e1 = [1,1,1]
 494             >>> e2 = [1,2,3]
 495             >>> SLG = SymmetricLinearGame(L, K, e1, e2)
 496             >>> print(SLG.dual())
 497             The linear game (L, K, e1, e2) where
 498               L = [  1  -1 -12]
 499                   [ -5   2 -15]
 500                   [-15  -3   1],
 501               K = Nonnegative orthant in the real 3-space,
 502               e1 = [ 1]
 503                    [ 2]
 504                    [ 3],
 505               e2 = [ 1]
 506                    [ 1]
 507                    [ 1].
 508
 509         """
 510         # We pass ``self._L`` right back into the constructor, because
 511         # it will be transposed there. And keep in mind that ``self._K``
 512         # is its own dual.
 513         return SymmetricLinearGame(self._L,
 514                                    self._K,
 515                                    self._e2,
 516                                    self._e1)