SlabLU: a two-level sparse direct solver for elliptic PDEs

Abstract

The paper describes a sparse direct solver for the linear systems that arise from the discretization of an elliptic PDE on a two-dimensional domain. The scheme decomposes the domain into thin subdomains, or “slabs” and uses a two-level approach that is designed with parallelization in mind. The scheme takes advantage of \(\varvec{\mathcal {H}}^\textbf{2}\)-matrix structure emerging during factorization and utilizes randomized algorithms to efficiently recover this structure. As opposed to multi-level nested dissection schemes that incorporate the use of \(\varvec{\mathcal {H}}\) or \(\varvec{\mathcal {H}}^\textbf{2}\) matrices for a hierarchy of front sizes, SlabLU is a two-level scheme which only uses \(\varvec{\mathcal {H}}^\textbf{2}\)-matrix algebra for fronts of roughly the same size. The simplicity allows the scheme to be easily tuned for performance on modern architectures and GPUs. The solver described is compatible with a range of different local discretizations, and numerical experiments demonstrate its performance for regular discretizations of rectangular and curved geometries. The technique becomes particularly efficient when combined with very high-order accurate multidomain spectral collocation schemes. With this discretization, a Helmholtz problem on a domain of size \(\textbf{1000} \varvec{\lambda } \times \textbf{1000} \varvec{\lambda }\) (for which \(\varvec{N}~\mathbf {=100} \textbf{M}\)) is solved in 15 min to 6 correct digits on a high-powered desktop with GPU acceleration.