An Optimized Multicolor Point-Implicit Solver for Unstructured Grid Applications on GPUs

Overview

Accurate and efficient simulations of complex aerodynamic flows are challenging and require significant computational resources. NASA Langley Research Center's FUN3D computational fluid dynamics (CFD) software is used to solve the Navier-Stokes (NS) equations for a broad range of aerodynamics and multidisciplinary applications across the speed range.

In this work, we have proposed an efficient implementation of the FUN3D solver for graphics processing units (GPUs). The solver was initially formulated to leverage two existing CUDA library functions. However, numerical experiments showed that the performance of these functions was suboptimal for matrices representative of those encountered in FUN3D simulations, such as those associated with unstructured grids over complex aerospace configurations. Instead, optimized versions of these functions using varying levels of parallelism have been developed. The new implementation shows performance gains of up to seven times that of the approach based on existing CUDA library functions.

Project Details

FUN3D uses a fully implicit strategy based on a finite-volume spatial discretization on mixed-element unstructured grids, which requires the frequent solution of large tightly-coupled systems of block-sparse linear equations of the form Ax=b. To develop an efficient GPU implementation of the multicolor point-implicit solver, functions provided by NVIDIA's CUDA Sparse Matrix (cuSPARSE) and CUDA Basic Linear Algebra Subroutines (cuBLAS) libraries were first considered. The function cusparseSbsrmv multiplies a block-sparse matrix with a vector, and the function cublasStrsmBatched solves block systems of equations by performing forward and backward substitutions using a previously factored coefficient matrix. These library functions proved suboptimal for linear systems representative of those encountered in FUN3D simulations.

Our proposed optimizations of these two functions include an algorithm that allocates a number of warps (or a group of 32 threads) to process a subset of the blocks in one row of the sparse matrix, in order to perform a sparse matrix-vector product. A separate kernel is also used to perform forward and backward substitutions. We encountered several challenges, including a variable extent of available parallelism, indirect memory addressing, low arithmetic intensity, and the need to accommodate different block sizes. To address these challenges, our latest implementation places particular emphasis on coalesced memory loads, the use of shared memory and pre-fetching, minimal thread divergence within warps, and strategic use of shuffle instructions available on recent hardware.

Results and Impact

All results were generated on NASA's Pleiades supercomputer, using an Intel Xeon E5-2670 dual socket, eight-core processor with an NVIDIA Tesla K40 GPU accelerator. Performance was evaluated using a series of block-sparse matrices based on tetrahedral meshes. For a given sparsity pattern, matrices were populated using random values over a range of block sizes. We identified several factors critical to the performance of the underlying matrix-vector multiplication and forward and backward substitution operations, and addressed several challenges specific to the target application. Speedups of as much as seven times that of the library-based implementation have been observed for a broad range of matrix and block sizes.

Why HPC Matters

Accurate and efficient simulations of complex aerodynamic flows require frequent solutions of billions of simultaneous equations. These intensive calculations can be efficiently performed using emerging high-performance computing hardware such as GPU accelerators.