src/dunshire/symmetric_linear_game.py

   1 """
   2 Symmetric linear games and their solutions.
   3
   4 This module contains the main SymmetricLinearGame class that knows how
   5 to solve a linear game.
   6 """
   7
   8 # These few are used only for tests.
   9 from math import sqrt
  10 from random import randint, uniform
  11 from unittest import TestCase
  12
  13 # These are mostly actually needed.
  14 from cvxopt import matrix, printing, solvers
  15 from cones import CartesianProduct, IceCream, NonnegativeOrthant
  16 from errors import GameUnsolvableException
  17 from matrices import append_col, append_row, identity, inner_product
  18 import options
  19
  20 printing.options['dformat'] = options.FLOAT_FORMAT
  21 solvers.options['show_progress'] = options.VERBOSE
  22
  23
  24 class Solution:
  25     """
  26     A representation of the solution of a linear game. It should contain
  27     the value of the game, and both players' strategies.
  28     """
  29     def __init__(self, game_value, p1_optimal, p2_optimal):
  30         """
  31         Create a new Solution object from a game value and two optimal
  32         strategies for the players.
  33         """
  34         self._game_value = game_value
  35         self._player1_optimal = p1_optimal
  36         self._player2_optimal = p2_optimal
  37
  38     def __str__(self):
  39         """
  40         Return a string describing the solution of a linear game.
  41
  42         The three data that are described are,
  43
  44           * The value of the game.
  45           * The optimal strategy of player one.
  46           * The optimal strategy of player two.
  47
  48         EXAMPLES:
  49
  50            >>> print(Solution(10, matrix([1,2]), matrix([3,4])))
  51            Game value: 10.0000000
  52            Player 1 optimal:
  53              [ 1]
  54              [ 2]
  55            Player 2 optimal:
  56              [ 3]
  57              [ 4]
  58
  59         """
  60         tpl = 'Game value: {:.7f}\n' \
  61               'Player 1 optimal:{:s}\n' \
  62               'Player 2 optimal:{:s}'
  63
  64         p1_str = '\n{!s}'.format(self.player1_optimal())
  65         p1_str = '\n  '.join(p1_str.splitlines())
  66         p2_str = '\n{!s}'.format(self.player2_optimal())
  67         p2_str = '\n  '.join(p2_str.splitlines())
  68
  69         return tpl.format(self.game_value(), p1_str, p2_str)
  70
  71
  72     def game_value(self):
  73         """
  74         Return the game value for this solution.
  75         """
  76         return self._game_value
  77
  78
  79     def player1_optimal(self):
  80         """
  81         Return player one's optimal strategy in this solution.
  82         """
  83         return self._player1_optimal
  84
  85
  86     def player2_optimal(self):
  87         """
  88         Return player two's optimal strategy in this solution.
  89         """
  90         return self._player2_optimal
  91
  92
  93 class SymmetricLinearGame:
  94     """
  95     A representation of a symmetric linear game.
  96
  97     The data for a linear game are,
  98
  99       * A "payoff" operator ``L``.
 100       * A cone ``K``.
 101       * A point ``e`` in the interior of ``K``.
 102       * A point ``f`` in the interior of the dual of ``K``.
 103
 104     In a symmetric game, the cone ``K`` is be self-dual. We therefore
 105     name the two interior points ``e1`` and ``e2`` to indicate that
 106     they come from the same cone but are "chosen" by players one and
 107     two respectively.
 108
 109     The ambient space is assumed to be the span of ``K``.
 110     """
 111     def __init__(self, L, K, e1, e2):
 112         """
 113         INPUT:
 114
 115           - ``L`` -- an square matrix represented as a list of lists
 116              of real numbers. ``L`` itself is interpreted as a list of
 117              ROWS, which agrees with (for example) SageMath and NumPy,
 118              but not with CVXOPT (whose matrix constructor accepts a
 119              list of columns).
 120
 121           - ``K`` -- a SymmetricCone instance.
 122
 123           - ``e1`` -- the interior point of ``K`` belonging to player one;
 124             it can be of any enumerable type having the correct length.
 125
 126           - ``e2`` -- the interior point of ``K`` belonging to player two;
 127             it can be of any enumerable type having the correct length.
 128
 129         EXAMPLES:
 130
 131         Lists can (and probably should) be used for every argument:
 132
 133             >>> from cones import NonnegativeOrthant
 134             >>> K = NonnegativeOrthant(2)
 135             >>> L = [[1,0],[0,1]]
 136             >>> e1 = [1,1]
 137             >>> e2 = [1,1]
 138             >>> G = SymmetricLinearGame(L, K, e1, e2)
 139             >>> print(G)
 140             The linear game (L, K, e1, e2) where
 141               L = [ 1  0]
 142                   [ 0  1],
 143               K = Nonnegative orthant in the real 2-space,
 144               e1 = [ 1]
 145                    [ 1],
 146               e2 = [ 1]
 147                    [ 1].
 148
 149         The points ``e1`` and ``e2`` can also be passed as some other
 150         enumerable type (of the correct length) without much harm, since
 151         there is no row/column ambiguity:
 152
 153             >>> import cvxopt
 154             >>> import numpy
 155             >>> from cones import NonnegativeOrthant
 156             >>> K = NonnegativeOrthant(2)
 157             >>> L = [[1,0],[0,1]]
 158             >>> e1 = cvxopt.matrix([1,1])
 159             >>> e2 = numpy.matrix([1,1])
 160             >>> G = SymmetricLinearGame(L, K, e1, e2)
 161             >>> print(G)
 162             The linear game (L, K, e1, e2) where
 163               L = [ 1  0]
 164                   [ 0  1],
 165               K = Nonnegative orthant in the real 2-space,
 166               e1 = [ 1]
 167                    [ 1],
 168               e2 = [ 1]
 169                    [ 1].
 170
 171         However, ``L`` will always be intepreted as a list of rows, even
 172         if it is passed as a ``cvxopt.base.matrix`` which is otherwise
 173         indexed by columns:
 174
 175             >>> import cvxopt
 176             >>> from cones import NonnegativeOrthant
 177             >>> K = NonnegativeOrthant(2)
 178             >>> L = [[1,2],[3,4]]
 179             >>> e1 = [1,1]
 180             >>> e2 = e1
 181             >>> G = SymmetricLinearGame(L, K, e1, e2)
 182             >>> print(G)
 183             The linear game (L, K, e1, e2) where
 184               L = [ 1  2]
 185                   [ 3  4],
 186               K = Nonnegative orthant in the real 2-space,
 187               e1 = [ 1]
 188                    [ 1],
 189               e2 = [ 1]
 190                    [ 1].
 191             >>> L = cvxopt.matrix(L)
 192             >>> print(L)
 193             [ 1  3]
 194             [ 2  4]
 195             <BLANKLINE>
 196             >>> G = SymmetricLinearGame(L, K, e1, e2)
 197             >>> print(G)
 198             The linear game (L, K, e1, e2) where
 199               L = [ 1  2]
 200                   [ 3  4],
 201               K = Nonnegative orthant in the real 2-space,
 202               e1 = [ 1]
 203                    [ 1],
 204               e2 = [ 1]
 205                    [ 1].
 206         """
 207         self._K = K
 208         self._e1 = matrix(e1, (K.dimension(), 1))
 209         self._e2 = matrix(e2, (K.dimension(), 1))
 210
 211         # Our input ``L`` is indexed by rows but CVXOPT matrices are
 212         # indexed by columns, so we need to transpose the input before
 213         # feeding it to CVXOPT.
 214         self._L = matrix(L, (K.dimension(), K.dimension())).trans()
 215
 216         if not K.contains_strict(self._e1):
 217             raise ValueError('the point e1 must lie in the interior of K')
 218
 219         if not K.contains_strict(self._e2):
 220             raise ValueError('the point e2 must lie in the interior of K')
 221
 222     def __str__(self):
 223         """
 224         Return a string representatoin of this game.
 225
 226         EXAMPLES:
 227
 228             >>> from cones import NonnegativeOrthant
 229             >>> K = NonnegativeOrthant(3)
 230             >>> L = [[1,-5,-15],[-1,2,-3],[-12,-15,1]]
 231             >>> e1 = [1,1,1]
 232             >>> e2 = [1,2,3]
 233             >>> SLG = SymmetricLinearGame(L, K, e1, e2)
 234             >>> print(SLG)
 235             The linear game (L, K, e1, e2) where
 236               L = [  1  -5 -15]
 237                   [ -1   2  -3]
 238                   [-12 -15   1],
 239               K = Nonnegative orthant in the real 3-space,
 240               e1 = [ 1]
 241                    [ 1]
 242                    [ 1],
 243               e2 = [ 1]
 244                    [ 2]
 245                    [ 3].
 246
 247         """
 248         tpl = 'The linear game (L, K, e1, e2) where\n' \
 249               '  L = {:s},\n' \
 250               '  K = {!s},\n' \
 251               '  e1 = {:s},\n' \
 252               '  e2 = {:s}.'
 253         indented_L = '\n      '.join(str(self._L).splitlines())
 254         indented_e1 = '\n       '.join(str(self._e1).splitlines())
 255         indented_e2 = '\n       '.join(str(self._e2).splitlines())
 256         return tpl.format(indented_L, str(self._K), indented_e1, indented_e2)
 257
 258
 259     def solution(self):
 260         """
 261         Solve this linear game and return a Solution object.
 262
 263         OUTPUT:
 264
 265         If the cone program associated with this game could be
 266         successfully solved, then a Solution object containing the
 267         game's value and optimal strategies is returned. If the game
 268         could *not* be solved -- which should never happen -- then a
 269         GameUnsolvableException is raised. It can be printed to get the
 270         raw output from CVXOPT.
 271
 272         EXAMPLES:
 273
 274         This example is computed in Gowda and Ravindran in the section
 275         "The value of a Z-transformation":
 276
 277             >>> from cones import NonnegativeOrthant
 278             >>> K = NonnegativeOrthant(3)
 279             >>> L = [[1,-5,-15],[-1,2,-3],[-12,-15,1]]
 280             >>> e1 = [1,1,1]
 281             >>> e2 = [1,1,1]
 282             >>> SLG = SymmetricLinearGame(L, K, e1, e2)
 283             >>> print(SLG.solution())
 284             Game value: -6.1724138
 285             Player 1 optimal:
 286               [ 0.5517241]
 287               [-0.0000000]
 288               [ 0.4482759]
 289             Player 2 optimal:
 290               [0.4482759]
 291               [0.0000000]
 292               [0.5517241]
 293
 294         The value of the following game can be computed using the fact
 295         that the identity is invertible:
 296
 297             >>> from cones import NonnegativeOrthant
 298             >>> K = NonnegativeOrthant(3)
 299             >>> L = [[1,0,0],[0,1,0],[0,0,1]]
 300             >>> e1 = [1,2,3]
 301             >>> e2 = [4,5,6]
 302             >>> SLG = SymmetricLinearGame(L, K, e1, e2)
 303             >>> print(SLG.solution())
 304             Game value: 0.0312500
 305             Player 1 optimal:
 306               [0.0312500]
 307               [0.0625000]
 308               [0.0937500]
 309             Player 2 optimal:
 310               [0.1250000]
 311               [0.1562500]
 312               [0.1875000]
 313
 314         """
 315         # The cone "C" that appears in the statement of the CVXOPT
 316         # conelp program.
 317         C = CartesianProduct(self._K, self._K)
 318
 319         # The column vector "b" that appears on the right-hand side of
 320         # Ax = b in the statement of the CVXOPT conelp program.
 321         b = matrix([1], tc='d')
 322
 323         # A column of zeros that fits K.
 324         zero = matrix(0, (self._K.dimension(), 1), tc='d')
 325
 326         # The column vector "h" that appears on the right-hand side of
 327         # Gx + s = h in the statement of the CVXOPT conelp program.
 328         h = matrix([zero, zero])
 329
 330         # The column vector "c" that appears in the objective function
 331         # value <c,x> in the statement of the CVXOPT conelp program.
 332         c = matrix([-1, zero])
 333
 334         # The matrix "G" that appears on the left-hand side of Gx + s = h
 335         # in the statement of the CVXOPT conelp program.
 336         G = append_row(append_col(zero, -identity(self._K.dimension())),
 337                        append_col(self._e1, -self._L))
 338
 339         # The matrix "A" that appears on the right-hand side of Ax = b
 340         # in the statement of the CVXOPT conelp program.
 341         A = matrix([0, self._e2], (1, self._K.dimension() + 1), 'd')
 342
 343         # Actually solve the thing and obtain a dictionary describing
 344         # what happened.
 345         soln_dict = solvers.conelp(c, G, h, C.cvxopt_dims(), A, b)
 346
 347         # The "status" field contains "optimal" if everything went
 348         # according to plan. Other possible values are "primal
 349         # infeasible", "dual infeasible", "unknown", all of which
 350         # mean we didn't get a solution. That should never happen,
 351         # because by construction our game has a solution, and thus
 352         # the cone program should too.
 353         if soln_dict['status'] != 'optimal':
 354             raise GameUnsolvableException(soln_dict)
 355
 356         p1_value = soln_dict['x'][0]
 357         p1_optimal = soln_dict['x'][1:]
 358         p2_optimal = soln_dict['z'][self._K.dimension():]
 359
 360         return Solution(p1_value, p1_optimal, p2_optimal)
 361
 362     def dual(self):
 363         """
 364         Return the dual game to this game.
 365
 366         EXAMPLES:
 367
 368             >>> from cones import NonnegativeOrthant
 369             >>> K = NonnegativeOrthant(3)
 370             >>> L = [[1,-5,-15],[-1,2,-3],[-12,-15,1]]
 371             >>> e1 = [1,1,1]
 372             >>> e2 = [1,2,3]
 373             >>> SLG = SymmetricLinearGame(L, K, e1, e2)
 374             >>> print(SLG.dual())
 375             The linear game (L, K, e1, e2) where
 376               L = [  1  -1 -12]
 377                   [ -5   2 -15]
 378                   [-15  -3   1],
 379               K = Nonnegative orthant in the real 3-space,
 380               e1 = [ 1]
 381                    [ 2]
 382                    [ 3],
 383               e2 = [ 1]
 384                    [ 1]
 385                    [ 1].
 386
 387         """
 388         return SymmetricLinearGame(self._L, # It will be transposed in __init__().
 389                                    self._K, # Since "K" is symmetric.
 390                                    self._e2,
 391                                    self._e1)
 392
 393
 394 class SymmetricLinearGameTest(TestCase):
 395     """
 396     Tests for the SymmetricLinearGame and Solution classes.
 397     """
 398
 399     def assert_within_tol(self, first, second):
 400         """
 401         Test that ``first`` and ``second`` are equal within our default
 402         tolerance.
 403         """
 404         self.assertTrue(abs(first - second) < options.ABS_TOL)
 405
 406
 407     def assert_solution_exists(self, L, K, e1, e2):
 408         """
 409         Given the parameters needed to construct a SymmetricLinearGame,
 410         ensure that that game has a solution.
 411         """
 412         G = SymmetricLinearGame(L, K, e1, e2)
 413         soln = G.solution()
 414         L_matrix = matrix(L).trans()
 415         expected = inner_product(L_matrix*soln.player1_optimal(),
 416                                  soln.player2_optimal())
 417         self.assert_within_tol(soln.game_value(), expected)
 418
 419     def test_solution_exists_nonnegative_orthant(self):
 420         """
 421         Every linear game has a solution, so we should be able to solve
 422         every symmetric linear game over the NonnegativeOrthant. Pick
 423         some parameters randomly and give it a shot. The resulting
 424         optimal solutions should give us the optimal game value when we
 425         apply the payoff operator to them.
 426         """
 427         ambient_dim = randint(1, 10)
 428         K = NonnegativeOrthant(ambient_dim)
 429         e1 = [uniform(0.1, 10) for idx in range(K.dimension())]
 430         e2 = [uniform(0.1, 10) for idx in range(K.dimension())]
 431         L = [[uniform(-10, 10) for i in range(K.dimension())]
 432              for j in range(K.dimension())]
 433         self.assert_solution_exists(L, K, e1, e2)
 434
 435     def test_solution_exists_ice_cream(self):
 436         """
 437         Like :meth:`test_solution_exists_nonnegative_orthant`, except
 438         over the ice cream cone.
 439         """
 440         # Use a minimum dimension of two to avoid divide-by-zero in
 441         # the fudge factor we make up later.
 442         ambient_dim = randint(2, 10)
 443         K = IceCream(ambient_dim)
 444         e1 = [1]
 445         e2 = [1]
 446         # If we choose the rest of the components of e1,e2 randomly
 447         # between 0 and 1, then the largest the squared norm of the
 448         # non-height part of e1,e2 could be is the 1*(dim(K) - 1). We
 449         # need to make it less than one (the height of the cone) so
 450         # that the whole thing is in the cone. The norm of the
 451         # non-height part is sqrt(dim(K) - 1), and we can divide by
 452         # twice that.
 453         fudge_factor = 1.0 / (2.0*sqrt(K.dimension() - 1.0))
 454         e1 += [fudge_factor*uniform(0, 1) for idx in range(K.dimension() - 1)]
 455         e2 += [fudge_factor*uniform(0, 1) for idx in range(K.dimension() - 1)]
 456         L = [[uniform(-10, 10) for i in range(K.dimension())]
 457              for j in range(K.dimension())]
 458         self.assert_solution_exists(L, K, e1, e2)
 459