This site is supported by donations to The OEIS Foundation.
User:Peter Luschny/AignerTriangles
Contents
- 1 Aigner Triangles
- 1.1 An operation ZN x ZN→ ZN
- 1.2 Modified Charlier polynomials
- 1.3 Orthogonal polynomials
- 1.4 The associated Hankel matrix
- 1.5 Motzkin numbers and shifted Chebyshev polynomials
- 1.6 Catalan numbers and Chebyshev polynomials
- 1.7 Set partitions without singletons
- 1.8 Feynman diagrams and modified Hermite polynomials
- 1.9 Schroeder's third problem and Hermite polynomials
- 1.10 Jacobi Fractions
- 1.11 Weighted Motzkin paths
- 1.12 Selected bibliography
Aigner Triangles
KEYWORDS: Orthogonal polynomials, Jacobi fractions, generalized Motzkin numbers, weighted Motzkin paths, binary sequence operation.
Concerned with sequences:
A000165, A000296, A000698, A001006, A023443, A026300, A049310, A049347, A053121, A056594, A064189, A066325, A098742, A099174, A104562, A123023, A123023, A126120, A137286, A137338, A216916, A217537, A217538.
An operation ZN x ZN→ ZN
Let us denote the nonnegative numbers by N and the integers by Z. Given two sequences
- s: N → Z
- t: N → Z
we associate a lower triangular matrix A, build recursively by
- A(0,0) = 1; A(0,k) = 0;
- A(n,k) = A(n-1,k-1)+s(k)*A(n-1,k)+t(k+1)*A(n-1,k+1).
This matrix will be called the Aigner triangle of s and t provided t(0) = 1 and t(n) <> 0 for all n.
Let dim denote the number of rows of the triangle. Than we can compute the Aigner triangle of s and t with Sage as
def ATri(s,t,dim):
T = matrix(SR,dim,dim)
for n in range(dim): T[n,n] = 1
for n in (1..dim-1):
for k in (0..n-1):
T[n,k] = T[n-1,k-1]+s(k)*T[n-1,k]+t(k+1)*T[n-1,k+1]
return T
Note that we have by definition always the diagonal equal to 1. What we are most interested in is the first column of this triangle.
def ASeq(s,t,dim): return ATri(s,t,dim).column(0)
Thus we have an operation ZN x ZN → ZN and we will call the resulting sequence the generalized Motzkin numbers corresponding to s, t.
Modified Charlier polynomials
... and indecomposable set partitions without singletons
A first example of this operation: Let
s = lambda n: n+1 t = lambda n: n+1
then ATri(s, t, 9) gives
1
1 1
3 3 1
9 12 6 1
33 51 34 10 1
135 237 193 79 15 1
609 1188 1132 584 160 21 1
2985 6381 6920 4268 1510 293 28 1
15747 36507 44213 31542 13576 3464 497 36 1
and ASeq(s, t, 12) gives the generalized Motzkin numbers corresponding to (s, t) = (n+1, n+1), which are
1, 1, 3, 9, 33, 135, 609, 2985, 15747, 88761, 531561, ...
Looking things up in the database we find the number of indecomposable set partitions of [1..n+2] without singletons (in our enumeration), Knuth's A098742, and the triangle A216916.
Since the Aigner triangle of s and t is invertible we may also consider the inverse matrix, which is again a lower triangular array with diagonal equal to 1.
def IATri(s,t,dim):
M = ATri(s,t,dim)
return M.inverse()
For example IATri(s, t, 9) with the above s and t gives
1
-1 1
0 -3 1
3 6 -6 1
-12 -9 26 -10 1
45 3 -109 71 -15 1
-198 81 501 -475 155 -21 1
1071 -786 -2663 3329 -1455 295 -28 1
-6984 6711 16510 -25495 13729 -3647 511 -36 1
This triangle has a nice interpretation: it gives the coefficients of the modified Charlier polynomials ordered by rising powers. In the database they were entered by R. L. Bagula, A137338.
By the way, the modified Charlier polynomials satisfy the recurrence
C(n, x) = if n > 0 then (x-n)*C(n-1, x) - n*C(n-2, x)
elif n = 0 then 1 else 0.
We call the first column of the inverse of an Aigner triangle the inverse generalized Motzkin numbers corresponding to (s, t). They can be computed by
def IASeq(s,t,dim): return IATri(s,t,dim).column(0)
In our example we find:
1, -1, 0, 3, -12, 45, -198, 1071, -6984, 53217, -462330, ..
This sequence is not yet in the database.
Orthogonal polynomials
The way we proceeded to compute the inverse triangle was laborious. First we computed the triangle ATri(s, t) and then computed the inverse. However, there is a much better way to compute IATri(s, t), namely again by recursion directly from s and t.
- A*(0,0) = 1; A*(0,k) = 0;
- A*(n+1,k) = A*(n,k-1)-s(n)A*(n,k)-t(n)A*(n-1,k).
Now we can replace the routine IATri by the more efficient one:
def IATri(s,t,dim):
T = matrix(SR,dim,dim)
for n in range(dim): T[n,n] = 1
for n in (0..dim-2):
for k in (0..n):
T[n+1,k] = T[n,k-1]-s(n)*T[n,k]-t(n)*T[n-1,k]
return T
We have already seen above that we could associate a family of polynomials (the modified Charlier polynomials) with the Aigner triangle corresponding to (s(n)=n+1, t(n)=n+1). Such an association is possible in general. To this end we define a recursive system
- p0(x) = 1;
- pn+1(x) = (x-sn)pn(x) - tnpn-1(x). (*)
Implemented with Sage:
def Poly(s, t, n, x):
if n < 0 : return 0
if n == 0 : return 1
p = (x-s(n-1))*Poly(s,t,n-1,x)-t(n-1)*Poly(s,t,n-2,x)
return expand(p)
Continuing our example:
for n in (0..7): Poly(s,t,n,x) 1 x - 1 x^2 - 3*x x^3 - 6*x^2 + 6*x + 3 x^4 - 10*x^3 + 26*x^2 - 9*x - 12 x^5 - 15*x^4 + 71*x^3 - 109*x^2 + 3*x + 45 x^6 - 21*x^5 + 155*x^4 - 475*x^3 + 501*x^2 + 81*x - 198
And these are again the modified Charlier polynomials, as we expected.
Now is there something special about polynomials as defined in (*) above? Indeed there is. It is described in a theorem due to Favard:
A sequence of polynomials pn(x) forms an orthogonal system if and only if (*) holds for pn(x) for some sequences s and t.
The associated Hankel matrix
To prove this theorem we have to introduce two more entities related with every Aigner triangle of s and t.
The first is trivial: We associate with t the sequence of partial products Tn = t0t1t1...tn and for notational convenience the diagonal matrix T, hinted by:
T_0
T_1 0
T_2
0 .
.
T can be easily computed with Sage:
def T(t, dim):
M = matrix(SR, dim, dim)
for n in range(dim):
M[n,n] = mul(t(i) for i in (0..n))
return M
Continuing our example we find:
T(t,7) [ 1 0 0 0 0 0 0] [ 0 2 0 0 0 0 0] [ 0 0 6 0 0 0 0] [ 0 0 0 24 0 0 0] [ 0 0 0 0 120 0 0] [ 0 0 0 0 0 720 0] [ 0 0 0 0 0 0 5040]
The second entity which we want to associate with an Aigner triangle is much more gripping: it is the Hankel matrix of ASeq(s, t), the generalized Motzkin numbers corresponding to s, t (see above).
def Hankel(f, dim):
H = matrix(SR, dim, dim)
for i in range(dim):
for j in range(dim):
H[i,j] = f[i+j]
return H
def Han(s, t, dim):
b = ASeq(s,t,2*dim-1) # generalized Motzkin numbers wrt s, t.
return Hankel(b, dim)
Now the technical bottom line can be simply stated as:
det Hn = det Tn for all n.
Let's check this with our example:
[Han(s,t,n).det() for n in (1..7)] [1, 2, 12, 288, 34560, 24883200, 125411328000] [T(t,n).det() for n in (1..7)] [1, 2, 12, 288, 34560, 24883200, 125411328000]
These are the superfactorials A000178.
The general bottom line can now be stated as:
A sequence of numbers are generalized Motzkin numbers corresponding to some s, t if and only if det Hn <> 0 for all n.
We will also look at the row sums and the alternating row sums of the triangles:
def Sum(L):
return add(L[k] for k in range(len(L)))
def AltSum(L):
return abs(add((-1)^k*abs(L[k]) for k in range(len(L))))
Let us summarize our findings so far:
| Triangle | Sequence | Sum | AlterSum | |
| AT | A216916 | A098742 | - | - |
| AT-1 | A137338 | - | - | - |
| Hankel: | A000178 | Poly: | modified Charlier |
Motzkin numbers and shifted Chebyshev polynomials
Now let us look at a second example, which is the most characteristic one: setting s(n) = t(n) = 1. This will give us the Urform of all the relations studied here.
s = lambda n: 1 t = lambda n: 1 ATri(s, t, 9) # A064189 A026300 [ 1 0 0 0 0 0 0 0 0] [ 1 1 0 0 0 0 0 0 0] [ 2 2 1 0 0 0 0 0 0] [ 4 5 3 1 0 0 0 0 0] [ 9 12 9 4 1 0 0 0 0] [ 21 30 25 14 5 1 0 0 0] [ 51 76 69 44 20 6 1 0 0] [127 196 189 133 70 27 7 1 0] [323 512 518 392 230 104 35 8 1] ASeq(s,t,12) # A001006 1, 1, 2, 4, 9, 21, 51, 127, 323, 835, 2188, 5798 ASum(s,t,12) # A005773 1, 2, 5, 13, 35, 96, 267, 750, 2123, 6046, 17303 IATri(s, t, 9) # A104562 [ 1 0 0 0 0 0 0 0 0] [ -1 1 0 0 0 0 0 0 0] [ 0 -2 1 0 0 0 0 0 0] [ 1 1 -3 1 0 0 0 0 0] [ -1 2 3 -4 1 0 0 0 0] [ 0 -4 2 6 -5 1 0 0 0] [ 1 2 -9 0 10 -6 1 0 0] [ -1 3 9 -15 -5 15 -7 1 0] [ 0 -6 3 24 -20 -14 21 -8 1] IASeq(s,t,12) # A049347 1, -1, 0, 1, -1, 0, 1, -1, 0, 1, -1, 0 IASum(s,t,12) # A056594 1, 0, -1, 0, 1, 0, -1, 0, 1, 0, -1, 0, [Han(s,t,n).det() for n in (1..7)] [1, 1, 1, 1, 1, 1, 1] # A000012
So we can summarize: The sequence pair s(n) = t(n) = 1 maps to the Motzkin triangle and the Motzkin numbers. The row sums are the number of directed animals of size n A005773 and the alternating sums the Riordan numbers A005043.
| Triangle | Sequence | Sum | AlterSum | |
| AT | A064189, A026300 | A001006 | A005773 | A005043 |
| AT-1 | A104562 | A049347 | A056594 | - |
| Hankel: | A000012 | Poly: | Chebyshev |
Catalan numbers and Chebyshev polynomials
Our third example gives the most important special case of the Aigner triangles. Here we set s(n) = 0 and t(n) = 1.
s = lambda n: 0 t = lambda n: 1 ATri(s, t, 9) # A053121 [ 1 0 0 0 0 0 0 0 0] [ 0 1 0 0 0 0 0 0 0] [ 1 0 1 0 0 0 0 0 0] [ 0 2 0 1 0 0 0 0 0] [ 2 0 3 0 1 0 0 0 0] [ 0 5 0 4 0 1 0 0 0] [ 5 0 9 0 5 0 1 0 0] [ 0 14 0 14 0 6 0 1 0] [14 0 28 0 20 0 7 0 1] ASeq(s,t,12) # A126120 1, 0, 1, 0, 2, 0, 5, 0, 14, 0, 42, 0 ASum(s,t,12) # A001405 1, 1, 2, 3, 6, 10, 20, 35, 70, 126, IATri(s, t, 9) # A049310 [ 1 0 0 0 0 0 0 0 0] [ 0 1 0 0 0 0 0 0 0] [ -1 0 1 0 0 0 0 0 0] [ 0 -2 0 1 0 0 0 0 0] [ 1 0 -3 0 1 0 0 0 0] [ 0 3 0 -4 0 1 0 0 0] [ -1 0 6 0 -5 0 1 0 0] [ 0 -4 0 10 0 -6 0 1 0] [ 1 0 -10 0 15 0 -7 0 1] IASeq(s,t,12) # A056594 1, 0, -1, 0, 1, 0, -1, 0, 1, 0, -1, 0, IASum(s,t,12) # A010892 1, 1, 0, -1, -1, 0, 1, 1, 0, -1, -1, 0, [Han(s,t,n).det() for n in (1..7)] [1, 1, 1, 1, 1, 1, 1] # A000012
Thus we find: The sequence pair s(n) = 0 and t(n) = 1 maps to the Catalan triangle and the aerated Catalan numbers. The alternating row sums of the inverse triangle are the Fibonacci numbers A000045.
| Triangle | Sequence | Sum | AlterSum | |
| AT | A053121 | A126120 | A001405 | A001405 |
| AT-1 | A049310 | A056594 | A010892 | A000045 |
| Hankel: | A000012 | Poly: | Chebyshev |
Set partitions without singletons
We come back to our first example, however with a small modification.
s = lambda n: n t = lambda n: n if n > 0 else 1 ATri(s, t, 9) # A217537 [ 1 0 0 0 0 0 0 0 0] [ 0 1 0 0 0 0 0 0 0] [ 1 1 1 0 0 0 0 0 0] [ 1 4 3 1 0 0 0 0 0] [ 4 11 13 6 1 0 0 0 0] [ 11 41 55 35 10 1 0 0 0] [ 41 162 256 200 80 15 1 0 0] [ 162 715 1274 1176 595 161 21 1 0] [ 715 3425 6791 7182 4361 1526 294 28 1] ASeq(s,t,12) # A000296 1, 0, 1, 1, 4, 11, 41, 162, 715, 3425, 17722, 98253, ASum(s,t,12) # A217924 1, 1, 3, 9, 35, 153, 755, 4105, 24323, 155513, IATri(s, t, 9) [ 1 0 0 0 0 0 0 0 0] [ 0 1 0 0 0 0 0 0 0] [ -1 -1 1 0 0 0 0 0 0] [ 2 -1 -3 1 0 0 0 0 0] [ -3 8 5 -6 1 0 0 0 0] [ 4 -31 0 25 -10 1 0 0 0] [ -5 119 -56 -95 70 -15 1 0 0] [ 6 -533 455 364 -455 154 -21 1 0] [ -7 2904 -3326 -1428 3059 -1428 294 -28 1] IASeq(s,t,12) # (-1)^(n+1)*A023443(n) 1, 0, -1, 2, -3, 4, -5, 6, -7, 8, -9, 10, -11, 12, -13, IASum(s,t,12) # A165900 1, 1, -1, -1, 5, -11, 19, -29, 41, -55, 71, -89, 109, [Han(s,t,n).det() for n in (1..7)] # A000178 [1, 1, 2, 12, 288, 34560, 24883200]
| Triangle | Sequence | Sum | AlterSum | |
| AT | A217537 | A000296 | A217924 | A000012 |
| AT-1 | - | A023443 | A165900 | - |
| Hankel: | A000178 | Poly: | P |
Feynman diagrams and modified Hermite polynomials
Next let's see what happens if we make the choice s(n) = 0 and t(n) = n+1.
s = lambda n: 0 t = lambda n: n+1 ATri(s, t, 9) # A217538 [ 1 0 0 0 0 0 0 0 0] [ 0 1 0 0 0 0 0 0 0] [ 2 0 1 0 0 0 0 0 0] [ 0 5 0 1 0 0 0 0 0] [ 10 0 9 0 1 0 0 0 0] [ 0 37 0 14 0 1 0 0 0] [ 74 0 93 0 20 0 1 0 0] [ 0 353 0 193 0 27 0 1 0] [ 706 0 1125 0 355 0 35 0 1] Aseq(s,t,13) # aerated A000698 1, 0, 2, 0, 10, 0, 74, 0, 706, 0, 8162, 0, 110410, Asum(s,t,13) # 1, 1, 3, 6, 20, 52, 188, 574, 2222, 7626, 31350, IATri(s, t, 9) # A137286 [ 1 0 0 0 0 0 0 0 0] [ 0 1 0 0 0 0 0 0 0] [ -2 0 1 0 0 0 0 0 0] [ 0 -5 0 1 0 0 0 0 0] [ 8 0 -9 0 1 0 0 0 0] [ 0 33 0 -14 0 1 0 0 0] [ -48 0 87 0 -20 0 1 0 0] [ 0 -279 0 185 0 -27 0 1 0] [ 384 0 -975 0 345 0 -35 0 1] IASeq(s,t,12) # aerated A000165 1, 0, -2, 0, 8, 0, -48, 0, 384, 0, -3840, 0 IASum(s,t,13) # 1, 1, -1, -4, 0, 20, 20, -120, -280, 800, 3600 [Han(s,t,n).det() for n in (1..7)] # A000178(n+1) [1, 2, 12, 288, 34560, 24883200, 125411328000]
The most intriguing observation is the connection with A000698. The inverse triangle is A137286, the coefficients of a variant of the Hermite polynomials. The before missing Aigner triangle is now A217538.
| Triangle | Sequence | Sum | AlterSum | |
| AT | A217538 | A000698 | - | - |
| AT-1 | A137286 | A000165 | - | A000932 |
| Hankel: | A000178 | Poly: | modified-Hermite |
Schroeder's third problem and Hermite polynomials
A variant of the foregoing is:
s = lambda n: 0 t = lambda n: n if n > 0 else 1 ATri(s, t, 9) # A099174 [ 1 0 0 0 0 0 0 0 0] [ 0 1 0 0 0 0 0 0 0] [ 1 0 1 0 0 0 0 0 0] [ 0 3 0 1 0 0 0 0 0] [ 3 0 6 0 1 0 0 0 0] [ 0 15 0 10 0 1 0 0 0] [ 15 0 45 0 15 0 1 0 0] [ 0 105 0 105 0 21 0 1 0] [105 0 420 0 210 0 28 0 1] ASeq(s,t,13) # A123023 1, 0, 1, 0, 3, 0, 15, 0, 105, 0, 945, 0, 10395, ASum(s,t,13) # A000085 1, 1, 2, 4, 10, 26, 76, 232, 764, 2620, 9496, IATri(s, t, 9) # A066325 [ 1 0 0 0 0 0 0 0 0] [ 0 1 0 0 0 0 0 0 0] [ -1 0 1 0 0 0 0 0 0] [ 0 -3 0 1 0 0 0 0 0] [ 3 0 -6 0 1 0 0 0 0] [ 0 15 0 -10 0 1 0 0 0] [ -15 0 45 0 -15 0 1 0 0] [ 0 -105 0 105 0 -21 0 1 0] [ 105 0 -420 0 210 0 -28 0 1] IASeq(s,t,12) # A123023 1, 0, -1, 0, 3, 0, -15, 0, 105, 0, -945, 0, 10395, IASum(s,t,13) # A001464 1, 1, 0, -2, -2, 6, 16, -20, -132, 28, 1216, 936, [Han(s,t,n).det() for n in (1..7)] # A000178 [1, 1, 2, 12, 288, 34560, 24883200, 125411328000]
| Triangle | Sequence | Sum | AlterSum | |
| AT | A099174 | A123023 | A000085 | A000085 |
| AT-1 | A066325 | A123023 | A001464 | A000085 |
| Hankel: | A000178 | Poly: | Hermite |
Jacobi Fractions
Time to switch back to some notions. The Jacobi fractions are continued fractions of the form
More on the Digital Library of Mathematical Functions: DLMF
Let's raise the question: What happens if we feed our sequence pair s and t into the Jacobi fractions?
# Jcf evaluates a continued fraction of Jacobi type
# using backwards recurrence.
def Jcf(s,t,n):
c = 0
for k in range(n,0,-1):
a = 1-s(k)*x
b = t(k)*x^2
c = b/(a-c)
c = 1/(1-s(0)*x-c)
return c.simplify_full()
- The Catalan example (s(n)=0, t(n)=1) gives us:
for n in (0..6): Jcf(s,t,n) 1 -1/(x^2 - 1) (x^2 - 1)/(2*x^2 - 1) -(2*x^2 - 1)/(x^4 - 3*x^2 + 1) (x^4 - 3*x^2 + 1)/(3*x^4 - 4*x^2 + 1) -(3*x^4 - 4*x^2 + 1)/(x^6 - 6*x^4 + 5*x^2 - 1) (x^6 - 6*x^4 + 5*x^2 - 1)/(4*x^6 - 10*x^4 + 6*x^2 - 1)
Next let's see what we can recover if we expand these rational functions.
for n in (0..6): n, taylor(Jcf(s,t,n),x,0,2*n) (0) 1 (1) 1 + 1*x^2 (2) 1 + 1*x^2 + 2*x^4 (3) 1 + 1*x^2 + 2*x^4 + 5*x^6 (4) 1 + 1*x^2 + 2*x^4 + 5*x^6 + 14*x^8 (5) 1 + 1*x^2 + 2*x^4 + 5*x^6 + 14*x^8 + 42*x^10 (6) 1 + 1*x^2 + 2*x^4 + 5*x^6 + 14*x^8 + 42*x^10 + 132*x^12
- In the Motzkin case (s(n)=1, t(n)=1) we find:
for n in (0..5): Jcf(s,t,n) 1/(x - 1) (x - 1)/(2*x - 1) -(2*x - 1)/(x^3 + x^2 - 3*x + 1) -(x^3+x^2-3*x+1)/(x^4-2*x^3-3*x^2+4*x-1) (x^4-2*x^3-3*x^2+4*x-1)/(4*x^4-2*x^3-6*x^2+5*x-1) -(4*x^4-2*x^3-6*x^2+5*x-1)/(x^6+2*x^5-9*x^4+10*x^2-6*x+1) and expanded: (0) 1 (1) 1+x+2*x^2 (2) 1+x+2*x^2+4*x^3+9*x^4 (3) 1+x+2*x^2+4*x^3+9*x^4+21*x^5+51*x^6 (4) 1+x+2*x^2+4*x^3+9*x^4+21*x^5+51*x^6+127*x^7+323*x^8
Thus we see: Going to infinity with n yields the infinite J-fraction expansion of Sum_{n≥0} ASeqn(s,t) zn. In other words, we found ordinary generating functions of the generalized Motzkin numbers associated with s and t (up to n).
Weighted Motzkin paths
We have seen that a great variety of combinatorial sequences arise from Aigner triangles, which can be computed recursively and have Jacobi fractions as generating functions. They give rise to orthogonal polynomials and indeed (the coefficients of) all possible orthogonal polynomials can be found from Aigner triangles.
But the icing on the cake is that there is a uniform combinatorial interpretation for all these objects. The name is weighted Motzkin paths.
A Motzkin path is a path in R2 from (0, 0) to (n, 0) consisting of a sequence of steps provided that it never falls below the x-axis. Three types of steps are permitted, denoted by (1,1), (1,0) and (1,-1), and called up-steps, horizontal steps and down-steps, respectively.
We will work here with Motzkin paths by letting 1’s represent steps in the direction (1,1) and 0’s represent steps in the direction (1,0) and -1’s represent steps in the direction (1,-1). At least when we are coding. When pretty printing a path we will use instead the letters 'H', 'D' and 'U' for '0', '-1' and '1', respectively.
def steps_to_string(S):
Z = ""
for s in S:
if s== 0: Z += 'H' # horizontal step
elif s==-1: Z += 'D' # down step
elif s== 1: Z += 'U' # up step
return Z
def PrintAllMotzkinPaths(n):
P = iter(MotzkinPaths(n))
for p in P:
print steps_to_string(p),
return
PrintAllMotzkinPaths(5)
HHHHH HHHUD HHUHD HHUDH HUHHD HUHDH HUDHH HUDUD
HUUDD UHHHD UHHDH UHDHH UHDUD UHUDD UDHHH UDHUD
UDUHD UDUDH UUHDD UUDHD UUDDH
Here we used the following elegant Sage function written by Dan Drake to generate the Motzkin paths.
# Copyright (C) 2011 Dan Drake. This program is free software:
# you can redistribute it and/or modify it under the terms of
# the GNU General Public License as published by the Free
# Software Foundation, either version 2 of the License, or (at
# your option) any later version. See: http://www.gnu.org/licenses
from sage.combinat.backtrack import GenericBacktracker
class MotzkinPaths(GenericBacktracker):
def __init__(self, n, h=None):
GenericBacktracker.__init__(self, [], (0, 0))
self._n = n
if h is not None:
self.max_ht = h
else:
self.max_ht = n
def _rec(self, path, state):
len, ht = state
if len < self._n:
# if length is less than n, we need to keep building
# the path, so new length will be 1 longer, and we
# don't yield the path yet.
newlen = len + 1
# can always tack on an E step, if the new height isn't
# more than the new number of steps remaining
if ht <= self._n - newlen:
yield path + [0], (newlen, ht), False
# if we're not touching the x-axis, we can yield a path
# with a downstep at the end
if ht > 0:
yield path + [-1], (newlen, ht - 1), False
# can take an upstep if the new height isn't more than
# the new number of steps remaining, and if we're not
# at the maximum height
if ht + 1 <= self._n - newlen and ht < self.max_ht:
yield path + [1], (newlen, ht + 1), False
else:
# if length is n, set state to None so we stop trying
# to make new paths, and yield what we've got
yield path, None, True
# - End of Dan Drake's MotzkinPaths.
Now comes the crucial step. We associate with a path a weight, depending on the function pair s and t of course.
def weight(s,t,P):
h = 0; w = 1
for p in P:
if p== 1: h += 1
elif p== 0: w *= s(h)
elif p==-1: w *= t(h); h -= 1;
return w
Basically this says: The weight of an up-step is 1, the weight of an horizontal step is given by the s-function and the weight of a down-step by the t-function where the functions are evaluated at the height of the path in this point. The weight of the path is the product of the weight of its steps.
For example we find for Motzkin paths with length 4 and
(I) s(n) = n+1 and t(n) = n+1 (first example above)
(II) s(n) = 1 and t(n) = 1 (Motzkin case)
(III) s(n) = 0 and t(n) = 1 (Catalan case)
[I] [II] [III]
--------------------------
w(HHHH) = 1 | 1 | 0
w(HHUD) = 2 | 1 | 0
w(HUHD) = 4 | 1 | 0
w(HUDH) = 2 | 1 | 0
w(UHHD) = 8 | 1 | 0
w(UHDH) = 4 | 1 | 0
w(UDHH) = 2 | 1 | 0
w(UDUD) = 4 | 1 | 1
w(UUDD) = 6 | 1 | 1
-- - -
33 9 2
The bottom line of this construction is:
The sums of the weights over all Motzkin paths of length n are the generalized Motzkin numbers associated with s and t.
Selected bibliography
- M. Aigner, Catalan and other numbers: a recurrent theme. In: Algebraic Combinatorics and Computer Science, Springer 2001. (What we called Aigner matrix here Aigner calls a 'recursive matrix'. However this term appears with different meanings in the literature.)
- P. Flajolet, Combinatorial aspects of continued fractions, Discrete Math. 32 (1980), 125-161
- C. G. J. Jacobi, Werke, vol. VI, p. 385-426.
- X. Viennot, Une théorie combinatoire des polynômes orthogonaux généraux, Lecture notes, Montréal (1984).
- H. Wronski, Philosophie de la Technie algorithmique: Loi suprême et universelle. Paris 1815-1817
