The Partition Transform
Peter Luschny, March 2016
We introduce a sequence-to-triangle transformation which we will call the partition transformation (or short P-transformation). The name derives from the main loop which runs over all partitions of n.
We say, $p$ is an integer partition of $n$ or $p \in P_{n}$ $\Longleftrightarrow$ $p = (p_{1}, \ldots, p_{n}), \ p_{i}$ integer, $\sum_{1\le i \le n} p = n$ and $p_{1} \ge p_{2} \ge \ldots \ge p_{n} \ge 0.$
def PartitionTransform(n, f, norm = None):
if n == 0: return [1]
R = [0]*(n+1)
for q in Partitions(n):
p = q + [0]
R[q[0]] += (-1)^q[0]*mul(binomial(p[j], p[j+1])*f(j+1)^p[j]
for j in range(len(q)))
if norm == None: return R
return [norm(n, k)*R[k] for k in (0..n)]
lah2 = lambda n: 1 if n == 1 else ((n-1)^2+1)/(n*(4*n-2))
norm = lambda n,k: (-1)^k*factorial(2*n)/factorial(2*k)
for n in range(8): print PartitionTransform(n, lah2, norm)
eul = lambda n: 1/((2*n-1)*(2*n))
nrm = lambda n,k: (-1)^k*factorial(2*n)
for n in range(8): print PartitionTransform(n, eul, nrm)
ward = lambda n: 1/(n+1)
norm = lambda n,k: (-1)^k*falling_factorial(n+k, n)
for n in range(8): print PartitionTransform(n, ward, norm)
def PMatrix(dim, f, norm = None, inverse = False, reduced = False):
i = 1; F = [1]
if reduced:
while i <= dim: F.append(f(i)); i += 1
else:
while i <= dim: F.append(F[i-1]*f(i)); i += 1
C = [[0 for k in range(m+1)] for m in range(dim)]
C[0][0] = 1
if inverse:
for m in (1..dim-1):
C[m][m] = -C[m-1][m-1]/F[1]
for k in range(m-1, 0, -1):
C[m][k] = -(C[m-1][k-1]+sum(F[i]*C[m][k+i-1]
for i in (2..m-k+1)))/F[1]
else:
for m in (1..dim-1):
C[m][m] = -C[m-1][m-1]*F[1]
for k in range(m-1, 0, -1):
C[m][k] = -sum(F[i]*C[m-i][k-1] for i in (1..m-k+1))
if norm == None: return C
for m in (1..dim-1):
for k in (1..m): C[m][k] *= norm(m,k)
return C
1 | 0 | 0 | 0 | 0 | 0 |
0 | -x1 | 0 | 0 | 0 | 0 |
0 | -x1x2 | x12 | 0 | 0 | 0 |
0 | -x1x2x3 | 2x12x2 | -x13 | 0 | 0 |
0 | -x1x2x3x4 | 2x12x2x3+x12x22 | -3x13x2 | x14 | 0 |
0 | -x1x2x3x4x5 | 2x12x2x3x4+2x12x22x3 | -3x13x2x3−3x13x22 | 4x14x2 | -x15 |
f = lambda n: var('x'+str(n))
PMatrix(5, f)
f = lambda n: var('x'+str(n))
P = PMatrix(5, f, inverse = True)
for p in P: print p
# The above output is a bit messy.
# Additionally most of the time we can require x1 = 1.
# So let's clean this up:
[[expand(p).substitute(x1=1) for p in L] for L in P]
PMatrix(8, lambda n: n)
def PEvaluateAt(n, f, val = None, norm = None):
P = PMatrix(n, f, norm)
if val == None: val = var('x')
p = [sum([val^k*P[j][k] for k in (0..j)]) for j in range(n)]
return p
# for linguistic reason only, equivalent to PEvaluateAt(n, f, norm)
PAsPolynomials = lambda n, f, norm=None: PEvaluateAt(n, f, norm = norm)
PAsPolynomials(5, lambda n: n)
PEvaluateAt(9, lambda n: n, val = 1)
def PAnalyser(Len, f, norm = None):
print "As polynomials:"
P = PAsPolynomials(Len-2, f, norm)
for p in P: print p
print "As triangle:"
P = PMatrix(Len-2, f, norm)
for p in P: print p
print "x=1: ", PEvaluateAt(Len,f, 1,norm)
print "x=-1: ", PEvaluateAt(Len,f, -1,norm)
print "x=1/2: ", PEvaluateAt(Len,f, 1/2,lambda n,k: 2^n*norm(n,k))
print "x=-1/2:", PEvaluateAt(Len,f,-1/2,lambda n,k: 2^n*norm(n,k))
\begin{equation} E_{2n} = (2n)! \sum_{p \in P_n} (-1)^{p_1} \binom{p_1}{p_2}\binom{p_2}{p_3} \ldots \binom{p_{n-1}}{p_{n}} \left(\frac{1}{1 \cdot 2}\right)^{p_1} \left(\frac{1}{3 \cdot 4}\right)^{p_2} \ldots \left(\frac{1}{(2n-1)2n}\right)^{p_n} \end{equation}
# Euler(2n)
eul = lambda n: 1/((2*n-1)*(2*n))
nrm = lambda n,k: factorial(2*n)
PAnalyser(8, eul, nrm)
\begin{equation} B_{2n} = \frac{(2n)!}{(2-2^{2n})} \sum_{p \in P_n} (-1)^{p_1} \binom{p_1}{p_2}\binom{p_2}{p_3} \ldots \binom{p_{n-1}}{p_{n}} \left(\frac{1}{2 \cdot 3}\right)^{p_1} \left(\frac{1}{4 \cdot 5}\right)^{p_2} \ldots \left(\frac{1}{2n (2n+1)}\right)^{p_n} \end{equation}
# Bernoulli(2n)
bern = lambda n: 1/((2*n)*(2*n+1))
norm = lambda n,k: factorial(2*n)/(2-2^(2*n))
PAnalyser(8, bern, norm)
def Pcoeffsum(len, f, norm = None):
R, C = [1], [1]+[0]*(len-1)
for n in (1..len-1):
for k in range(n, 0, -1):
C[k] = C[k-1] * f(k)
C[0] = -sum(C[k] for k in (1..n))
if norm == None: R.append(C[0])
else: R.append(C[0]*norm(n))
return R
print Pcoeffsum(9, lambda n: n)
print PEvaluateAt(9, lambda n: n, 1)
def A001896_list(len):
R, C = [1], [1]+[0]*(len-1)
for n in (1..len-1):
for k in range(n, 0, -1):
C[k] = C[k-1] / (8*k*(2*k+1))
C[0] = -sum(C[k] for k in (1..n))
R.append((C[0]*factorial(2*n)).numerator())
return R
print A001896_list(8)
def bernoulli_list(len):
f, R, C = 1, [1], [1]+[0]*(len-1)
for n in (1..len-1):
f *= n
for k in range(n, 0, -1):
C[k] = C[k-1] / (k+1)
C[0] = -sum(C[k] for k in (1..n))
R.append(C[0]*f)
return R
print bernoulli_list(17)
P = plot([1/2*x^2 - 1/3*x, -3/4*x^3 + x^2 - 1/4*x, 3/2*x^4 - 3*x^3 + 5/3*x^2
-1/5*x, -15/4*x^5 + 10*x^4 - 35/4*x^3 + 8/3*x^2 - 1/6*x, 45/4*x^6 - 75/2*x^5
+ 45*x^4 -137/6*x^3 + 17/4*x^2 - 1/7*x], (x,0,1))
show(P, frame=True, title=
'Polynomials with p_n(0) = 0 and p_n(1) = B(n) = BernoulliNumber(n)')
def PInverseSequence(n, f, norm = None):
M = PMatrix(n, f, norm, inverse = true)
return [M[j][1] for j in (1..n-1)]
print PInverseSequence(8, lambda n: 1)
print PInverseSequence(8, lambda n: n)
print PInverseSequence(8, lambda n: n^2)
print PInverseSequence(8, lambda n: factorial(n))
eul = lambda n: 1/((2*n-1)*(2*n))
nrm = lambda n,k: factorial(2*n)/4^k
PMatrix(7, eul, nrm, inverse=True)
import time
def timeit(statement):
durations = []
for i in range(3):
start = time.time()
eval(statement)
end = time.time()
t = end - start
durations.append(t)
return min(durations)
def benchmark(stmt):
t = timeit(stmt)
r = round(t*1000, 3)
print "%s best time: %s ms" % (stmt, r)
for n in range(10,50,10):
benchmark("PartitionTransform(%d, lambda n: n)" % n)
for n in range(10,50,10):
benchmark("PMatrix(%d, lambda n: n)" % n)
def PMultiCoefficients(dim, norm = None, inverse = False):
def coefficient(p): # admittedly somewhat obscure
c = SR(p).fraction(ZZ).numerator().coefficients()
return [0] if not c else c
f = lambda n: var('x'+str(n))
P = PMatrix(dim, f, norm, inverse)
return [[coefficient(p) for p in L] for L in P]
PMultiCoefficients(8)
PMultiCoefficients(6, inverse = True)
stirset2 = lambda n: 1 if n == 1 else 1/(n*(4*n-2))
stircycle2 = lambda n: 1 if n == 1 else (n-1)^2/(n*(4*n-2))
lah2 = lambda n: 1 if n == 1 else ((n-1)^2+1)/(n*(4*n-2))
norm = lambda n,k: (-1)^k*factorial(2*n)/factorial(2*k)
M = PMatrix(7, stirset2, norm)
for m in M: print m
M = PMatrix(7, stirset2, norm, inverse = True)
for m in M: print m
M = PMatrix(7, stircycle2, norm) # A204579
for m in M: print m
M = PMatrix(7, stircycle2, norm, inverse = True)
for m in M: print m
M = PMatrix(7, lah2, norm)
for m in M: print m
M = PMatrix(7, lah2, norm, inverse = True)
for m in M: print m
def T(n, k, w):
if n==k: return 1
if k<0 or k>n: return 0
return T(n-1, k-1, w) + w(n,k)*T(n-1, k, w)
stirset3 = lambda n,k: k^3
stircycle3 = lambda n,k: (n-1)^3
lah3 = lambda n,k: stircycle3(n,k)+stirset3(n,k)
for n in range(8): print [T(n, k, stirset3) for k in (0..n)]
for n in range(8): print [T(n, k, stircycle3) for k in (0..n)]
for n in range(8): print [T(n, k, lah3) for k in (0..n)]
StirlingSetStar = lambda n: 1/(n+1)
StirlingCycleStar = lambda n: n/(n+1)
norm = lambda n, k: (-1)^k*factorial(2*n)
\begin{equation} \genfrac\{\}{0pt}{}{n}{k}^{\star} = \frac{(2n)!}{(n+k)^{\underline{n}}} \sum_{m=0}^{k} (-1)^{m+k} \binom{n+k}{n+m}\genfrac\{\}{0pt}{}{n+m}{m} \end{equation}
\begin{equation} \genfrac{ [ }{ ] }{0pt}{}{n}{k}^{\star} = \frac{(2n)!}{(n+k)^{\underline{n}}} \sum_{m=0}^{k} (-1)^{m+k} \binom{n+k}{n+m}\genfrac{[ }{ ] }{0pt}{}{n+m}{m} \end{equation}
PMatrix(7, StirlingSetStar, norm)
T = lambda n, k: sum((-1)^(m+k)*binomial(n+k, n+m)*stirling_number2(n+m, m) for m in (0..k))
N = lambda n, k: factorial(2*n)/falling_factorial(n+k,n)
for n in (0..6): print [T(n,k)*N(n,k) for k in (0..n)]
PMatrix(7, StirlingCycleStar, norm)
T = lambda n, k: sum((-1)^(m+k)*binomial(n+k, n+m)*stirling_number1(n+m, m) for m in (0..k))
N = lambda n, k: factorial(2*n)/falling_factorial(n+k,n)
for n in (0..6): print [T(n,k)*N(n,k) for k in (0..n)]
\begin{gather} (n-k+1)^{\overline{n-k}} \, \genfrac\{\}{0pt}{}{n}{k} &= \binom{n}{k} \sum_{i=0}^{k} \binom{k}{i} \genfrac\{\}{0pt}{}{n-k}{i}^{\star} &= \binom{-k}{-n} \sum_{i=0}^{n-k} \binom{-n}{i} \genfrac{[ }{ ] }{0pt}{}{n-k}{i}^{\star} \\ (n-k+1)^{\overline{n-k}} \, \genfrac{ [ }{ ] }{0pt}{}{n}{k} &= \binom{n}{k} \sum_{i=0}^{k} \binom{k}{i} \genfrac{ [ }{ ] }{0pt}{}{n-k}{i}^{\star} &= \binom{-k}{-n} \sum_{i=0}^{n-k} \binom{-n}{i} \genfrac\{\}{0pt}{}{n-k}{i}^{\star} \end{gather}