# Independent check (my own code, not the bundle's): exact Wald coverage on n=5..100, p=0.001..0.999
from math import comb
from fractions import Fraction
import math
A=1959963984540054; B=10**15; G=1000
b93=b95=0; zf=1.959963984540054; b93f=b95f=0; ties=0
for n in range(5,101):
    ch=[comb(n,k) for k in range(n+1)]
    den=G**n
    for i in range(1,G):
        num=0; numf=0.0
        for k in range(1,n):
            lhs=(i*n-k*G)**2*n*B*B; rhs=A*A*k*(n-k)*G*G
            if lhs==rhs: ties+=1
            if lhs<=rhs: num+=ch[k]*i**k*(G-i)**(n-k)
        if 100*num<93*den: b93+=1
        if 100*num<95*den: b95+=1
print('pairs',96*999,'below93',b93,'below95',b95,'ties',ties)
print('shares',round(b93/95904,6),round(b95/95904,6))
