Parallel GPU-Accelerated Randomized Construction of Approximate Cholesky Preconditioners

TL;DR

The paper addresses efficient preconditioning for large sparse Laplacian systems by introducing ParAC, a parallel, randomized incomplete Cholesky method. It combines dynamic dependency tracking with a randomized fill-in strategy to expose substantial parallelism on both CPU and GPU, while avoiding costly nested-dissection preprocessing. The CPU implementation uses a left-looking approach with a preallocated triangular buffer, whereas the GPU employs right-looking updates with a persistent kernel and a block-based hash map to manage irregular fill-ins. Experimental results show competitive performance against HyPre, AmgX, and cuSPARSE ichol across diverse matrices, with notable speedups and robustness to changing sparsity patterns, highlighting practical impact for PDEs, spectral graph tasks, and graph-based learning.

Abstract

We introduce a parallel algorithm to construct a preconditioner for solving a large, sparse linear system where the coefficient matrix is a Laplacian matrix (a.k.a., graph Laplacian). Such a linear system arises from applications such as discretization of a partial differential equation, spectral graph partitioning, and learning problems on graphs. The preconditioner belongs to the family of incomplete factorizations and is purely algebraic. Unlike traditional incomplete factorizations, the new method employs randomization to determine whether or not to keep fill-ins, i.e., newly generated nonzero elements during Gaussian elimination. Since the sparsity pattern of the randomized factorization is unknown, computing such a factorization in parallel is extremely challenging, especially on many-core architectures such as GPUs. Our parallel algorithm dynamically computes the dependency among row/column indices of the Laplacian matrix to be factorized and processes the independent indices in parallel. Furthermore, unlike previous approaches, our method requires little pre-processing time. We implemented the parallel algorithm for multi-core CPUs and GPUs, and we compare their performance to other state-of-the-art methods.

Figures (4)

• Figure 1: An overview linking matrix sparsity, graph representations, and factorization. Top left: sparse pattern of a sample Laplacian matrix. Top right: its Cholesky factor, highlighting fill-ins from elimination. Bottom left: the corresponding graph with edges connecting vertices. Bottom right: the e-tree from the factorization order. Each vertex’s removal and fill-ins create new connections among neighbors, as captured by the tree governing factorization dependencies.
• Figure 2: A possible graph after eliminating vertex 0 in \ref{['fig:mat_pattern']} using randomized Cholesky. Instead of forming a clique around the neighbors of 0, a spanning tree is formed.
• Figure 3: Factor scaling time of three different orderings on CPU, all in seconds. We show the scaling results for all three orderings.
• Figure 4: Top figure shows e-tree depth using the classical e-tree computation vs. actual e-tree height vs. triangular solve critical path length. Bottom figure shows the corresponding time usage by each ordering on GPU, and the ratio of fill-in in the resulting lower triangular factor. The ratio is defined as $\frac{2 * \texttt{nnz($G$)}}{\texttt{nnz($L$)}}$, where $G$ is the resulting factor and $L$ is the input.

Theorems & Definitions (3)

• definition 1: Graph Laplacian
• definition 2
• lemma 1