H-LU factorization on many-core systems.
The H-LU factorization is a very important algorithm
for the construction of robust H-matrix preconditioners
for BEM and FEM problems.
This paper describes how the construction of the factorization can
be implemented efficiently on a many-core system by using a
S. Börm, K. Reimer:
Efficient arithmetic operations for rank-structured matrices based on
hierarchical low-rank updates.
to appear in Comp. Vis. Sci.
Most arithmetic algorithms for H-matrices are based on a sequence
of low-rank updates that are applied to submatrices of a given
This paper introduces an algorithm that performs similar updates on
H²-matrices in linear time and shows how it can
be used to develop efficient arithmetic routines for H²-matrices.
Numerical experiments demonstrate that the new approach leads to
preconditioners for FEM and BEM problems that require
O(n log n) operations for setup and, at least for two-dimensional
problems, O(n) units of storage.
M. Bebendorf, C. Kuske, R. Venn:
Wideband nested cross approximation for Helmholtz problems.
Num. Math. (2014)
Treating medium- to high-frequency Helmholtz equations by standard
H-matrix techniques leads to very high local ranks and makes
this approach rather unattractive.
This paper combines the idea of directional decompositions introduced
in [EngquistYing2007] with the
nested ACA algorithm presented in
factoring out plane waves leads to kernel functions that are smooth
in a cone originating in a cluster, so ACA can be expected to be
appropriate for these functions.
Since plane waves are separable, the same holds for the original
The resulting algorithm uses a heuristic to choose representative
columns of the total cluster basis for the cones and reaches
a complexity of O(n log n).
S. Börm, S. Christophersen:
Approximation of integral operators by Green quadrature and nested cross
The Green quadrature method (cf.
[BörmGördes2014]) can lead to relatively high ranks,
particularly in three-dimensional applications.
This paper combines this approach with a cross approximation technique
in order to obtain a fast algorithm that construct an approximation
of a boundary integral element matrix in linear complexity.
S. Börm, J. Gördes:
Low-rank approximation of integral operators using the Green formula and quadrature.
Num. Alg. 64(3):567-592 (2013)
Boundary element methods typically rely on a representation formula
that expresses the solution of a partial differential equation by
This property can also be used to construct low-rank approximations
by applying quadrature rules to the surfaces of suitable subdomains.
This paper introduces the corresponding algorithm and proves that
it converges exponentially.
M. Bebendorf, R. Venn:
Constructing nested bases approximations from the entries of
Num. Math. (2012), 121(4):609-635
This paper presents an extension of the ACA algorithm
to H²-matrices: the total cluster bases
(cf. [Börm2005]) are replace by
a cross approximation relying on a heuristic choice of
Efficient Numerical Methods for Non-local Operators:
H²-Matrix Compression, Algorithms and Analysis.
EMS Tracts in Mathematics, Vol. 14 (2010)
H²-matrices can be used as data-sparse representations of
dense matrices resulting, e.g., from the discretization of an
integral or elliptic partial differential operator.
This book offers a comprehensive introduction into the relevant
algorithms, the techniques for the corresponding complexity
and error analysis, and a collection of numerical examples
illustrating the properties of the H²-matrix method.
L. Grasedyck, R. Kriemann, S. Le Borne:
Domain-decomposition based H-LU preconditioning.
Num. Math. 112(4):565-600 (2009)
The performance of triangular factorizations can be improved
significantly by using a nested dissection ordering.
This paper describes how to construct a nested dissection cluster
tree and use it to obtain efficient H-LU factorizations
that can be used as preconditioners for partial differential
Hierarchische Matrizen: Algorithmen und Analysis.
This book gives a comprehensive overview of the current state
of the art in the field of hierarchical matrices.
It covers algebraic and analytic compression techniques,
matrix arithmetic operations, treatment of matrix functions and
a new type of multilevel domain decomposition method based on
approximations of local Poincaré-Steklov operators.
B. Engquist, L. Ying:
Fast directional multilevel algorithms for oscillatory kernels
SIAM J. Sci. Comput., 29(4):1710-1737 (2007)
Applying H- or H²-methods directly to Helmholtz
problems with large wave numbers leads to relatively high local ranks,
which results in poor compression results.
This paper presents an alternative approach:
the Helmholtz kernel function is split into a plane wave and a factor
that is smooth in a cone originating from a given cluster.
This smooth factor can be approximated by standard techniques, and
combining the resulting expansions on different levels of the cluster
tree leads to an algorithm of complexity O(n log n).
Construction of data-sparse H2-matrices by
SIAM J. Sci. Comput., 31:1820-1839 (2009)
From an algorithmic point of view, H-matrices have a relatively
simple structure: the matrix is split into a hierarchy of submatrices,
and each submatrix is approximated by low rank.
Since all submatrices are independent, simple recursive algorithms can
be used to perform various operations.
H²-matrices introduce an additional multilevel
structure that can significantly improve the efficiency, but until now
also meant that algorithms had to preserve the connections between
different levels and blocks.
This paper introduces a simple procedure that turns an H-matrix
into an H²-matrix and has the advantage that it
can compress each submatrix as soon as it becomes available, thus
reducing the storage requirements to those of native
H²-matrix algorithms without sacrificing the
simplicity of H-matrix algorithms.
Approximation of solution operators of elliptic partial differential
equations by H- and H²-matrices
Num. Math. 115:165-193 (2010)
A major advantage of hierarchical matrix methods compared to multigrid
techniques is that they can handle jumping and anisotropic coefficients
without the need of specially adapted algorithms.
This paper improves the first approximation estimate of
and extends it to H²-matrices.
Adaptive variable-rank approximation of general dense matrices
SIAM J. Sci. Comput. 30:148-168 (2007)
The goal of variable-order schemes
[Sauter2000] is to find data-sparse
approximations of optimal complexity O(n) which ensure
that the approximation error is proportional to the discretization
The analysis of these schemes is usually quite complicated and closely
connected to the underlying continuous problem.
This paper presents an algorithm which finds variable-rank
approximations of general matrices automatically: if the given matrix
can be represented by an H²-matrix of variable rank,
the algorithmus will find such a representation without any additional
input from the user.
Data-sparse approximation of non-local operators by H²-matrices
Lin. Alg. Appl. 422:380-403 (2007)
H²-matrices reach their high efficiency by exploiting
a nested structure of the expansion systems used for the approximation of
admissible blocks. For classical integral operators, the existence of
suitable nested systems is obvious, but for more general situations, a
careful analysis is required.
This paper presents the necessary theoretical results and uses them
to prove that classical integral operators and solution operators of
elliptic partial differential equations can be approximated by
Steffen Börm, Lars Grasedyck:
Hybrid cross approximation of integral operators.
Num. Math. 101:221-249 (2005)
A popular approach to constructing low-rank approximations of discretized
integral operators is based on the application of cross approximation
to entries of the discrete matrix. This has the advantage that a complete
H-matrix approximation can be created based only on a routine for
evaluating integrals and some geometric information.
Unfortunately, this approach can be shown to fail in a number of important
situations, e.g., for the classical double layer operator on polygonal
domains. This is due to the fact that too much information is lost in the
transition from the continuous to the discrete problem.
[Börm/Grasedyck2004] proposes a different approach: the cross
approximation is applied to the original kernel function, and the resulting
kernel expansion is then discretized. The paper gives a rigorous proof
for the convergence of this approach and numerical experiments which
indicate that the hybrid technique is both faster and more reliable
than previous methods.
Steffen Börm, Lars Grasedyck, Wolfgang Hackbusch:
These are the lecture notes of the winterschool on hierarchical matrices.
We cover low-rank approximation schemes, clustering techniques, truncated
arithmetics, complexity estimates, H²-matrices and
the application of these techniques to integral operators, elliptic
partial differential equations and control theory. The lectures notes
also contain a description of the basic and some of the more advanced
functions and structures of the HLib package and a number of
theoretical and practical exercises.
Approximation of integral operators by H²-matrices
with adaptive bases
Computing 74:249-271 (2005)
The H²-matrices constructed by interpolation
typically contain a degree of redundancy, since the general polynomial
approximation scheme cannot take special properties of the geometry
of the kernel function into account. This paper presents two algorithms
that can be used to improve the efficiency: Orthogonalization detects
expansion functions that have become irrelevant during the discretization,
while recompression takes the kernel function into account.
H²-matrix arithmetics in linear complexity
Computing 77:1-28 (2006)
One of the main advantages of hierarchical matrices lies in the fact
that arithmetic operations like matrix addition, multiplication and
inversion can be performed in almost linear complexity.
"Almost" linear complexity means that additional logarithmic factors
appear in theory and practice. While this is natural for hierarchical
matrices, it is not for H²-matrices. This paper
presents algorithms that compute best approximations of sums and products
of H²-matrices in linear complexity. Numerical
experiments demonstrate the practical advantages of the new methods
compared to general H-matrix arithmetics.
Adaptive recompression of H-matrices for BEM
Computing 74(3):205-223 (2005)
The approximation constructed by ACA will typically not be the optimal
representation of a given integral operator in the hierarchical matrix
format. This paper introduces an improved ACA algorithm and a coarsening
technique that lead to a significant reduction of the storage complexity.
Mario Bebendorf, Wolfgang Hackbusch:
Existence of H-matrix approximants to the inverse FE-matrix
of elliptic operators with Loo-coefficients
Num. Math. (2003), 95:1-28
The problem of finding a good approximation of the inverse of a
finite element stiffness matrix is similar to that of finding local
low-rank approximations of the integral operator defined by the
corresponding Green function. In this paper, these approximations
are constructed for elliptic operators with non-smooth coefficients.
Lars Grasedyck, Wolfgang Hackbusch:
Construction and Arithmetics of H-Matrices
Computing (2003), 70:295-334
The complexity of the H-matrix arithmetics (addition,
multiplication, inversion) is estimated
in a general and rigorous framework. The proofs are based on standard
geometric (regular) clustering and cover the field of finite
element as well as boundary element discretisations on
non-regular (locally refined) grids in arbitrary dimension.
Lars Grasedyck, Wolfgang Hackbusch, Boris Khoromskij:
Solution of large-scale algebraic matrix Riccati equations by
use of hierarchical matrices
Computing (2003), 70:121-165
The H-matrix arithmetic can be applied to matrix equations
from control theory if the solution and all intermediate results can
be represented in the hierarchical matrix format. This paper gives
a rigorous proof of the approximability of the solutions of
algebraic matrix Riccati equations.
Mario Bebendorf, Sergej Rjasanow:
Adaptive low-rank approximation of collocation matrices
Computing (2003), 70:1-24
This paper describes a variant of the ACA algorithm and investigates
its convergence in the case of integral operators that have been
discretized by a collocation method. The paper also introduces an
partitioning strategy for matrices that reduces the number of blocks
but is no longer compatible with the hierarchical matrix structure.
Steffen Börm, Wolfgang Hackbusch:
Data-sparse approximation by adaptive H²-matrices
Computing (2002), 69:1-35
While finding the best approximation of an arbitrary matrix in a suitable
hierarchical matrix format is a straightforward application of the singular
value decomposition, the situation is much more complicated for
H²-matrices, since we can not approximate all blocks
independently, and have to ensure that the approximation spaces are nested.
This paper describes and analyzes an algorithm that finds a quasi-optimal
approximation of an arbitrary matrix in the H²-matrix
Steffen Börm, Wolfgang Hackbusch:
H²-matrix approximation of integral operators by
Appl. Num. Math. (2002), 43:129-143
The construction introduced in [Giebermann2001]
can be used to create H²-matrices. The resulting
algorithms are very fast and can be applied to arbitrary geometries and
arbitrary asymptotically smooth kernel functions.
Theorie und Anwendungen Hierarchischer Matrizen
Ph.D. thesis, accepted 2001 by the University of Kiel
This thesis gives an introduction to the theory of hierarchical
matrices, covering the approximation of integral operators, the
arithmetic operations and very general complexity estimates based
on the concepts of sparsity and idempotence of block cluster trees.
It also provides error estimates for arithmetic operations.
PDF (2822 KB)
Multilevel approximation of boundary integral operators
Computing, 67:183-207, 2001
In this paper, the author introduces a multilevel strategy for approximating
integral operators by nested interpolation. This approach is closely
related to the interpolation methods used to provide initial approximations
for current H²-matrix techniques.
Wolfgang Hackbusch, Boris Khoromskij:
A sparse H-matrix arithmetic. Part 2: Application to multi-dimensional
Computing (2000), 64:21-47
The topic of this paper is the generalization of the basic H-matrix
concepts, in particular of the complexity estimates, to the multi-dimensional
Incomplete cross approximation in the mosaic-skeleton method
Computing (2000), 64:367-380
This paper describes a method for finding low-rank approximations of
admissible blocks by looking for subspaces corresponding to maximal
volumes. The resulting approximations are closely related to the adaptive
cross approximation technique
Variable order panel clustering
Computing (2000), 64:223-261
This paper introduces and analyzes an algorithm for the approximation
of certain integral operators which ensures the optimal order of
convergence of the discretization scheme with the optimal order of
This remarkable result is due to the fact that a low-order expansion
is sufficient on geometrically small domains, while higher orders are
only required on large domains.
Wolfgang Hackbusch, Boris Khoromskij, Stefan Sauter:
Lectures on Applied Mathematics (2002), page 9-29, Springer Berlin
This is the first paper on H²-matrices, a specialized
variant of hierarchical matrices that combines the hierarchical block
structure with a hierarchy of subspaces in order to improve the efficiency
of the representation.
A sparse matrix arithmetic based on H-matrices. Part I:
Introduction to H-matrices.
Computing (1999), 62:89-108
This is the first paper ever written on hierarchical matrices. All the
basic concepts of the H-matrix arithmetic are introduced for
a simple one-dimensional model problem with rank-one approximations of
the admissible blocks.
PDF (415 KB),
S.A. Goreinov, E.E. Tyrtyshnikov, N.L. Zamarashkin:
A theory of pseudoskeleton approximations
Linear Algebra and its Applications (1997), 261:1-21
This paper demonstrates that if a matrix can be approximated by low
rank k, it can also be approximated by a matrix constructed
from k of its rows and columns (called a pseudoskeleton
This results indicates that the related cross approximation technique
can work if the correct rows and columns are chosen.
S.A. Goreinov, E.E. Tyrtyshnikov, N.L. Zamarashkin:
Pseudo-skeleton approximation by matrices of maximal volume
Mathematical Notes (1997), 62:515-519
This paper suggests a strategy for choosing the rows and columns
required for the pseudo-skeleton approximation