High-order matrix-free incompressible flow solvers with GPU acceleration and low-order refined preconditioners

Abstract

We present a matrix-free flow solver for high-order finite element discretizations of the incompressible Navier-Stokes and Stokes equations with GPU acceleration. For high polynomial degrees, assembling the matrix for the linear systems resulting from the finite element discretization can be prohibitively expensive, both in terms of computational complexity and memory. For this reason, it is necessary to develop matrix-free operators and preconditioners, which can be used to efficiently solve these linear systems without access to the matrix entries themselves. The matrix-free operator evaluations utilize GPU-accelerated sum-factorization techniques to minimize memory movement and maximize throughput. The preconditioners developed in this work are based on a low-order refined methodology with parallel subspace corrections, as described for diffusion problems in [39]. The saddle-point Stokes system is solved using block-preconditioning techniques, which are robust in mesh size, polynomial degree, time step, and viscosity. For the incompressible Navier-Stokes equations, we make use of projection (fractional step) methods, which require Helmholtz and Poisson solves at each time step. The performance of our flow solvers is assessed on several benchmark problems in two and three spatial dimensions.

Conservative and accurate solution transfer between high-order and low-order refined finite element spaces

“…As a consequence, the matrix-free algorithm has significantly higher arithmetic intensity than the matrix-based algorithm. On GPU-based platforms, memory transfer is typically the bottleneck, and the matrix-free algorithms can be expected to outperform algorithms requiring fully assembled matrices [23,14,13]. The appropriate choice of algorithm will depend on both polynomial degrees p and q.…”

Accelerating High-Order Mesh Optimization Using Finite Element Partial Assembly on GPUs

“…Our default choice is the Minimum Residual (MINRES) method, as ∂ 2 F is symmetric but not necessarily positivedefinite. Preconditioning for matrix-free inversion is a substantial challenge and an active area of research [5]. We have the option to perform Jacobi preconditioning, as the diagonal of ∂ 2 F can be computed through tensor contractions without having the global matrix; these algorithms can be foung in files fem/tmop/tmop pa h2d.cpp and fem/tmop/tmop pa h3d.cpp for 2D and 3D, respectively.…”

“…Obtaining the above PA complexities, however, requires that the finite element basis functions are tensor products of 1D basis functions, e.g., quadrilaterals in 2D and hexahedra in 3D. Partial assembly has become even more relevant in recent years [4,5,6] owing to its efficient use of GPU-based architectures, which are desirable for arithmetically intensive applications that do not require a large amount of data to be moved between the CPU and GPU [1].…”

Matrix-free approaches for GPU acceleration of a high-order finite element hydrodynamics application using MFEM, Umpire, and RAJA

“…MARBL is built on modular physics and computer science components and makes extensive use of high-order finite element numerical methods. Compared to standard low-order finite volume schemes, high-order numerical methods have more resolution/accuracy per unknown and have higher FLOP/ byte ratios meaning that more floating-point operations are performed for each piece of data retrieved from memory (Dobrev et al, 2012); (Franco et al, 2020). This leads to improved strong parallel scalability, better throughput on GPU platforms, and increased computational efficiency.…”