Cholesky Decomposition, Combinatorics, and Fractals
Published:
Cholesky factorization is a core concept in linear algebra and frequently used to solve scientific computing problems. While seemingly the furthest thing from Pascal’s triangle and binomial coefficients there is an interesting link using results from a path finding problem. Furthermore, this can be extended to create fractals and other lovely patterns.
So first, a background in each of these areas…
Pascal’s Triangle
Pascal’s triangle is generated by setting each number in a triangle to the sum of its two parents. Beginning with ones this gives something like below.
1
1 1
1 2 1
1 3 3 1
1 4 6 4 1
1 5 10 10 5 1
Probably the most important property of Pascal’s triangle is that the number at \((n,k)\) is equal to the binomial coefficient \(\binom{n}{k}\). This pattern repeats on forever.
From here on out let \(P_n\) be the Pascal matrix with Pascal’s triangle in its lower triangle.
Cholesky Decomposition
In a completely opposite direction Cholesky Decomposition is very common linear algebra technique used in scientific computing and data analytics. For a real, symmetric matrix \(A\) Cholesky decomposition will decompose \(A\) into lower and upper triangular matrices.
\[A = LL^{T}\]As it turns out real, symmetric matrices are quite common in applications and Cholesky decomposition is faster than other techniques, such as LU decomposition.
Delannoy Numbers
Delannoy Numbers have their origin in path finding and enumeration. Their wikipedia page describes them quite well along with having helpful pictures.
Let their be a grid of points with \((0, 0)\) in the lower-left corner. Then \(D(n,k)\) is the total number of paths from \((0, 0)\) to \((n, k)\) where the only direction moved is up, right, or up-right diagonal.
This exposes a quite nice recurrence relation \(D(a,b) = D(a-1,b) + D(a,b-1) + D(a-1,b-1)\). However, there are more practical methods for calculating these numbers:
\[D(n,k) = \sum_{d=0}^{n} \binom{k}{d} \binom{n+k-d}{k} = \sum_{d=0}^{n} 2^d \binom{k}{d} \binom{n}{d}\]For convenience purposes, we’ll define the matrix \(D_n\) as an \(n\times n\) matrix where \([D_n]_{ij} = D(i,j)\). Essentially, \(D_n\) is a table containing all the Delannoy Numbers up to \(n\).
The Connection
Apart from Pascal’s triangle and Delannoy Numbers both having a relation to combinatorics, all three of these are very different from each other. Yet there is a pretty cool, in my opinion, tie-in from Delannoy Numbers to Pascal’s triangle using Cholesky decomposition.
First, we will take the Delannoy matrix \(D_n\) for some \(n\) and decompose it using Cholesky decomposition. This will give us
\[D_n = L L^T\]Now take the diagonal matrix \(B_n = \text{diag} \left( 2^{-0/2}, 2^{-1/2}, \ldots, 2^{-n/2} \right)\). Using this matrix and the upper triangular one from above we can construct the Pascal matrix:
\[P_n = L^T B_n\]This can easily be verified for a given \(n\) using a computer.
\[D_8 = \left( \begin{array}{ccccccccc} 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 \\ 1 & 3 & 5 & 7 & 9 & 11 & 13 & 15 & 17 \\ 1 & 5 & 13 & 25 & 41 & 61 & 85 & 113 & 145 \\ 1 & 7 & 25 & 63 & 129 & 231 & 377 & 575 & 833 \\ 1 & 9 & 41 & 129 & 321 & 681 & 1289 & 2241 & 3649 \\ 1 & 11 & 61 & 231 & 681 & 1683 & 3653 & 7183 & 13073 \\ 1 & 13 & 85 & 377 & 1289 & 3653 & 8989 & 19825 & 40081 \\ 1 & 15 & 113 & 575 & 2241 & 7183 & 19825 & 48639 & 108545 \\ 1 & 17 & 145 & 833 & 3649 & 13073 & 40081 & 108545 & 265729 \\ \end{array} \right)\] \[L^T = \left( \begin{array}{ccccccccc} 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & \sqrt{2} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 2 \sqrt{2} & 2 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 3 \sqrt{2} & 6 & 2 \sqrt{2} & 0 & 0 & 0 & 0 & 0 \\ 1 & 4 \sqrt{2} & 12 & 8 \sqrt{2} & 4 & 0 & 0 & 0 & 0 \\ 1 & 5 \sqrt{2} & 20 & 20 \sqrt{2} & 20 & 4 \sqrt{2} & 0 & 0 & 0 \\ 1 & 6 \sqrt{2} & 30 & 40 \sqrt{2} & 60 & 24 \sqrt{2} & 8 & 0 & 0 \\ 1 & 7 \sqrt{2} & 42 & 70 \sqrt{2} & 140 & 84 \sqrt{2} & 56 & 8 \sqrt{2} & 0 \\ 1 & 8 \sqrt{2} & 56 & 112 \sqrt{2} & 280 & 224 \sqrt{2} & 224 & 64 \sqrt{2} & 16 \\ \end{array} \right)\] \[B_n = \left( \begin{array}{ccccccccc} 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & \frac{1}{\sqrt{2}} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & \frac{1}{2} & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{2 \sqrt{2}} & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{4} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{4 \sqrt{2}} & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & \frac{1}{8} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & \frac{1}{8 \sqrt{2}} & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & \frac{1}{16} \\ \end{array} \right)\] \[P_n = L^T B_n = \left( \begin{array}{ccccccccc} 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 2 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 3 & 3 & 1 & 0 & 0 & 0 & 0 & 0 \\ 1 & 4 & 6 & 4 & 1 & 0 & 0 & 0 & 0 \\ 1 & 5 & 10 & 10 & 5 & 1 & 0 & 0 & 0 \\ 1 & 6 & 15 & 20 & 15 & 6 & 1 & 0 & 0 \\ 1 & 7 & 21 & 35 & 35 & 21 & 7 & 1 & 0 \\ 1 & 8 & 28 & 56 & 70 & 56 & 28 & 8 & 1 \\ \end{array} \right)\]While the numeric verification is not a proof, it is a cool feeling to see that Pascal triangle pop-out at the end.
From Matrices to Fractals
Now that was already interesting enough, but we can make it more interesting by turning the result into a fractal. Let’s define a new matrix \(F\) such that
\[[F]_{ij} = [L^T B_n]_{ij} \mod 2\]Now the lower triangle of \(F\) is equivalent to the Sierpinski sieve, which forms a beautiful triangular fractal.