mgritter · January 28, 2024 08:16 · Jan 28, 2024
diff --git a/compute.py b/compute.py
@@ -0,0 +1,96 @@
+outcomes = ["AS", "AF", "NS", "NF", "DS", "DF"]
+
+# Game played with N-sided dice
+# Player may reroll R times
+# Success threshold is max of D dice.
+def prob(N, R, D):
+    P = {}
+    EV = {}
+    advantage = range( 2 * N // 3 + 1, N+1 )
+    disadvantage = range( 1, N // 3 + 1 )
+    neutral = range( N // 3 + 1, 2 * N // 3 + 1 )
+
+    for a in range( 1, N+1 ):
+        pS = (a * 1.0 / N) ** D
+        pF = 1.0 - pS
+        for b in range( 1, N+1 ):
+            for o in outcomes:
+                P[(o,a,b,0)] = 0.0
+
+            if b in advantage:
+                P[("AS",a,b,0)] = pS
+                P[("AF",a,b,0)] = pF                
+            elif b in disadvantage:
+                P[("DS",a,b,0)] = pS
+                P[("DF",a,b,0)] = pF                
+            else:    
+                P[("NS",a,b,0)] = pS
+                P[("NF",a,b,0)] = pF
+
+    for r in range( 1, R+1 ):
+        for a in range(1,N+1):
+            for b in range(1,N+1):
+                EV[(a,b,r-1)] = \
+                    4 * P[("AS",a,b,r-1)] + \
+                    3 * P[("NS",a,b,r-1)] + \
+                    2 * P[("DS",a,b,r-1)] + \
+                    2 * P[("AF",a,b,r-1)] + \
+                    1 * P[("NF",a,b,r-1)]
+
+        for b in range( 1, N+1 ):
+            #print( f"EV[*,{b},{r-1}] = 1/{N} * (", end='' )
+            #for a in range(1, N+1):
+            #    print( f"{EV[(a,b,r-1)]} + ", end='' )
+            #print( ")")
+            EV[("*",b,r-1)] = sum(EV[(a,b,r-1)] for a in range(1,N+1)) / N
+
+        for a in range(1,N+1):
+            for b in range(1,N+1):
+                # if EV( a, b, 0 ) < EV( *, a, r-1 ) then reroll
+                if EV[(a,b,0)] < EV[("*",a,r-1)]:
+                    #print( f'{a},{b} {r} {EV[(a,b,0)]} vs {EV[("*",a,r-1)]} reroll' )
+                    # Reroll moves a to b and replaces a with unform outcome [1,N]
+                    for o in outcomes:
+                        P[(o,a,b,r)] = sum( P[(o,k,a,r-1)] for k in range(1,N+1) ) / N
+                else:
+                    #print( f'{a},{b} {r} {EV[(a,b,0)]} vs {EV[("*",a,r-1)]} stop' )
+                    for o in outcomes:
+                        P[(o,a,b,r)] = P[(o,a,b,0)]
+
+    # Initial state probabilities
+    init = {}
+    for o in outcomes:
+        init[o] = sum( P[(o,a,b,R)]
+                       for a in range(1,N+1)
+                       for b in range(1,N+1) ) / (N ** 2)
+
+    return init, P, EV
+
+def show_decision(N, r, EV):
+    print( f"[{r} rerolls left]" )
+    print( "a\\b " + " ".join( str(x) for x in range(1,N+1) ) )
+    for a in range( 1, N+1):
+        print( f"{a:3} ", end='') 
+        for b in range(1, N+1):
+            if EV[(a,b,0)] < EV[("*",a,r-1)]:
+                print( "R ", end='')
+            else:
+                print( "x ", end='')
+        print()
+
+def main():
+    N = 6
+    R = 10
+    for D in range(1,11):
+        print( f"### Difficulty {D}")
+        init, p, ev = prob(N, R, D)
+        for o in outcomes:
+            print( f"{o} {init[o]}" )
+
+        print()
+        for r in range(1, R+1):
+            show_decision(N, r, ev)
+        print()
+
+if __name__ == "__main__":
+    main()