Hari Sundar
Geometric Multigrid for Higher-order Discretizations

Institute for Computational Engineering & Sciences
University of Texas at Austin
201 East 24th St
Stop C0200
Austin
Texas 78712-1229
hari@ices.utexas.edu
Georg Stadler
Omar Ghattas
George Biros

We are interested on asymptotically optimal-- $\mathcal{O}(N)$ --complexity solvers for approximating the solution of elliptic partial differential equations (PDEs), where $N$ is the number of unknowns. Multigrid is such a solver. In practice however, multigrid performs best for low-order uniform discretizations with smooth coefficients.

Our goal is to develop a parallel geometric multigrid for solving systems arising from higher-order discretizations of variable-coefficient elliptic partial differential equations on arbitrary geometries using highly adapted meshes. High-order discretizations offer several advantages. According to standard isoparametric polynomial approximation theory, by using a finite element basis of at least degree $p$ , we can achieve very fast $\mathcal{O}(N^{-(p+1)})$ convergence for sufficiently smooth problems while improving the locality and thus the CPU efficiency of the calculations.

Our method is designed for meshes that are built from an unstructured hexahedral macro mesh, in which each macro element is adaptively refined as an octree. This forest-of-octrees approach enables us to generate meshes for complex geometries with arbitrary levels of local refinement. We use geometric multigrid (GMG) for each of the octrees and algebraic multigrid (AMG) as the coarse grid solver. We designed our GMG sweeps to entirely avoid collectives, thus minimizing communication cost. Recently, we presented weak and strong scaling results for the 3D variable-coefficient Poisson problem using linear discretization that demonstrate high parallel scalability. Here we explore various approaches for extending our geometric multigrid solver to support higher-order discretizations.

For higher-order finite-element discretizations, the following approaches are commonly used,

Schwarz-based methods

This is the most common approach for solving systems arising from higher-order discretizations, which consists of local block solves and a coarse-grid solve. The main challenge with these approaches is they require solving dense local blocks with direct methods. Additionally, the coarse-grid solve can become fairly expensive and is not straightforward to achieve good parallel scalability.

$p$ -multigrid

These methods are a direct extension of multigrid concepts to higher-order discretization. The usual approach has been to use coarse grids based on lower-order polynomials. Starting from a fine grid with order-$p$ polynomial basis, the coarser grids correspond to polynomials of order $p/2, p/4,\ldots,1$ , followed by geometric coarsening of the $p=1$ grid. The main shortcoming of this approach has been the dependence of the convergence factor on the order of the polynomial basis.

precondition using lower-order operator

This approach preconditions the higher-order operator using a lower-order operator obtained by overlaying the higher-order nodes with a lower-order (typically linear) finite element mesh. Multigrid is used to solve the lower-order operator. Although this approach is nearly independent of $p$ , and is relatively straightforward to parallelize, it is not work optimal and the convergence factors are lower than multigrid applied directly to the higher-order operator.

direct

Directly apply multigrid to iteratively solve the linear systems that result from the higher-order discretizations. This approach is more difficult to implement and the cost per iteration increases. However, our preliminary results suggest that such an approach is the most general.

In summary, the overall theme of existing work appears to use low-order approximations as preconditioners. The advantages of doing this are mainly in the simplicity of the approach and the availability of parallel multigrid solvers capable of solving such lower-order operators. The sparsity of the lower-order operators also permits the use of AMG for solving the lower-order operators, possibly obtained via discretizations on unstructured meshes. Although there are examples of using Algebraic Multigrid directly on operators resulting from higher-order discretizations, limited work has been done on using geometric multigrid with higher-order discretizations. To the best of our knowledge, no prior work on using geometric multigrid for solving systems arising from higher-order discretizations on arbitrary geometries using highly adapted meshes. In this work, we develop geometric multigrid methods to support higher-order discretizations ( $1\le p\le 8$ ) and compare compare against preconditioning using the co-located linear operator. We evaluate using variable-coefficient Poisson problems on $2D$ and $3D$ domains. We demonstrate that by using appropriate inter-grid transfer operators and smoothers, mesh-independent convergence is possible ( $1\le p\le 8$ ) for the direct approach. For the direct approach, best results are obtained using the symmetric successive over-relaxation (SSOR) smoother. We conclude with thoughts on the parallelization of the proposed approach.