#OK to post homework
#Victoria Chayes, 1/30/2022, Assignment 3

Help:=proc(): print(`AllGames(a,b),FullStat(a,b),SampleStat(a,b,K) `) end:

#############################################

#########

#3 Using (a) G:=RandDisGame(4,4) ; (b) G:=RandDisGame(5,7); Find (BY HAND!)
#   (i) The outcome of the sequential (aka dynamical) version of the game where Player Row goes first
#   (ii) The outcome of the sequential (aka dynamical) version of the game where Player Column goes first
#   (iii) The set (possibly empty) of (pure) Nash Equilibria
#Check these answers using procedures DynRC, DynCR, and your own PureNashEqui(G) that hopefully you wrote for hw2
 
# G:=RandDisGame(4,4): matrix(G)
#  [[1, 10], [3, 14], [16, 1],  [12, 4]]
#  [[6, 16], [9, 8],  [11, 7],  [8, 15]]
#  [[13, 9], [2, 3],  [5, 6],   [10, 2]]
#  [[7, 13], [14, 5], [15, 12], [4, 11]]]

# (i) if Player Row goes first, their considerations are: they choose 1, PC chooses 2 for a payoff of 3; they choose 2, PC chooses 1 for a payoff of 6; they choose 3, PC chooses 1 for a payoff of 13; or they choose 4 and PC chooses 1 for a payoff of 7. Therefore, PR should choose 3 if they go first for the maximum payoff. We check with:
# DynRC(G)
#     [3, 1], [13, 9]

# (ii) if Player Column goes first, their considerations are: they choose 1, PR chooses 3 for payoff 9; they choose 2, PR chooses 4 for payoff 5; they choose 3, PR chooses 1 for payoff 1; or they choose 4 and PR chooses 1 for payoff 4. Therefore, PC should choose 1 for maximum payoff. We check with
# DynCR(G)
#     [3, 1], [13, 9]

# (iii) [3, 1] is clearly the only Nash Equilibrium from the argument above.
# PureNashEqui(G)
#     {[3, 1]}

 
# G:=RandDisGame(5,7): matrix(G)
# [[33, 25], [23, 30], [26, 26], [11, 16], [20, 27], [35, 23], [32, 12]] 
# [[19, 10], [4, 13],  [17, 18], [25, 29], [5, 31],  [8, 7],   [31, 5]] 
# [[10, 4],  [2, 33],  [18, 15], [1, 6],   [3, 21],  [15, 24], [9, 17]] 
# [[22, 20], [16, 9],  [12, 14], [7, 3],   [6, 32],  [29, 11], [28, 34]] 
# [[27, 8],  [21, 35], [14, 1],  [30, 2],  [13, 22], [34, 28], [24, 19]]]

# (i) if Player Row goes first, their considerations are: they choose 1, PC chooses 2 for a payoff of 23; they choose 2, PC chooses 5 for payoff of 5; they choose 3, PC chooses 2 for payoff of 2; they choose 4, PC chooses 7 for payoff of 28; or they choose 5, PC chooses 2 for a payoff of 21. Thus they are best choosing 4.
# DynRC(G)
#     [4, 7], [28, 34]

# (ii) if Player Column goes first, their considerations are: they choose 1, PR chooses 1 for payoff 25; they choose 2, PR chooses 1 for payoff 30; they choose 3, PR chooses 1 for payoff 26; they choose 4, PR chooses 5 for payoff 2; they choose 6, PR chooses 1 for payoff 23; or they choose 7, PR chooses 1 for payoff 12. Thus they are best choosing 2.
# DynCR(G)
#     [1, 2], [23, 30]

# (iii) it is fairly simple to spot that R1 is always Player 1's best choice, except in Column 4, where it is R5. As Column 4 is not Player 2's best choice in Row 5, it can be discarded. In R1, Column 2 is Player 2's best choice. Therefore, [1,2] is the only Nash Equilibria 
# PureNashEqui(G)
#     {[1, 2]}


#########

