]> gitweb.michael.orlitzky.com - dunshire.git/blobdiff - src/dunshire/games.py
Rename "xi" and "gamma" to "p" and "q" to avoid name clashes with CVXOPT.
[dunshire.git] / src / dunshire / games.py
index 9f0b0ced50628d387ef5fbe2b8d166aefc9bb6a8..43fa007c61077c90ef627e9738afbaa09dcde9c1 100644 (file)
@@ -5,17 +5,11 @@ This module contains the main :class:`SymmetricLinearGame` class that
 knows how to solve a linear game.
 """
 
-# These few are used only for tests.
-from math import sqrt
-from random import randint, uniform
-from unittest import TestCase
-
-# These are mostly actually needed.
 from cvxopt import matrix, printing, solvers
-from cones import CartesianProduct, IceCream, NonnegativeOrthant
-from errors import GameUnsolvableException
-from matrices import append_col, append_row, identity, inner_product
-import options
+from .cones import CartesianProduct
+from .errors import GameUnsolvableException
+from .matrices import append_col, append_row, identity
+from . import options
 
 printing.options['dformat'] = options.FLOAT_FORMAT
 solvers.options['show_progress'] = options.VERBOSE
@@ -210,7 +204,7 @@ class SymmetricLinearGame:
     Examples
     --------
 
-        >>> from cones import NonnegativeOrthant
+        >>> from dunshire import *
         >>> K = NonnegativeOrthant(3)
         >>> L = [[1,-5,-15],[-1,2,-3],[-12,-15,1]]
         >>> e1 = [1,1,1]
@@ -231,7 +225,7 @@ class SymmetricLinearGame:
 
     Lists can (and probably should) be used for every argument::
 
-        >>> from cones import NonnegativeOrthant
+        >>> from dunshire import *
         >>> K = NonnegativeOrthant(2)
         >>> L = [[1,0],[0,1]]
         >>> e1 = [1,1]
@@ -253,7 +247,7 @@ class SymmetricLinearGame:
 
         >>> import cvxopt
         >>> import numpy
-        >>> from cones import NonnegativeOrthant
+        >>> from dunshire import *
         >>> K = NonnegativeOrthant(2)
         >>> L = [[1,0],[0,1]]
         >>> e1 = cvxopt.matrix([1,1])
@@ -274,7 +268,7 @@ class SymmetricLinearGame:
     otherwise indexed by columns::
 
         >>> import cvxopt
-        >>> from cones import NonnegativeOrthant
+        >>> from dunshire import *
         >>> K = NonnegativeOrthant(2)
         >>> L = [[1,2],[3,4]]
         >>> e1 = [1,1]
@@ -319,10 +313,10 @@ class SymmetricLinearGame:
         # feeding it to CVXOPT.
         self._L = matrix(L, (K.dimension(), K.dimension())).trans()
 
-        if not K.contains_strict(self._e1):
+        if not self._e1 in K:
             raise ValueError('the point e1 must lie in the interior of K')
 
-        if not K.contains_strict(self._e2):
+        if not self._e2 in K:
             raise ValueError('the point e2 must lie in the interior of K')
 
     def __str__(self):
@@ -363,7 +357,7 @@ class SymmetricLinearGame:
         This example is computed in Gowda and Ravindran in the section
         "The value of a Z-transformation"::
 
-            >>> from cones import NonnegativeOrthant
+            >>> from dunshire import *
             >>> K = NonnegativeOrthant(3)
             >>> L = [[1,-5,-15],[-1,2,-3],[-12,-15,1]]
             >>> e1 = [1,1,1]
@@ -383,7 +377,7 @@ class SymmetricLinearGame:
         The value of the following game can be computed using the fact
         that the identity is invertible::
 
-            >>> from cones import NonnegativeOrthant
+            >>> from dunshire import *
             >>> K = NonnegativeOrthant(3)
             >>> L = [[1,0,0],[0,1,0],[0,0,1]]
             >>> e1 = [1,2,3]
@@ -433,21 +427,41 @@ class SymmetricLinearGame:
         # what happened.
         soln_dict = solvers.conelp(c, G, h, C.cvxopt_dims(), A, b)
 
+        # The optimal strategies are named ``p`` and ``q`` in the
+        # background documentation, and we need to extract them from
+        # the CVXOPT ``x`` and ``z`` variables. The objective values
+        # :math:`nu` and :math:`omega` can also be found in the CVXOPT
+        # ``x`` and ``y`` variables; however, they're stored
+        # conveniently as separate entries in the solution dictionary.
+        p1_value = -soln_dict['primal objective']
+        p2_value = -soln_dict['dual objective']
+        p1_optimal = soln_dict['x'][1:]
+        p2_optimal = soln_dict['z'][self._K.dimension():]
+
         # The "status" field contains "optimal" if everything went
         # according to plan. Other possible values are "primal
-        # infeasible", "dual infeasible", "unknown", all of which
-        # mean we didn't get a solution. That should never happen,
-        # because by construction our game has a solution, and thus
-        # the cone program should too.
-        if soln_dict['status'] != 'optimal':
+        # infeasible", "dual infeasible", "unknown", all of which mean
+        # we didn't get a solution. The "infeasible" ones are the
+        # worst, since they indicate that CVXOPT is convinced the
+        # problem is infeasible (and that cannot happen).
+        if soln_dict['status'] in ['primal infeasible', 'dual infeasible']:
             raise GameUnsolvableException(soln_dict)
-
-        p1_value = soln_dict['x'][0]
-        p1_optimal = soln_dict['x'][1:]
-        p2_optimal = soln_dict['z'][self._K.dimension():]
+        elif soln_dict['status'] == 'unknown':
+            # When we get a status of "unknown", we may still be able
+            # to salvage a solution out of the returned
+            # dictionary. Often this is the result of numerical
+            # difficulty and we can simply check that the primal/dual
+            # objectives match (within a tolerance) and that the
+            # primal/dual optimal solutions are within the cone (to a
+            # tolerance as well).
+            if abs(p1_value - p2_value) > options.ABS_TOL:
+                raise GameUnsolvableException(soln_dict)
+            if (p1_optimal not in self._K) or (p2_optimal not in self._K):
+                raise GameUnsolvableException(soln_dict)
 
         return Solution(p1_value, p1_optimal, p2_optimal)
 
+
     def dual(self):
         r"""
         Return the dual game to this game.
