I still think that the given implementation could be improved.
An alternative would be to work in the Symbolic Ring.
def Bell_partial_polynomial(n,k,v):
@cached_function
def T(n, k):
if k==0: return v[0]^n
return sum(binomial(n-1,j-1)*v[j]*T(n-j,k-1) for j in (0..n-k+1))
return T(n,k)
If one wants to see the representation as a multivariate
polynomial one can define the variables x_0,.._x_{n+1}
R = range(5)
v = var(['x_'+str(i) for i in (0..n+1)])
for n in R:
[Bell_partial_polynomial(n,k,v).expand() for k in (0..n)]
[1]
[x_0, x_1]
[x_0^2, x_0*x_1 + x_2, x_1^2]
[x_0^3, x_0^2*x_1 + 2*x_0*x_2 + x_3, x_0*x_1^2 + 3*x_1*x_2, x_1^3]
[x_0^4, x_0^3*x_1 + 3*x_0^2*x_2 + 3*x_0*x_3 + x_4, x_0^2*x_1^2 +
5*x_0*x_1*x_2 + 3*x_2^2 + 4*x_1*x_3, x_0*x_1^3 + 6*x_1^2*x_2, x_1^4]
Substituting x_0=0 one gets the more familiar form starting with x_1.
for n in R:
[Bell_partial_polynomial(n,k,v).subs(x_0=0) for k in (0..n)]
[1]
[0, x_1]
[0, x_2, x_1^2]
[0, x_3, 3*x_1*x_2, x_1^3]
[0, x_4, 3*x_2^2 + 4*x_1*x_3, 6*x_1^2*x_2, x_1^4]
This is the form best suited for the Stirling/Lah family of
combinatorial numbers. For example to compute the unsigned
Stirling numbers first kind A132393 one can write
R = range(6)
s = [0]+[factorial(i) for i in R]
for n in R:
[Bell_partial_polynomial(n,k,s) for k in (0..n)]
[1]
[0, 1]
[0, 1, 1]
[0, 2, 3, 1]
[0, 6, 11, 6, 1]
[0, 24, 50, 35, 10, 1]
In any case this implementation seems to be more
efficient than the loop in the current implementation
for p in Partitions(n, length=k):
for part, count in p.to_exp_dict().iteritems():
We give a last example which shows that there are other
important examples which do /not/ reduce to integers.
R = range(9)
t = [2^(1-j) if is_odd(j) else 0 for j in (0..n+2)]
for n in R:
[Bell_partial_polynomial(n,k,t) for k in (0..n)]
[1]
[0, 1]
[0, 0, 1]
[0, 1/4, 0, 1]
[0, 0, 1, 0, 1]
[0, 1/16, 0, 5/2, 0, 1]
[0, 0, 1, 0, 5, 0, 1]
[0, 1/64, 0, 91/16, 0, 35/4, 0, 1]
[0, 0, 1, 0, 21, 0, 14, 0, 1]
These are the coefficients of the central factorials.
http://mathworld.wolfram.com/CentralFactorial.html
///////////////////////////////////////////////////////////
///////////////////////////////////////////////////////////
Bell polynomials reconsidered.
=============================
Consider the transformation
T(n, k, x) = Sum_{j=0..n-k+1} K(n,j)*x[j]*T(n-j,k-1,f)
if k != 0 else T(n, k, x) = x[0]^n.
As kernels K(n,j) we consider here the two cases
K(n,j) = binomial(n-1,j-1) and
K(n,j) = number_of_partitions(n,j),
where number_of_partitions(n,j) is the number of j-partition of n.
If the kernel K is fixed than T can be seen as a
sequence-to-triangle transformation, mapping the sequence
x to the triangle T(x), defined by
T(x)_{n,k} = T(n, k, x).
Let us look at the first few rows of the generated triangles.
Case K(n,j) = binomial(n-1,j-1):
x[0]^0
x[0]^1, x[1]
x[0]^2, x[1]*x[0]+x[2], x[1]^2
x[0]^3, x[1]*x[0]^2+2*x[2]*x[0]+x[3], x[1]^2*x[0]+3*x[2]*x[1], x[1]^3
x[0]^4, x[1]*x[0]^3+3*x[2]*x[0]^2+3*x[3]*x[0]+x[4], x[1]^2*x[0]^2+5*x[1]*x[2]*x[0]
+4*x[3]*x[1]+3*x[2]^2, x[1]^3*x[0]+6*x[2]*x[1]^2, x[1]^4
Case K(n,j) = number_of_partitions(n,j):
x[0]^0
x[0]^1, x[1]
x[0]^2, x[1]*x[0]+2*x[2], x[1]^2
x[0]^3, x[1]*x[0]^2+2*x[2]*x[0]+3*x[3], x[1]^2*x[0]+4*x[2]*x[1], x[1]^3
x[0]^4, x[1]*x[0]^3+3*x[2]*x[0]^2+4*x[3]*x[0]+5*x[4], x[1]^2*x[0]^2
+5*x[1]*x[2]*x[0]+7*x[3]*x[1]+6*x[2]^2, x[1]^3*x[0]+7*x[2]*x[1]^2, x[1]^4
Note that these triangles have the offset (0,0), i. e. the enumeration is
T(0,0)
T(1,0) T(1,1)
T(2,0) T(2,1) T(2,2)
T(3,0) T(3,1) T(3,2) T(3,3)
T(4,0) T(4,1) T(4,2) T(4,3) T(4,4)
This is important to note because in the older literature
the first column T(n,0) often is ignored and the triangle
starts at T(1,1). This happens in particular in the important
special cases which arise when x[0] = 0. Then the above
two triangles reduce to, in the
case K(n,j) = binomial(n-1,j-1):
1
0, x[1]
0, x[2], x[1]^2
0, x[3], 3*x[2]*x[1], x[1]^3
0, x[4], 4*x[3]*x[1] + 3*x[2]^2 , 6*x[2]*x[1]^2 , x[1]^4
and in the case K(n,j) = number_of_partitions(n,j) to:
1
0, x[1]
0, 2*x[2], x[1]^2
0, 3*x[3], 4*x[2]*x[1], x[1]^3
0, 5*x[4], 7*x[3]*x[1] + 6*x[2]^2 , 7*x[2]*x[1]^2 , x[1]^4
Let's get numerical now and look for instance at the indicator
sequence of the positive numbers, 1 - 0^n (in other words,
let's set x[0] = 0 and x[k] = 1 for all k>0).
Then the first triangle becomes
1;
0, 1;
0, 1, 1;
0, 1, 3, 1;
0, 1, 7, 6, 1;
0, 1, 15, 25, 10, 1;
0, 1, 31, 90, 65, 15, 1;
which is A048993, the triangle of the Stirling set numbers.
The second triangle becomes
1;
0, 1;
0, 2, 1;
0, 3, 4, 1;
0, 5, 13, 7, 1;
0, 7, 30, 30, 10, 1;
0, 11, 76, 119, 65, 14, 1;
a triangle which has not yet been explored.
But let's now switch our point of view. Until now we
understood x as a sequence x[0], x[1], x[2],....
Now we fix the kernel binomial(n-1,j-1), write x_n instead of x[n]
and understand x_n as a variable in some ring. Then an expression
like x_1*x_0^3+3*x_2*x_0^2+3*x_3*x_0+x_4 (see the fifth line
in the first triangle above) becomes an element of a polynomial
ring in n-k+1 variables.
These polynomials are well known and famous: they are called
the Bell polynomials in the special case where x_0 is
substituted by 0 and the first column is disregarded.
x_1
x_2, x_1^2
x_3, 3*x_2*x_1, x_1^3
x_4, 4*x_3*x_1 + 3*x_2^2, 6*x_2*x_1^2, x_1^4
This is the definition as given for example by L. Comtet
in Advanced Combinatorics, page 135. However this definition
is not satisfactory. For example substituting 1 for the x_n
in this beheaded triangle we get
1;
1, 1;
1, 3, 1;
1, 7, 6, 1;
1, 15, 25, 10, 1;
1, 31, 90, 65, 15, 1;
which is A008277, the triangle of the trimmed Stirling numbers
of the second kind. (In passing note that the first comment
in A008277 is not correct. It says: "Also known as Stirling set
numbers and written {n,k}." However the notation {n,k} and
the term 'Stirling set numbers' were introduced by D. E. Knuth
explicitly to denote the full triangle A048993, not the
trimmed Stirling numbers.)
Of course there are many other instances where the reduced
Bell polynomials as defined by Comtet do not comply anymore
to the way important numbers nowadays are derived from the
Bell polynomials. The reader might check this for example
for the case x[n] = 1!,2!,3!,.. (Lah numbers) or
x[n] = 0!,1!,2!,3!,.. (signless Stirling numbers of the
first kind).
To sum up: We strongly suggest to use the polynomials with
the full set of variables x_0, x_1, ... and the enumeration
of the generated number triangles starting at (0,0).
This gives, among other important numbers, the Stirling set
numbers, the Stirling cycle numbers and the signless Lah
numbers as used today by the majority of writers.
---
With Maple one can nicely generate the above triangles
with four lines of code:
T := proc(n,k,x) option remember; `if`(k=0,x[0]^n,
add(binomial(n-1,j-1)*x[j]*T(n-j,k-1,x),j=0..n-k+1)) end:
for n from 0 to 6 do seq(simplify(T(n,k,x)),k=0..n) od;
for n from 0 to 6 do seq(simplify(subs(x[0]=0,T(n,k,x))), k=0..n) od;