#OK to post homework #Blair Seidler, 2/13/22, Assignment 6 Help:=proc(): print(`Only one procedure: TragedyOfTheCommons(v,G,c,n)`): end: #1. Using procedure RG(a,b,K) of Lecture 5 (also contained in C6.txt ), Create five random games # with a=2, b=2, K=1000, and for each of these games, pick three random choices of p1 and p2 and # run SimulateMG(G,p1,p2,K) with K=10000, each five times. Make sure that you get similar-looking # answers, and compare them to the output of evalf(PayOffG(G,p1,p2)) (* ***** Game #1 ***** > G := RG(2, 2, 1000); G := [[[860, 758], [750, 889]], [[300, 991], [5, 993]]] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 19 p1 := -- 50 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 69 p2 := --- 100 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [445.5940000, 917.1584000] [442.2690000, 918.9353000] [443.2070000, 918.1132000] [440.6175000, 918.7794000] [443.5795000, 918.8297000] > evalf(PayOffG(G, p1, p2)); [443.1430000, 918.2762000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 7 p1 := --- 100 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 27 p2 := -- 50 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [211.2060000, 979.7989000] [208.3465000, 980.1328000] [209.0110000, 979.7008000] [210.7765000, 979.5621000] [209.5665000, 979.5698000] > evalf(PayOffG(G, p1, p2)); [209.4570000, 979.7638000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 93 p1 := --- 100 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 3 p2 := - 5 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [770.2060000, 823.8369000] [773.3840000, 821.7731000] [773.9770000, 822.7253000] [774.3885000, 821.1641000] [770.9080000, 822.9965000] > evalf(PayOffG(G, p1, p2)); [771.6200000, 823.0980000] ***** Game #2 ***** > G := RG(2, 2, 1000); G := [[[967, 465], [194, 211]], [[555, 688], [938, 296]]] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 3 p1 := -- 25 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 11 p2 := -- 50 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [795.1889000, 369.0831000] [792.9347000, 366.9929000] [794.9052000, 368.6315000] [798.6635000, 366.8864000] [792.8251000, 367.0049000] > evalf(PayOffG(G, p1, p2)); [794.9784000, 368.3968000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 19 p1 := -- 50 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 61 p2 := --- 100 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [689.2530000, 471.4596000] [686.8019000, 470.3740000] [691.5465000, 469.0708000] [688.1559000, 471.0907000] [688.2527000, 472.7467000] > evalf(PayOffG(G, p1, p2)); [771.6200000, 823.0980000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 31 p1 := --- 100 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 1 p2 := - 2 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [697.9403000, 443.9856000] [695.0133000, 441.4749000] [697.3548000, 443.6382000] [699.4435000, 444.4960000] [692.2643000, 444.7463000] > evalf(PayOffG(G, p1, p2)); [695.0400000, 444.2600000] ***** Game #3 ***** > G := RG(2, 2, 1000); G := [[[261, 551], [977, 132]], [[356, 636], [117, 780]]] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 13 p1 := -- 20 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 9 p2 := -- 25 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [532.7690000, 442.3163000] [533.1746000, 442.1480000] [534.0412000, 442.1528000] [537.4828000, 439.8419000] [542.6639000, 435.3370000] > evalf(PayOffG(G, p1, p2)); [538.5700000, 438.7020000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 73 p1 := --- 100 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 17 p2 := --- 100 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [663.8895000, 354.6470000] [665.4949000, 352.6853000] [666.0278000, 352.6778000] [666.8171000, 351.9252000] [661.6614000, 355.8570000] > evalf(PayOffG(G, p1, p2)); [666.9145000, 352.3483000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 16 p1 := -- 25 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 19 p2 := -- 25 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [388.5562000, 527.2916000] [385.9390000, 528.1479000] [380.9591000, 532.2700000] [385.3477000, 529.3007000] [385.1527000, 529.5004000] > evalf(PayOffG(G, p1, p2)); [384.5280000, 529.6832000] ***** Game #4 ***** > G := RG(2, 2, 1000); G := [[[616, 329], [225, 227]], [[135, 759], [750, 763]]] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 17 p1 := -- 20 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 67 p2 := --- 100 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [465.7064000, 365.7412000] [464.9416000, 363.0826000] [464.5919000, 367.7148000] [464.9772000, 361.1790000] [465.4746000, 366.3310000] > evalf(PayOffG(G, p1, p2)); [464.6170000, 365.0870000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 1 p1 := -- 10 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 31 p2 := -- 50 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [375.1972000, 712.7764000] [382.0539000, 714.6196000] [378.9847000, 716.3052000] [377.6215000, 714.8776000] [382.9124000, 714.0464000] > evalf(PayOffG(G, p1, p2)); [378.5720000, 713.4920000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 14 p1 := -- 25 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 93 p2 := --- 100 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [409.9126000, 515.1058000] [407.3084000, 514.7706000] [409.1617000, 513.3498000] [403.2051000, 518.9568000] [406.9373000, 516.2636000] > evalf(PayOffG(G, p1, p2)); [407.9748000, 514.3248000] ***** Game #5 ***** > G := RG(2, 2, 1000); G := [[[83, 26], [394, 954]], [[298, 480], [165, 105]]] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 63 p1 := --- 100 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 19 p2 := --- 100 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [281.5436000, 554.4475000] [281.3569000, 554.0073000] [280.8850000, 553.3260000] [280.4871000, 552.3771000] [279.5245000, 547.4326000] > evalf(PayOffG(G, p1, p2)); [281.3932000, 555.1509000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 1 p1 := -- 20 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 61 p2 := --- 100 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [244.3855000, 337.3723000] [243.4998000, 334.9660000] [243.4778000, 335.0054000] [243.4763000, 334.6920000] [243.8357000, 336.0661000] > evalf(PayOffG(G, p1, p2)); [244.0380000, 336.4585000] > p1 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 1 p1 := - 5 > p2 := RandomTools[Generate](rational(denominator = 100, range = 0 .. 1)); 13 p2 := -- 20 > for i to 5 do > SimulateMG(G, p1, p2, 10000); > end do; [238.8412000, 347.6795000] [239.1114000, 348.8330000] [239.8867000, 350.1612000] [239.4086000, 349.1328000] [240.8485000, 351.9137000] > evalf(PayOffG(G, p1, p2)); [239.5300000, 349.1600000] *) #3. [A little bit challenging] Confirm the "Tragedy of the commons by confirming that the solution # of Eq. (1.2.6) for n >1 is larger than the solution of Eq. (1.2.7) (which is the former equation # with n=1), and hence that too many goats are grazed in the Nash equilibrium when more than one # farmer shares the common pasture. #TragedyOfTheCommons(v,G,c,n): inputs # (i) An expression,v, for a function that is positive in the variable G, in 0 ≤ G ≤ 1, # but whose first and second derivatives are negative there (like in Fig. 1.2.4 in the book with G_max=1); # (ii) a number c, the cost of buying and caring for the goat, (between 0 and 1); #(iii) a positive integer n, the number of farmers; #and outputs the Nash Equilibrium G*. TragedyOfTheCommons:=proc(v,G,c,n) local f: fsolve(v+G*diff(v,G)/n-c=0,G): end: (* Output > seq(TragedyOfTheCommons(sqrt(1 - G), G, 0.5, n), n = 1 .. 10); 0.4108027069, 0.5240408206, 0.5790991062, 0.6120715420, 0.6341606640, 0.6500444109, 0.6620396765, 0.6714307918, 0.6789893673, 0.6852079369 > seq(TragedyOfTheCommons((1 - G)^(1/4), G, 0.5, n), n = 1 .. 10); 0.5980909164, 0.7159713145, 0.7695811633, 0.8009276122, 0.8217454810, 0.8366883478, 0.8479928247, 0.8568756419, 0.8640589838, 0.8700002405 > seq(TragedyOfTheCommons((1 - G)^(3/4), G, 0.5, n), n = 1 .. 10); 0.3111336476, 0.4083516515, 0.4565628023, 0.4855049838, 0.5048474518, 0.5187020389, 0.5291205867, 0.5372433987, 0.5437555690, 0.5490938253 As expected, as the number of farmers get larger, the total number of goats in the Nash equilibrium gets farther and farther away from the social optimum. Too many goats! *) #### From C6.txt #### #C6.txt: Feb. 7, 2022 Help6:=proc(): print(` SimulateG(G,p1,p2), SimulateMG(G,p1,p2,K) , PayOffG(G,p1,p2) `): end: #SimulateG(G,p1,p2) #Simulates ONE 2-person game with 2 strategies given by a bimatrix G each # game where the stochastic Strategy of Player is to play Odd with Pron. p1 #and Player 2 with prob. p2 SimulateG:=proc(G,p1,p2) local c1,c2: #c1:=What player 1 played #1 is ODD #2 is EVEN c1:=LoadedCoin(p1): c2:=LoadedCoin(p2): G[c1][c2]: end: #SimulateMG(G,p1,p2,K): Plays the game G with mixed strategy (p1,p2) K times and outputs #the pair of average payoffs SimulateMG:=proc(G,p1,p2,K) local i: evalf(add(SimulateG(G,p1,p2),i=1..K)/K): end: #PayOffG(G,p1,p2): Inputs a 2-player game with two strategies for each player, given in terms of its bi-matrix, and rational numbers p1 and p2 between #0 and 1 (inclusive) (of left as symbols), outputs the pair [PayOffOfPlayer1,PayOffOfPlayer2] for the expected pay-off if #Player 1 adopts strategy R1 with prob. p1 (and hence Strategy R2 with prob. 1-p1) and #Player 2 adopts strategy C1 with prob. p2 (and hence Strategy C2 with prob. 1-p2) PayOffG:=proc(G,p1,p2): expand(p1*p2*G[1][1]+ p1*(1-p2)*G[1][2]+(1-p1)*p2*G[2][1]+ (1-p1)*(1-p2)*G[2][2]): end: #OLD STUFF #C5.txt Maple Code for Lecture 5 Help5:=proc(): print(` IsNE(G,a1,a2), NE(G), BestTot(G) , BetterForBoth(G,a1,a2), RG(a,b,K) `):end: #BetterForBoth(G,a1,a2): Given a game G and a strategy choice (a1,a2) finds all the strategy choices that are better for BOTH players BetterForBoth:=proc(G,a1,a2) local b1,b2,S: S:={}: for b1 from 1 to nops(G) do for b2 from 1 to nops(G[1]) do if G[b1][b2][1]>G[a1][a2][1] and G[b1][b2][2]>G[a1][a2][2] then S:=S union {[b1,b2]}: fi: od: od: S: end: #BestTot(G): Given a 2-player game G, outputs the places where the total pay-off is the best, together with the pay-off BestTot:=proc(G) local a1,a2,S,rec: S:={[1,1]}: rec:=G[1][1][1]+G[1][1][2]: for a1 from 1 to nops(G) do for a2 from 1 to nops(G[1]) do if G[a1][a2][1]+G[a1][a2][2]>rec then S:={[a1,a2]}: rec:= G[a1][a2][1]+G[a1][a2][2]: elif G[a1][a2][1]+G[a1][a2][2]=rec then S:=S union {[a1,a2]}: fi: od: od: S,rec: end: #IsNE(G,a1,a2): Is [a1,a2] a Nash Equilibrium of the game G? IsNE:=proc(G,a1,a2) member(a1,BR1(G,a2)) and member(a2,BR2(G,a1)): end: #NE(G): Given a 2-player game G given by a bi-matrix G finds the set of Nash Equilibria NE:=proc(G) local a1,a2,S: S:={}: for a1 from 1 to nops(G) do for a2 from 1 to nops(G[1]) do if IsNE(G,a1,a2) then S:=S union {[a1,a2]}: fi: od: od: S: end: #RG(a,b,K): A random 2-player (static) game where the Row player has a strategies (Row 1, .., Row a), the Column player has b strategies (Col. 1, ..., Crol. b) #and the pay-offs are from 0 to K. For example, to see a random game where Player Row has four strategy choices, and Player Column has five strategy choices #and the pay-offs are integers from 0 to 20 type: RG:=proc(a,b,K) local ra,i,j: ra:=rand(0..K): [seq([seq([ra(),ra()],j=1..b)],i=1..a)]: end: #Start FROM C4.txt #C4.txt: Maple Code for Lecture 4 of Math640 (Spring 2022) Help4:=proc(): print(` FP(F), BR12(G), BR21(G) `): end: #FP(f): Given a discrete function f:= {1,...,n} ->{1,...,n} described by a list of length n, find the set of fixed points FP:=proc(L) local S,i: S:={}: for i from 1 to nops(L) do if L[i]=i then S:=S union {i}: fi: od: S: end: #BR12(G): Given a game G outputs the discrete function f: A1->A1 such that f(a1)= Best Response of Player 1 to the Best Response of Player 2 to Player 1's strategy a1. #Try: G:=RandDisGame(G,10,20): BR12(G); BR12:=proc(G) local L1,L2,a1: L1:=BR1dv(G): L2:=BR2dv(G): [seq(L1[L2[a1]],a1=1..nops(L2))]: end: #BR21(G): Given a game G outputs the discrete function g: A2->A2 such that g(a2)= Best Response of Player 2 to the Best Response of Player 1 to Player 2's strategy a2. #Try: G:=RandDisGame(G,10,20): BR21(G); BR21:=proc(G) local L1,L2,a2: L1:=BR1dv(G): L2:=BR2dv(G): [seq(L2[L1[a2]],a2=1..nops(L1))]: end: #FROM C3.txt Help3:=proc(): print(`RandDisGame(a,b), BR1d(G,a2), BR2d(G,a1), BR1dv(G), BR2dv(G), DynRC(G), DynCR(G) `):end: with(combinat): #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: #End C3.txt #FROM C2.txt Help2:=proc(): print(` GameDB(), Rand2PlayerGame(a,b,K), IsStictDom(v1,v2), FindR(G), FindC(G), ShrinkGame(G), ReducedGame(G) , MyMaxIndex(L), LoadedCoin(p), BR1(G,a2), BR2(G,a1), BR1v(G), BR2v(G) `):end: ###PREPARED BEFORE CLASS #GameDB(): A list of length 5 consisiting of a a "data base" of famous gamess (and less famous ones): #Prisoner's dillema, Boattle of the Sexes, Matching pennies, Figure 1.1.1. in Gibbons, Figure 1.1.4 in Gibbons #For example, to see the Prisoner's dillema, type #GameDB()[1]) GameDB:=proc(): [ [ [ [[-1,-1],[-9,0]], [[0,-9],[-6,-6]]], [`Mum`, `Fink`], [`Mum`, `Fink`]], [ [ [[2,1],[0,0]], [[0,0],[1,2]]], [`Box`, `Opera`],[`Box`, `Opera`]], [ [ [[-1,1],[1,-1]], [[1,-1],[-1,1]]], [`Odd`, `Even`], [`Odd`, `Even`]], [ [ [ [1,0],[1,2],[0,1]], [ [0,3],[0,1],[2,0]]],[`Up`, `Down`], [`Left`, `Middle`, `Right`]], [[ [ [0,4],[4,0],[5,3]],[ [4,0],[0,4],[5,3]], [ [3,5],[3,5],[6,6]]], [`T`, `M`, `B`],[`L`,`C`,`R`]], [ [ [[0,0],[-1,1],[1,-1]], [[1,-1],[0,0],[-1,1]], [[-1,1],[1,-1],[0,0]] ], [`Scissors`, `Rock`, `Paper`], [`Scissors`, `Rock`, `Paper`] ] ]: end: PrintGame:=proc(G) local i: matrix( [ [" ", op(G[3])], seq([G[2][i],op(G[1][i])],i=1..nops(G[2])) ]): end: #Rand2PlayerGame(a,b,K): A random 2-player (static) game where the Row player has a strategies (Row 1, .., Row a), the Column player has b strategies (Col. 1, ..., Crol. b) #and the pay-offs are from 0 to K. For example, to see a random game where Player Row has four strategy choices, and Player Column has five strategy choices #and the pay-offs are integers from 0 to 20 type: #matrix(Rand2PlayerGame(4,5,20)); Rand2PlayerGame:=proc(a,b,K) local ra,i,j: ra:=rand(0..K): [ [seq([seq([ra(),ra()],j=1..b)],i=1..a)],[seq(i,i=1..a)],[seq(j,j=1..b)]]: end: #IsStictDom(v1,v2): Given two lists of numbers is v1[i]<=v2[i] for all i. For example #IsStrictDom([1,3,2],[2,4,3]); should return true but IsDom([1,3,2],[2,4,1]); should return false IsStrictDom:=proc(v1,v2) local i: for i from 1 to nops(v1) do if v1[i]>=v2[i] then RETURN(false): fi: od: true: end: #FindR(G): Finds a strictly dominated row stategy FindR:=proc(G) local G1,RowS,i1,i2,j,v1,v2: G1:=G[1]: RowS:=G[2]: for i1 from 1 to nops(RowS) do for i2 from i1+1 to nops(RowS) do v1:=[seq(G1[i1][j][1],j=1..nops(G1[i1]))]: v2:=[seq(G1[i2][j][1],j=1..nops(G1[i2]))]: if IsStrictDom(v1,v2) then RETURN(i1): elif IsStrictDom(v2,v1) then RETURN(i2): fi: od: od: FAIL: end: #FindC(G): Finds a strictly dominated row stategy FindC:=proc(G) local G1,ColS,j1,j2,i,v1,v2: G1:=G[1]: ColS:=G[3]: for j1 from 1 to nops(ColS) do for j2 from j1+1 to nops(ColS) do v1:=[seq(G1[i][j1][2],i=1..nops(G1))]: v2:=[seq(G1[i][j2][2],i=1..nops(G1))]: if IsStrictDom(v1,v2) then RETURN(j1): elif IsStrictDom(v2,v1) then RETURN(j2): fi: od: od: FAIL: end: #ShrinkGame(G): Given a game G, tries to shrink it, if it can't it returs it ShrinkGame:=proc(G) local i,j,G1,RowS,ColS,i1: G1:=G[1]: RowS:=G[2]; ColS:=G[3]: i:=FindR(G): if i<>FAIL then RETURN([ [op(1..i-1,G1),op(i+1..nops(G1),G1)], [op(1..i-1,RowS),op(i+1..nops(RowS),RowS)],ColS ]): fi: j:=FindC(G): if j<>FAIL then RETURN([ [seq( [op(1.. j-1,G1[i1]),op(j+1..nops(G1[i1]),G1[i1])],i1=1..nops(G1))],RowS, [op(1..j-1,ColS),op(j+1..nops(ColS),ColS)]]): fi: FAIL: end: #ReducedGame(G): The reduced game of G after all the possible Elimination of dominated strategy ReducedGame:=proc(G) local G1,G2: G1:=G: G2:=ShrinkGame(G1): while G2<>FAIL do G1:=G2: G2:=ShrinkGame(G1): od: G1: end: ####END PREPARED BEFORE CLASS #DONE DURING CLASS WITH STUDENTS' HELP #MyMaxIndex(L): [Suggested by Victoria Chayes]. #Inputs a list L of numbers, output the SET of places (indices) where it is maximum. For example #MyMaxIndex([5,6,7,7,1,7]); should give {3,4,6} 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: #Clarification (Jan. 26, 2022) thanks to Blair Seidler: Note that G is the actual BI-MATRIX, so to experiment with it use, for example #G:=Rand2PlayerGame(10,20,100)[1] #BR1(G,a2): Inputs a bi-matrix of a 2-player game G and a member a2 of A2 outputs the SUBSET of A1 that consists of the best response BR1:=proc(G,a2) local a1: MyMaxIndex([seq(G[a1][a2][1],a1=1..nops(G))]): end: #BR2(G,a1): Inputs a bimatrix of a game G and a member a1 of A1 outputs the SUBSET of A2 that consists of the best response BR2:=proc(G,a1) local a2: MyMaxIndex([seq(G[a1][a2][2],a2=1..nops(G[a1]))]): end: #BR1v(G): Inputs a bimatrix of a game G, outputs The list of size A2 with all the set of best responses # BR1v:=proc(G) local a2: [seq(BR1(G,a2),a2=1..nops(G[1]))]:end: #BR2v(G): Inputs a bimatrix of a game G, outputs The list of size A1 with all the set of best responses BR2v:=proc(G) local a1: [seq(BR2(G,a1),a1=1..nops(G))]:end: #end DURING CLASS WITH STUDENTS' HELP #Dne after class #LoadedCoin(p): outputs 1 or 2 with probability p and 1-p respecively. Try: #[seq(LoadedCon(1/3)),i=1..300)]; LoadedCoin:=proc(p) local m,n,ra: m:=numer(p): n:=denom(p): ra:=rand(1..n)(): if ra<=m then 1: else 2: fi: end: #End FROM C2.txt