Merge d633619 into 113fea9

src/repeated_game.jl

@@ -83,7 +83,7 @@ best_dev_payoff_2(rpd::RepGame2, a1::Int) =
 """
     unitcircle(npts)
 
-Places `npts` equally spaced points along the 2 dimensional circle and returns 
+Places `npts` equally spaced points along the 2 dimensional circle and returns
 the points with x coordinates in first column and y coordinates in second
 column.
 """
@@ -105,7 +105,7 @@ end
 """
     initialize_sg_hpl(nH, o, r)
 
-Initializes subgradients, extreme points and hyperplane levels for the 
+Initializes subgradients, extreme points and hyperplane levels for the
 approximation of the convex value set of a 2 player repeated game.
 
 # Arguments
@@ -118,7 +118,7 @@ approximation of the convex value set of a 2 player repeated game.
 
 - `C::Array{Float64}(nH, 1)` : The array containing the hyperplane levels.
 - `H::Array{Float64}(nH, 2)` : The array containing the subgradients.
-- `Z::Array{Float64}(nH, 2)` : The array containing the extreme points of the 
+- `Z::Array{Float64}(nH, 2)` : The array containing the extreme points of the
   value set.
 """
 function initialize_sg_hpl(nH::Int, o::Vector{Float64}, r::Float64)
@@ -144,7 +144,7 @@ end
 """
     initialize_sg_hpl(rpd, nH)
 
-Initializes subgradients, extreme points and hyperplane levels for the 
+Initializes subgradients, extreme points and hyperplane levels for the
 approximation of the convex value set of a 2 player repeated game by choosing
 an appropriate origin and radius.
 
@@ -157,7 +157,7 @@ an appropriate origin and radius.
 
 - `C::Array{Float64}(nH, 1)` : The array containing the hyperplane levels.
 - `H::Array{Float64}(nH, 2)` : The array containing the subgradients.
-- `Z::Array{Float64}(nH, 2)` : The array containing the extreme points of the 
+- `Z::Array{Float64}(nH, 2)` : The array containing the extreme points of the
   value set.
 """
 function initialize_sg_hpl(rpd::RepeatedGame, nH::Int)
@@ -179,7 +179,7 @@ end
 """
     initialize_LP_matrices(rpd, H)
 
-Initialize matrices for the linear programming problems. 
+Initialize matrices for the linear programming problems.
 
 # Arguments
 
@@ -188,11 +188,11 @@ Initialize matrices for the linear programming problems.
 
 # Returns
 
-- `c::Array{Float64}(nH, 1)` : Vector used to determine which subgradient is 
+- `c::Array{Float64}(nH, 1)` : Vector used to determine which subgradient is
   being used.
-- `A::Array{Float64}(nH, 2)` : Matrix with nH set constraints and to be filled 
+- `A::Array{Float64}(nH, 2)` : Matrix with nH set constraints and to be filled
   with 2 additional incentive compatibility constraints.
-- `b::Array{Float64}(nH, 1)` : Vector to be filled with the values for the 
+- `b::Array{Float64}(nH, 1)` : Vector to be filled with the values for the
   constraints.
 """
 function initialize_LP_matrices(rpd::RepGame2, H)
@@ -221,19 +221,23 @@ end
 
 Given a constraint w ∈ W, this finds the worst possible payoff for agent i.
 
-# Arugments 
+# Arugments
 
 - `rpd::RepGame2` : Two player repeated game.
 - `H::Array{Float64, 2}` : The subgradients used to approximate the value set.
 - `C::Array{Float64, 1}` : The array containing the hyperplane levels.
 - `i::Int` : The player of interest.
+- `lp_solver` : Allows users to choose a particular solver for linear
+   programming problems. Options include ClpSolver(), CbcSolver(),
+   GLPKSolverLP() and GurobiSolver(). By default, it choooses ClpSolver().
+
 
 # Returns
 
 - `out::Float64` : Worst possible payoff for player i.
 """
 function worst_value_i(rpd::RepGame2, H::Array{Float64, 2},
-                       C::Array{Float64, 1}, i::Int)
+                       C::Array{Float64, 1}, i::Int, lp_solver=ClpSolver())
     # Objective depends on which player we are minimizing
     c = zeros(2)
     c[i] = 1.0
@@ -242,7 +246,7 @@ function worst_value_i(rpd::RepGame2, H::Array{Float64, 2},
     lb = [-Inf, -Inf]
     ub = [Inf, Inf]
 
-    lpout = linprog(c, H, '<', C, lb, ub, ClpSolver())
+    lpout = linprog(c, H, '<', C, lb, ub, lp_solver)
     if lpout.status == :Optimal
         out = lpout.sol[i]
     else
@@ -253,22 +257,24 @@ function worst_value_i(rpd::RepGame2, H::Array{Float64, 2},
 end
 
 "See worst_value_i for documentation"
-worst_value_1(rpd::RepGame2, H::Array{Float64, 2}, C::Array{Float64, 1}) =
-    worst_value_i(rpd, H, C, 1)
+worst_value_1(rpd::RepGame2, H::Array{Float64, 2}, C::Array{Float64, 1},
+              lp_solver=ClpSolver()) = worst_value_i(rpd, H, C, 1)
 "See worst_value_i for documentation"
-worst_value_2(rpd::RepGame2, H::Array{Float64, 2}, C::Array{Float64, 1}) =
-    worst_value_i(rpd, H, C, 2)
+worst_value_2(rpd::RepGame2, H::Array{Float64, 2}, C::Array{Float64, 1},
+              lp_solver=ClpSolver()) = worst_value_i(rpd, H, C, 2)
 "See worst_value_i for documentation"
-worst_values(rpd::RepGame2, H::Array{Float64, 2}, C::Array{Float64, 1}) =
-    (worst_value_1(rpd, H, C), worst_value_2(rpd, H, C))
+worst_values(rpd::RepGame2, H::Array{Float64, 2}, C::Array{Float64, 1},
+             lp_solver=ClpSolver()) = (worst_value_1(rpd, H, C),
+                                       worst_value_2(rpd, H, C))
 
 #
 # Outer Hyper Plane Approximation
 #
 """
     outerapproximation(rpd; nH=32, tol=1e-8, maxiter=500, check_pure_nash=true,
-                       verbose=false, nskipprint=50, 
-                       plib=getlibraryfor(2, Float64))
+                       verbose=false, nskipprint=50,
+                       plib=getlibraryfor(2, Float64),
+                       lp_solver=ClpSolver())
 
 Approximates the set of equilibrium value set for a repeated game with the
 outer hyperplane approximation described by Judd, Yeltekin, Conklin 2002.
@@ -280,24 +286,28 @@ outer hyperplane approximation described by Judd, Yeltekin, Conklin 2002.
   - `tol` : Tolerance in differences of set.
   - `maxiter` : Maximum number of iterations.
   - `verbose` : Whether to display updates about iterations and distance.
-  - `nskipprint` : Number of iterations between printing information 
+  - `nskipprint` : Number of iterations between printing information
     (assuming verbose=true).
   - `check_pure_nash`: Whether to perform a check about whether a pure Nash
     equilibrium exists.
-  - `plib`: Allows users to choose a particular package for the geometry 
+  - `plib`: Allows users to choose a particular package for the geometry
           computations.
           (See [Polyhedra.jl](https://github.com/JuliaPolyhedra/Polyhedra.jl)
-          docs for more info). By default, it chooses to use 
+          docs for more info). By default, it chooses to use
           [CDDLib.jl](https://github.com/JuliaPolyhedra/CDDLib.jl)
+  - `lp_solver` : Allows users to choose a particular solver for linear
+     programming problems. Options include ClpSolver(), CbcSolver(),
+     GLPKSolverLP() and GurobiSolver(). By default, it choooses ClpSolver().
 
-# Returns 
+# Returns
 
   - `vertices::Array{Float64}` : Vertices of the outer approximation of the
     value set.
 """
 function outerapproximation(rpd::RepGame2; nH=32, tol=1e-8, maxiter=500,
                             check_pure_nash=true, verbose=false, nskipprint=50,
-                            plib=getlibraryfor(2, Float64))
+                            plib=getlibraryfor(2, Float64),
+                            lp_solver=ClpSolver())
     # Long unpacking of stuff
     sg, delta = unpack(rpd)
     p1, p2 = sg.players
@@ -332,8 +342,8 @@ function outerapproximation(rpd::RepGame2; nH=32, tol=1e-8, maxiter=500,
     iter, dist = 0, 10.0
     while (iter < maxiter) & (dist > tol)
         # Compute the current worst values for each agent
-        _w1 = worst_value_1(rpd, H, C)
-        _w2 = worst_value_2(rpd, H, C)
+        _w1 = worst_value_1(rpd, H, C, lp_solver)
+        _w2 = worst_value_2(rpd, H, C, lp_solver)
 
         # Update all set constraints -- Copies elements 1:nH of C into b
         copy!(b, 1, C, 1, nH)
@@ -365,7 +375,7 @@ function outerapproximation(rpd::RepGame2; nH=32, tol=1e-8, maxiter=500,
                           (1-delta)*best_dev_payoff_2(rpd, a1) - delta*_w2
 
                 # Solve corresponding linear program
-                lpout = linprog(c, A, '<', b, lb, ub, ClpSolver())
+                lpout = linprog(c, A, '<', b, lb, ub, lp_solver)
                 if lpout.status == :Optimal
                     # Pull out optimal value and compute
                     w_sol = lpout.sol
@@ -393,7 +403,7 @@ function outerapproximation(rpd::RepGame2; nH=32, tol=1e-8, maxiter=500,
             end
 
             # Update the points
-            Z[:, ih] = (1-delta)*flow_u(rpd, a1star, a2star) + delta*[Wstar[1], 
+            Z[:, ih] = (1-delta)*flow_u(rpd, a1star, a2star) + delta*[Wstar[1],
                         Wstar[2]]
         end