Programming a fast, efficient, parallel multigrid solver for Helmholtz and other equations on dynamically adaptive meshes can't be indefinitely hard

Tobias Weinzierl (The University of Durham)


Multigrid algorithms are among the fastest linear equation system solvers in
the world. Geometric matrix-free multigrid, i.e. multigrid that derives its
coarse grids from a mesh hierarchy but never assembles
any matrix, seems so be a almost a silver bullet: no assembly overhead,
structured memory access and built-in domain decomposition all help to make it
perform, while, on Cartesian grids, the operators' tensor product structure
and straightforward adaptive mesh refinement help us to tackle moderately
high-dimensional problems. Unfortunately, it is not robust--even simple
convection-diffusion equations, jumping coefficients or small non-diffusive
operator parts as they arise with Helmholtz equations make it
struggle - and its complex data access patterns render the programming of an
efficient implementation challenging.

We introduce a family of spacetree-based multigrid realisations that rely on
the tree's multiscale nature to derive coarse grids. They fall into the class
of matrix-free geometric multigrid solvers. The most sophisticated realisation variants from the family pick up
the concept of BoxMG, which defines operator-dependent prolongation and
restriction in combination with Galerkin/Petrov-Galerkin coarse-grid operators.
This yields robust solvers. We propose to
embed the algebraic, problem- and grid-dependent multigrid operators as
stencils into the grids, and to evaluate all matrix-vector products in-situ
throughout grid traversals. While such a realisation idiom is not literally
matrix-free - the grid carries the matrix - we propose to switch to a
hierarchical representation of all operators. Only differences of algebraic
operators to their geometric counterparts are held. These hierarchical
differences can be stored and exchanged with small memory footprint. Our
realisations support arbitrary dynamically adaptive grids that can be rotated in
the complex plane, while they vertically
integrate the multilevel operations through spacetree linearisation. This
yields good memory access characteristics while standard colouring of mesh
entities with a multiscale domain decomposition allows us to use
manycore clusters; notably if we fuse coupled multi-parameter runs into one

Our techniques yield reasonably mature and robust implementation blueprints for
convection-diffusion problems with jumping coefficients. In this talk, we
however emphasise their application to a set of Helmholtz equations of moderate
dimensionality in combination with a scaling of the mesh with complex
We show that we still can construct, in most cases, fast and robust solvers.

This is joint work with Bram Reps and Marion Weinzierl.

Import this event to your Outlook calendar
▲ Up to the top