#4 Using Maple's command permute(n), write a procedure AllGames(a,b) that inputs positive integers a and b and outputs the set of all (a*b)!^2 possible games where Player Row has a strategies, Player Column has b strategies, each player's payoffs are distinct and drawn from {1, ..., a*b}

AllGames:=proc(a,b) local S, L, i:
if a=1 and b=1 then
return({[[1,1]]}):
fi:
L:=permute(permute(a*b,2)):
S:={}:
for i from 1 to nops(L) do
 S:=S union {[seq([seq(L[i][(k-1)*b+j],j=1..b)],k=1..a)]}
od:
S:
end:

#########

#5 Write a procedure FullStat(a,b) that inputs positive integers a and b (not too big!) and outputs the following list [Number of Games with No Nash Equilibria, Number of Games with One Nash Equilibria, Number of Games with Two Nash Equilibria, Number of Games where the Row Player was better off playing second, Number of Games where the Column Player was better off playing second, Number of Games where DynRC(G)=DynCR(G)] What are (i) FullStat(2,2), (ii) FullStat(3,2). 

FullStat:=proc(a,b) local AG,L,i,x1,C1,R1:
AG:=AllGames(a,b):
L:=[0,0,0,0,0,0]:
for i from 1 to nops(AG) do
 x1:=nops(PureNashEqui(AG[i])):
  if x1=0 then
   L[1]:=L[1]+1:
  elif x1=1 then
   L[2]:=L[2]+1:
  elif x1=2 then
   L[3]:=L[3]+1:
  fi:
 R1:=[DynRC(AG[i])][2]:
 C1:=[DynCR(AG[i])][2]:
 if C1[1]>=R1[1] then
   L[4]:=L[4]+1:
 fi:
 if R1[2]>=C1[2] then
   L[5]:=L[5]+1:
 fi:
 if R1=C1 then    # as we know that payoffs are unique, we can do this
   L[6]:=L[6]+1:
 fi:
od: 
L:
end:

#NOTE: this would all be a lot faster if we didn't have to reconstruct the games in AllGames, but instead kept them as a list and used [a,b] to indicate the dimensions. In that case, we use:

FullStatFast:=proc(a,b) local AG,L,i,x1,C1,R1:
AG:=AllGamesFast(a,b):
L:=[0,0,0,0,0,0]:
for i from 1 to nops(AG) do
 x1:=nops(PureNashEquiFast(AG[i],a,b)):
  if x1=0 then
   L[1]:=L[1]+1:
  elif x1=1 then
   L[2]:=L[2]+1:
  elif x1=2 then
   L[3]:=L[3]+1:
  fi:
 R1:=[DynRCFast(AG[i],a,b)][2]:
 C1:=[DynCRFast(AG[i],a,b)][2]:
 if C1[1]>=R1[1] then
   L[4]:=L[4]+1:
 fi:
 if R1[2]>=C1[2] then
   L[5]:=L[5]+1:
 fi:
 if R1=C1 then 
   L[6]:=L[6]+1:
 fi:
od: 
L:
end:

AllGamesFast:=proc(a,b):
permute(permute(a*b,2)):
end: 

BR1dFast:=proc(G,a,b,a2) local a1:
max[index]([seq(G[(a1-1)*b+a2][1],a1=1..a)]):
end:

BR2dFast:=proc(G,a,b,a1) local a2:
max[index]([seq(G[(a1-1)*b+a2][2],a2=1..b)]):
end:

BR1dvFast:=proc(G,a,b) local a2: [seq(BR1dFast(G,a,b,a2),a2=1..b)]:end:
BR2dv:=proc(G) local a1: [seq(BR2dFast(G,a,b,a1),a1=1..a)]:end:

DynRCFast:=proc(G,a,b) local L,i,iBest, jBest:
L:=BR2dvFast(G,a,b):
iBest:=max[index]([seq(G[(i-1)*b+L[i]][1],i=1..a)]):
jBest:=L[iBest]:
[iBest,jBest], G[(iBest-1)*b+jBest]:
end:

