This is working beautifully. I went ahead and added some checks to save time on solutions with pure strategies, but even without this its still faster than LRSLib at the nxn matrices tested.
function strat_vec(l::Int64, i::Int64)
vec_to_return = vec(zeros(l))
@inbounds vec_to_return[i] = 1.0
return vec_to_return
end
minmax(R::Array{Float64, 2}, m::Int64) = @inbounds mapreduce(x -> maximum(R[:, x]), min, 1:m)
maxmin(R::Array{Float64, 2}, n::Int64) = @inbounds mapreduce(x -> minimum(R[x, :]), max, 1:n)
findminmax(R::Array{Float64, 2}, n::Int64) = @inbounds strat_vec(n, argmax(map(x -> minimum(R[x, :]), 1:n)))
findmaxmin(R::Array{Float64, 2}, m::Int64) = @inbounds strat_vec(m, argmin(map(x -> maximum(R[:, x]), 1:m)))
function nash(R::Array{Float64, 2})
n, m = size(R)
# Check if we have to do linear programming
n == 1 && return minimum(R), vec([1.0]), strat_vec(m, argmin(vec(R)))
m == 1 && return maximum(R), strat_vec(n, argmin(vec(R))), vec([1.0])
minmax(R, m) == maxmin(R, n) && return minmax(R, m), findminmax(R, n), findmaxmin(R, m)
# Set up model and payoff
model = direct_model(GLPK.Optimizer())
@variable(model, z)
@objective(model, Max, 1.0 * z)
# Solve for row player
@variable(model, x[1:n], lower_bound = 0.0)
@constraint(model, c1, x' * R .>= z)
@constraint(model, sum(x) == 1.0)
optimize!(model)
return value(z), value.(x), shadow_price.(c1)
end
