# given binsizes=b_1, b_2, ..., b_n, compute sum(b_i)! / (b_1! b_2! ... b_n!)
# or 0 if sum(b_i) = 0. This is the number of ways to arrange sum(b_i) objects
# in a row,
# where group i has b_i indistinguishable elements, and objects from different
# groups are distinguishable
# s is an optional argument, which if supplied should be equal to sum(binsizes)
def multichoose(binsizes, s=-1):
if s == -1:
s = sum(binsizes)
if s == 0:
return 0
i = 0
result = 1
while i < len(binsizes):
result *= choose(s, binsizes[i])
s -= binsizes[i]
i += 1
return result
# Given a combination, written as a list of numbers like C = [0,1,1,0,0], assigns
# a number N(C) to this combination. N is invertible, mapping any n, 0 <= n <
# max_C {N(C)} to some C, where the maximum is taken over all C matching
# the pattern of the given C (in this case, the pattern is 3 0's and 2 1's). The
# numbers appearing in C must be in range(0,multichoose(binsizes)), with no gaps
def combcount(C):
binsizes = [0]*len(set(C))
for i in xrange(len(C)):
binsizes[C[i]] += 1
s = sum(binsizes)
result = 0
for i in xrange(len(C)):
# add to result the number of valid combinations V whose prefixes (up to index i-1)
# are as given, and where V[i] < C[i]
for j in xrange(C[i]):
if binsizes[j] == 0:
continue
binsizes[j] -= 1
s -= 1
result += multichoose(binsizes, s)
binsizes[j] += 1
s += 1
binsizes[C[i]] -= 1
s -= 1
return result
# given a number N and a list of bin sizes (as in the argument of multichoose),
# where 0 <= N < multichoose(binsizes)
# produce the combination C matching binsizes such that combcount(C)=N.
def combcreate(N,binsizes):
s = sum(binsizes)
l = len(binsizes)
C = [0]*s
for i in xrange(s):
mchoose = 0
j = -1
# subtract from N the number of valid combinations V whose prefixes (up to
# index i-1) are as already found, and where V[i] < j, where j is the largest
# integer such that N remains nonnegative during this process
while N >= 0:
newj = j+1
while newj != l and binsizes[newj] == 0:
newj += 1
mchoose = 0
if newj == l:
break
j = newj
binsizes[j] -= 1
s -= 1
mchoose = multichoose(binsizes,s)
N -= mchoose
binsizes[j] += 1
s += 1
C[i] = j
N += mchoose
binsizes[j] -= 1
s -= 1
return C
def choose(n, k):
# iterative solution, dividing first when possible
nmink = n - k
if nmink < k:
nmink, k = k, nmink
result = 1
while n > 0:
if n == nmink:
break
result *= n
n -= 1
if k > 1 and result % k == 0:
result /= k
k -= 1
for i in xrange(2,k+1):
result /= i
return result
