Understanding the Sparse Cholesky Elimination Tree

Sparse Cholesky factorization is a cornerstone of numerical linear algebra, used extensively in solving symmetric positive definite systems. While the dense Cholesky algorithm is straightforward, applying it to sparse matrices—where most elements are zero—leads to a phenomenon known as "fill-in," where zeros in the original matrix $A$ become non-zero in the resulting lower triangular matrix $L$.

To manage this efficiently, software must predict where these non-zeros will appear and organize the computation to avoid unnecessary work. This is where the Elimination Tree comes into play. It serves as a compressed representation of the task dependency graph, allowing for $O(n)$ construction and efficient symbolic factorization.

From Dense to Sparse: The Problem of Fill-in

To understand why the elimination tree is necessary, we must first look at the right-looking dense Cholesky algorithm. The core process involves factoring the pivot, scaling the column below it, and performing a rank-1 update on the trailing matrix:

for (int k = 0; k < n; ++k) {
    L[k][k] = sqrt(A[k][k]);
    for (int i = k + 1; i < n; ++i) {
        L[i][k] /= L[k][k];
    }
    for (int j = k + 1; j < n; ++j) {
        for (int i = j; i < n; ++i) {
            L[i][j] -= L[i][k] * L[j][k];
        }
    }
}

In a sparse context, the third step—the rank-1 update—is the only operation that creates fill-in. Specifically, if $k < j \le i$, and both $L[i,k]$ and $L[j,k]$ are non-zero, then $L[i,j]$ must also be non-zero. This structural rule implies that many of the operations in a dense factorization are redundant when $A$ is sparse. If we prune the task dependency graph of all unnecessary operations, we are left with a streamlined DAG (Directed Acyclic Graph).

The Column Elimination Tree

While the pruned task graph is a DAG, the structural rule mentioned above allows us to simplify it further. If we have dependencies $0 \to 1$ and $0 \to 2$, and we know that $0 \to 1$ implies $L[1,0] \neq 0$ and $0 \to 2$ implies $L[2,0] \neq 0$, the structural rule dictates that $L[2,1]$ must be non-zero, creating an implied edge $1 \to 2$. Therefore, the edge $0 \to 2$ is redundant.

By removing all such redundant edges, the DAG collapses into a tree: the Column Elimination Tree. This tree tells us exactly how to generate the fill-pattern of $L$ from the initial pattern of $A$.

Practical Implementation: Symbolic and Numeric Factorization

Sparse factorization is typically split into two phases: symbolic and numeric.

1. Symbolic Factorization

Symbolic factorization determines the non-zero structure of $L$ without performing any floating-point arithmetic. Using the elimination tree (stored as a parents array), we can traverse the paths from the non-zeros of $A$ up to the root of the tree to identify all fill-in entries.

void symbolic_cholesky(int n, int** A_row, int* A_row_len, int* parent, int* mark, vector_int* L_col) {
    for (int k = 0; k < n; ++k) {
        mark[k] = -1;
        vector_clear(&L_col[k]);
        vector_push(&L_col[k], k);
    }

    for (int i = 0; i < n; ++i) {
        for (int p = 0; p < A_row_len[i]; ++p) {
            int k = A_row[i][p];
            while (k != -1 && k < i && mark[k] != i) {
                mark[k] = i;
                vector_push(&L_col[k], i);
                k = parent[k];
            }
        }
    }
}

2. Numeric Factorization

With the non-zero pattern preallocated, the numeric phase performs the actual calculations. Because the structure is known, the algorithm no longer needs to check if an entry exists before updating it, significantly increasing performance.

Computing the Elimination Tree

To compute the elimination tree efficiently in $O(n)$ time, we process the rows of $A$ in increasing order. For each non-zero entry $A[r,c]$ (where $c < r$), we walk upward from $c$ through the previously discovered ancestors until we find a column that does not yet have a parent assigned for the current row $r$. This column's parent then becomes $r$.

This approach ensures that the tree is built incrementally, grounding the graph theory in the actual mechanics of the linear algebra algorithm.

Conclusion

The elimination tree is more than a theoretical construct; it is a practical tool that transforms the complex problem of fill-in prediction into a simple tree traversal. By distilling the task DAG of a Cholesky factorization into a tree, developers can implement high-performance sparse solvers that minimize memory overhead and maximize computational efficiency.