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

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.