33 research outputs found
Randomized Strong Recursive Skeletonization: Simultaneous compression and factorization of -matrices in the Black-Box Setting
The hierarchical matrix (-matrix) formalism provides a way
to reinterpret the Fast Multipole Method and related fast summation schemes in
linear algebraic terms. The idea is to tessellate a matrix into blocks in such
as way that each block is either small or of numerically low rank; this enables
the storage of the matrix and the application of it to a vector in linear or
close to linear complexity. A key motivation for the reformulation is to extend
the range of dense matrices that can be represented. Additionally,
-matrices in principle also extend the range of operations
that can be executed to include matrix inversion and factorization. While such
algorithms can be highly efficient for certain specialized formats (such as
HBS/HSS matrices based on ``weak admissibility''), inversion algorithms for
general -matrices tend to be based on nested recursions and
recompressions, making them challenging to implement efficiently. An exception
is the \textit{strong recursive skeletonization (SRS)} algorithm by Minden, Ho,
Damle, and Ying, which involves a simpler algorithmic flow. However, SRS
greatly increases the number of blocks of the matrix that need to be stored
explicitly, leading to high memory requirements. This manuscript presents the
\textit{randomized strong recursive skeletonization (RSRS)} algorithm, which is
a reformulation of SRS that incorporates the randomized SVD (RSVD) to
simultaneously compress and factorize an -matrix. RSRS is a
``black box'' algorithm that interacts with the matrix to be compressed only
via its action on vectors; this extends the range of the SRS algorithm (which
relied on the ``proxy source'' compression technique) to include dense matrices
that arise in sparse direct solvers
A distributed-memory parallel algorithm for discretized integral equations using Julia
Boundary value problems involving elliptic PDEs such as the Laplace and the
Helmholtz equations are ubiquitous in physics and engineering. Many such
problems have alternative formulations as integral equations that are
mathematically more tractable than their PDE counterparts. However, the
integral equation formulation poses a challenge in solving the dense linear
systems that arise upon discretization. In cases where iterative methods
converge rapidly, existing methods that draw on fast summation schemes such as
the Fast Multipole Method are highly efficient and well established. More
recently, linear complexity direct solvers that sidestep convergence issues by
directly computing an invertible factorization have been developed. However,
storage and compute costs are high, which limits their ability to solve
large-scale problems in practice. In this work, we introduce a
distributed-memory parallel algorithm based on an existing direct solver named
``strong recursive skeletonization factorization.'' The analysis of its
parallel scalability applies generally to a class of existing methods that
exploit the so-called strong admissibility. Specifically, we apply low-rank
compression to certain off-diagonal matrix blocks in a way that minimizes data
movement. Given a compression tolerance, our method constructs an approximate
factorization of a discretized integral operator (dense matrix), which can be
used to solve linear systems efficiently in parallel. Compared to iterative
algorithms, our method is particularly suitable for problems involving
ill-conditioned matrices or multiple right-hand sides. Large-scale numerical
experiments are presented to demonstrate the performance of our implementation
using the Julia language