91 lines
1.8 KiB
Python
91 lines
1.8 KiB
Python
# Apply Burnside's Lemma to a group
|
|
# G = S_r x S_c. G is a group action on
|
|
# the set of R x C matrices with elements
|
|
# 1 to n. (S_r, S_c can be thought to
|
|
# permute the rows and columns).
|
|
# We need only count the total number
|
|
# matrices unchanged by a group operation.
|
|
# This is the cycle index of S_r x S_c.
|
|
|
|
from itertools import permutations
|
|
|
|
def fac(n):
|
|
a = 1
|
|
for i in range(2, n+1): a *= i
|
|
return a
|
|
|
|
|
|
def gcd(a, b):
|
|
while b > 0:
|
|
a, b = b, a % b
|
|
return a
|
|
|
|
|
|
def mult(l):
|
|
v = 1
|
|
for x in l: v *= x
|
|
return v
|
|
|
|
|
|
def gen_struct(full, l, n, i=0, v=0):
|
|
if v >= n:
|
|
if v == n: full.append([x for x in l])
|
|
return
|
|
if i == n: return
|
|
l[i] = (n-v) / (i+1)
|
|
for c in range(l[i]):
|
|
gen_struct(full, l, n, i+1, v+l[i]*(i+1))
|
|
l[i] -= 1
|
|
gen_struct(full, l, n, i+1, v+l[i]*(i+1))
|
|
|
|
|
|
def cycle_monomials(n):
|
|
# Generates the cycle structures of distinct
|
|
# lengths of S_n and their number of occurances.
|
|
a = range(1, n+1)
|
|
perms = permutations(range(1, n+1))
|
|
valid = {}
|
|
structs = []
|
|
gen_struct(structs, [0]*n, n)
|
|
|
|
for s in structs:
|
|
valid[tuple(s)] = distinct_disjoint(n, s)
|
|
return valid
|
|
|
|
|
|
def distinct_disjoint(n, s):
|
|
# Counts distinct number of cycle structures
|
|
# given a structure.
|
|
a = mult([fac(x)*(i+1)**x for i,x in enumerate(s)])
|
|
return fac(n) / a
|
|
|
|
|
|
def precompute_gcds(n):
|
|
gcds = {}
|
|
for i in range(1, n+1):
|
|
for j in range(1, n+1):
|
|
gcds[(i,j)] = gcd(i, j)
|
|
return gcds
|
|
|
|
|
|
def count_pow(w, h, gcds):
|
|
p = 0
|
|
for i,x in enumerate(w):
|
|
for j,y in enumerate(h):
|
|
p += x*y*gcds[(i+1, j+1)]
|
|
return p
|
|
|
|
|
|
def solution(w, h, s):
|
|
gcds = precompute_gcds(w if w > h else h)
|
|
cycle_w = cycle_monomials(w)
|
|
cycle_h = cycle_monomials(h)
|
|
|
|
a = 0
|
|
for m_w in cycle_w:
|
|
for m_h in cycle_h:
|
|
power = count_pow(m_w, m_h, gcds)
|
|
a += cycle_w[m_w]*cycle_h[m_h]*s**(power)
|
|
return str(a/(fac(w)*fac(h)))
|
|
|