DynCRFast:=proc(G,a,b) local L,j,iBest, jBest:
L:=BR1dvFast(G,a,b):
jBest:=max[index]([seq(G[(L[j]-1)*a+j][2],j=1..b)]):
iBest:=L[jBest]:
[iBest,jBest], G[(iBest-1)*b+jBest]:
end:

BR1Fast:=proc(G,a,b,a2) local a1:
MyMaxIndex([seq(G[(a1-1)*b+a2][1],a1=1..a)]):
end:
BR2Fast:=proc(G,a,b,a1) local a2:
MyMaxIndex([seq(G[(a1-1)*b+a2][2],a2=1..b)]):
end:
BR1vFast:=proc(G,a,b) local a2: [seq(BR1(G,a,b,a2),a2=1..b)]: end:
BR2vFast:=proc(G,a,b) local a1: [seq(BR2(G,a,b,a1),a1=1..a)]: end:

PureNashEquiFast:=proc(G,a,b) local N1,N2,S,i,j,k:
N1:=BR1vFast(G,a,b):
N2:=BR2vFast(G,a,b):
S:={}:
for i from 1 to nops(N1) do
for j from 1 to nops(N1[i]) do
k:=N1[i][j]:
if member(i,N2[k]) then
S:=S union {[k,i]}:
fi:
od:
od:
S:
end:


#########

#6 Write a procedure SampleStat(a,b,K) that inputs positive integers a and b, and a fairly large positive integer K, and generates K random games where Row has a strategies and Col has b strategies and outputs the following list of numbers (in floating-point, i.e. decimals) evalf( [(Number of Games with No Nash Equilibria)/K, (Number of Games with One Nash Equilibria)/K, (Number of Games with Two Nash Equilibria)/K, (Number of Games where the Row Player was better off playing second)/K, (Number of Games where the Column Player was better off playing second)/K, (Number of Games where DynRC(G)=DynCR(G))/K]): Run the following five times (i) SampleStat(10,10,10000) and see if the answers are similar.

SampleStat:=proc(a,b,K) local G,L,i,x1,C1,R1:
L:=[0,0,0,0,0,0]:
for i from 1 to K do
G:=RandDisGame(a,b):
 x1:=nops(PureNashEqui(G)):
  if x1=0 then
   L[1]:=L[1]+1:
  elif x1=1 then
   L[2]:=L[2]+1:
  elif x1=2 then
   L[3]:=L[3]+1:
  fi:
 R1:=[DynRC(G)][2]:
 C1:=[DynCR(G)][2]:
 if C1[1]>=R1[1] then
   L[4]:=L[4]+1:
 fi:
 if R1[2]>=C1[2] then
   L[5]:=L[5]+1:
 fi:
 if R1=C1 then 
   L[6]:=L[6]+1:
 fi:
od: 
evalf(L/~K):
end:

# We see:
# SampleStat(10,10,10000)
#   [0.3282000000, 0.4075000000, 0.2081000000, 0.6333000000, 0.6316000000, 0.3372000000]
# SampleStat(10,10,10000)
#   [0.3330000000, 0.4079000000, 0.2012000000, 0.6346000000, 0.6336000000, 0.3367000000]
# SampleStat(10,10,10000)
#   [0.3258000000, 0.4073000000, 0.2055000000, 0.6282000000, 0.6258000000, 0.3361000000]
# SampleStat(10,10,10000)
#   [0.3230000000, 0.4099000000, 0.2012000000, 0.6239000000, 0.6253000000, 0.3363000000]
# SampleStat(10,10,10000)
#   [0.3278000000, 0.4059000000, 0.2059000000, 0.6325000000, 0.6332000000, 0.3377000000]


#########


#############################################

# Programs from HW 2

IsNash:=proc(G,a1,a2):
if member(a1, BR1(G,a2)) and member(a2, BR2(G,a1)) then
return(true):
fi:
false:
end:

PureNashEqui:=proc(G) local N1,N2,S,i,j,k:
N1:=BR1v(G):
N2:=BR2v(G):
S:={}:
for i from 1 to nops(N1) do
  for j from 1 to nops(N1[i]) do
    k:=N1[i][j]:
    if member(i,N2[k]) then
     S:=S union {[k,i]}:
    fi:
 od:
od:
S:
end:

#############################################

#  Programs From C3.txt  



#C3.txt: Maple Code for Lecture 2 of Math640 (Spring 2022)

with(combinat):

Help3:=proc(): print(`RandDisGame(a,b), BR1d(G,a2), BR2d(G,a1),  BR1dv(G), BR2dv(G), DynRC(G), DynCR(G)  `):end:


#RandDisGame(a,b): A random game where Player Row has a strategies R1, R2, ..., Ra; and Player Col. has b strategies: C1, C2, ..., Cb and all pay-offs are DISCTINCT
#and consist of the set {1,2,...,ab}. Try:
#RandDisGame(4,6);
RandDisGame:=proc(a,b) local pi1,pi2,i1,j1:
pi1:=randperm(a*b):
pi2:=randperm(a*b):
[seq([seq([pi1[b*i1+j1],pi2[b*i1+j1]],j1=1..b)],i1=0..a-1)]:
end:



#BR1d(G,a2): Inputs a bimatrix of a game G and a member a2 of A2 outputs the UNIQUE member of A1 that consists of the best response  
BR1d:=proc(G,a2) local a1:
max[index]([seq(G[a1][a2][1],a1=1..nops(G))]):
end:



#BR2d(G,a1): Inputs a bimatrix of a game G and a member a1 of A1 outputs the UNIQUE member of A2 that consists of the best response  
BR2d:=proc(G,a1) local a2:
max[index]([seq(G[a1][a2][2],a2=1..nops(G[a1]))]):
end:

#BR1dv(G): Inputs a bimatrix of a game G, outputs The list of size A2 with all the  best responses
#
BR1dv:=proc(G) local a2: [seq(BR1d(G,a2),a2=1..nops(G[1]))]:end:

#BR2dv(G): Inputs a bimatrix of a game G, outputs The list of size A1 with all the  best responses
BR2dv:=proc(G) local a1: [seq(BR2d(G,a1),a1=1..nops(G))]:end:


#DynRC(G): The outcome of the Dynamical version of the game G if Row goes first and Column goes next
DynRC:=proc(G) local L,i,iBest, jBest:
#L is the list of size A1 where L[i] is the best response of Player Column to Strategy i and Row gets G[i][L[i]][1]
L:=BR2dv(G):
iBest:=max[index]([seq(G[i][L[i]][1],i=1..nops(G))]):
jBest:=L[iBest]:
[iBest,jBest], G[iBest][jBest]:
end:


#DynCR(G): The outcome of the Dynamical version of the game G if Column goes first and Row goes next
DynCR:=proc(G) local L,j,iBest, jBest:
#L is the list of size A2 where L[j] is the best response of Player Row to Strategy j and Col gets G[L[j]][j]][2]
L:=BR1dv(G):
jBest:=max[index]([seq(G[L[j]][j][2],j=1..nops(L))]):
iBest:=L[jBest]:
[iBest,jBest], G[iBest][jBest]:
end:


#  Programs From C2.txt

MyMaxIndex:=proc(L) local i,m,S:
m:=max(L):
S:={}:
for i from 1 to nops(L) do
 if L[i]=m then
    S:=S union {i}:
 fi:
od:
S:
end:

BR1:=proc(G,a2) local a1:
MyMaxIndex([seq(G[a1][a2][1],a1=1..nops(G))]):
end:

BR2:=proc(G,a1) local a2:
MyMaxIndex([seq(G[a1][a2][2],a2=1..nops(G[a1]))]):
end:

BR1v:=proc(G) local a2: [seq(BR1(G,a2),a2=1..nops(G[1]))]:end:

BR2v:=proc(G) local a1: [seq(BR2(G,a1),a1=1..nops(G))]:end: