A two-scale approach for efficient on-the-fly operator assembly in massively parallel high performance multigrid codes

TL;DR

This paper introduces a two-scale method for efficient on-the-fly operator assembly in massively parallel multigrid codes, improving performance for non-polyhedral domains by approximating operators with surrogate polynomials.

Contribution

The paper proposes a novel two-scale approach that leverages polynomial approximations for operator assembly, enhancing efficiency over traditional matrix-free methods in complex geometries.

Findings

Significant speedup in operator assembly for non-polyhedral domains.

Effective use of surrogate polynomials for on-the-fly stencil computation.

Strong and weak scaling results demonstrate improved performance in large-scale PDE solvers.

Abstract

Matrix-free finite element implementations of massively parallel geometric multigrid save memory and are often significantly faster than implementations using classical sparse matrix techniques. They are especially well suited for hierarchical hybrid grids on polyhedral domains. In the case of constant coefficients all fine grid node stencils in the interior of a coarse macro element are equal. However, for non-polyhedral domains the situation changes. Then even for the Laplace operator, the non-linear element mapping leads to fine grid stencils that can vary from grid point to grid point. This observation motivates a new two-scale approach that exploits a piecewise polynomial approximation of the fine grid operator with respect to the coarse mesh size. The low-cost evaluation of these surrogate polynomials results in an efficient stencil assembly on-the-fly for non-polyhedral domains that can be significantly more efficient than matrix-free techniques that are based on an element-wise assembly. The performance analysis and additional hardware-aware code optimizations are based on the Execution-Cache-Memory model. Several aspects such as two-scale a priori error bounds and double discretization techniques are presented. Weak and strong scaling results illustrate the benefits of the new technique when used within large scale PDE solvers.