@@ -460,7 +474,7 @@ class SymmetricLinearGame:
         Examples
         --------
 
-            >>> from cones import NonnegativeOrthant
+            >>> from dunshire import *
             >>> K = NonnegativeOrthant(3)
             >>> L = [[1,-5,-15],[-1,2,-3],[-12,-15,1]]
             >>> e1 = [1,1,1]
@@ -487,72 +501,3 @@ class SymmetricLinearGame:
                                    self._K,
                                    self._e2,
                                    self._e1)
-
-
-class SymmetricLinearGameTest(TestCase):
-    """
-    Tests for the SymmetricLinearGame and Solution classes.
-    """
-
-    def assert_within_tol(self, first, second):
-        """
-        Test that ``first`` and ``second`` are equal within our default
-        tolerance.
-        """
-        self.assertTrue(abs(first - second) < options.ABS_TOL)
-
-
-    def assert_solution_exists(self, L, K, e1, e2):
-        """
-        Given the parameters needed to construct a SymmetricLinearGame,
-        ensure that that game has a solution.
-        """
-        G = SymmetricLinearGame(L, K, e1, e2)
-        soln = G.solution()
-        L_matrix = matrix(L).trans()
-        expected = inner_product(L_matrix*soln.player1_optimal(),
-                                 soln.player2_optimal())
-        self.assert_within_tol(soln.game_value(), expected)
-
-    def test_solution_exists_nonnegative_orthant(self):
-        """
-        Every linear game has a solution, so we should be able to solve
-        every symmetric linear game over the NonnegativeOrthant. Pick
-        some parameters randomly and give it a shot. The resulting
-        optimal solutions should give us the optimal game value when we
-        apply the payoff operator to them.
-        """
-        ambient_dim = randint(1, 10)
-        K = NonnegativeOrthant(ambient_dim)
-        e1 = [uniform(0.1, 10) for idx in range(K.dimension())]
-        e2 = [uniform(0.1, 10) for idx in range(K.dimension())]
-        L = [[uniform(-10, 10) for i in range(K.dimension())]
-             for j in range(K.dimension())]
-        self.assert_solution_exists(L, K, e1, e2)
-
-    def test_solution_exists_ice_cream(self):
-        """
-        Like :meth:`test_solution_exists_nonnegative_orthant`, except
-        over the ice cream cone.
-        """
-        # Use a minimum dimension of two to avoid divide-by-zero in
-        # the fudge factor we make up later.
-        ambient_dim = randint(2, 10)
-        K = IceCream(ambient_dim)
-        e1 = [1] # Set the "height" of e1 to one
-        e2 = [1] # And the same for e2
-
-        # If we choose the rest of the components of e1,e2 randomly
-        # between 0 and 1, then the largest the squared norm of the
-        # non-height part of e1,e2 could be is the 1*(dim(K) - 1). We
-        # need to make it less than one (the height of the cone) so
-        # that the whole thing is in the cone. The norm of the
-        # non-height part is sqrt(dim(K) - 1), and we can divide by
-        # twice that.
-        fudge_factor = 1.0 / (2.0*sqrt(K.dimension() - 1.0))
-        e1 += [fudge_factor*uniform(0, 1) for idx in range(K.dimension() - 1)]
-        e2 += [fudge_factor*uniform(0, 1) for idx in range(K.dimension() - 1)]
-        L = [[uniform(-10, 10) for i in range(K.dimension())]
-             for j in range(K.dimension())]
-        self.assert_solution_exists(L, K, e1, e2)
-