def max_power(k, n):
if n % k > 0:
return 0
return max_power(k, n / k) + 1
def list_of_p_groups(n):
primes = n.prime_divisors()
degs = []
for i in xrange(len(primes)):
degs.append(max_power(primes[i], n))
return primes, degs
def all_abelian_groups(ord):
primes, degs = list_of_p_groups(ord)
grouplist = []
for i in xrange(len(primes)):
list = []
prime = primes[i]
def all_partitions(parti, n, i):
for j in xrange(1, min(n, i) + 1):
all_partitions(parti + [prime**j], n - j, j)
if n <= 0:
list.append(parti)
all_partitions([], degs[i], degs[i])
grouplist.append(list)
return grouplist
i = 10000
print all_abelian_groups(i)