The Spectral Element Libraries in Fortran is a library that provides data structures and type bound procedure that can be used to implement PDE solvers using Spectral Element Methods.
SELF follows a pattern driven design where routines that are provided are determined by the types of data and the types of operations we usually need to do to implement a PDE solver. Additionally, the SELF specification indicates that routines are provided for serial CPU only, parallel CPU only, single GPU, and multi-GPU platforms. Thus, for every algorithm that must be implemented, back-end implementations exist for each target system.
In the current SELF specification, the following data types are
Scalars in 1, 2, and 3 dimensions
Vectors in 2 and 3 dimensions
Vectors in 3 dimensions
Tensors in 2 and 3 dimensions
By "Scalar", "Vector", and "Tensor", we mean their mathematical definition; all data types are stored as multi-dimensional arrays
When developing PDE solvers using spectral element methods, there are a number of common routines that are needed
Calculus Operations
Derivative (1-D)
Divergence (2-D,3-D)
Gradient (2-D, 3-D)
Curl (2-D, 3-D)
Grid Interpolation
Boundary Interpolation
For the Calculus Operations, SELF currently offers two discrete formulations that offer spectral accuracy
Collocation
Discontinuous Galerkin
*Continuous Galerkin implementations are planned for Winter 2021.
SELF stores mesh connectivity information following the HOPr Mesh Format. Because of this, meshes in SELF are always expressed as unstructured meshes. Further, elements in the mesh have isoparametric geometry; each element's internal geometry is approximated with high order polynomials
Computational coordinates always vary between [-1,1]
Nodal (not "Modal") Spectral Element Methods are assumed
All elements are "tensor-product" elements.
Elements in 2-D have a computational grid spanning [-1,1]x[-1,1] (Logically quadrilaterals)
Elements in 3-D have a computational grid spanning [-1,1]x[-1,1]x[-1,1] (Logically hexahedrals)
At the core of SELF is a memory management module (SELF_Memory.F90) that defines data structures for various sized multi-dimensional arrays. Each data structure provides a host (CPU-side) array, a device (GPU-side) pointer, constructor routine, destructor routine, host-to-device copy, and device-to-host. Other types implemented in SELF use the SELF_Memory module to manage memory.
In spectral element methods, data is organized into elements on an unstructured grid. Within each element, functions are approximated by Lagrange interpolating polynomials. The nodal points, where the function values are assigned, are arranged in a structured grid. The size of the grid is determined by the degree of the polynomial.
We use the following indexing notation conventions
i, j, k - Indices that vary from 0 to N, where N is the polynomial degree in an element.
iVar - Index used to index over multiple variables/functions.
iEl - Index used to index over elements in the spectral element mesh
iSide - Index used to index over the sides/boundaries of each element
row,col - Indices used to index over the rows and columns of tensor functions
idir - Index used to index over spatial dimensions in vector functions
The mathematical data types are stored in arrays/pointers using the following ordering. In this description, we use Fortran-syntax, meaning that the first index varies in order in memory.
Scalars 1-D : f(i,iVar,iEl), bf(iVar,iSide,iEl)
Scalars 2-D : f(i,j,iVar,iEl), bf(i,iVar,iSide,iEl)
Scalars 3-D : f(i,j,k,iVar,iEl), bf(i,j,iVar,iSide,iEl)
Vectors 2-D : f(idir,i,j,iVar,iEl), bf(idir,i,iVar,iSide,iEl)
Vectors 3-D : f(idir,i,j,k,iVar,iEl), bf(idir,i,j,iVar,iSide,iEl)
Tensors 2-D : f(row,col,i,j,iVar,iEl), bf(row,col,i,iVar,iSide,iEl)
Tensors 3-D : f(row,col,i,j,k,iVar,iEl), bf(row,col,i,j,iVar,iSide,iEl)
The Lagrange_Class.F90 module provides the core set of routines that implement interpolation and calculus algorithms. The data-structure itself stores interpolation and derivative matrices that can be applied to data to achieve the desired result.
The Lagrange Class stores matrices for
Derivatives
Grid Interpolation - Interpolate from a "source grid" to a "target grid".
Boundary Interpolation - Interpolate from the a "source grid" to element boundaries. (Usually Required for Gauss Quadrature).
As an, a derivative matrix can be applied to a discrete array of values of a scalar function to return a spectrally accurate approximation of a derivative.
DO iEl = 1,nElements
DO iVar = 1,nVariables
DO i = 0, N
df(i,iVar,iEl) = 0.0_prec
DO ii = 0,N
df(i,iVar,iEl) = df(i,iVar,iEl) + myPoly % dMatrix % hostData(ii,i)*f(ii,iVar,iEl)
END DO
END DO
END DO
END DO
Additionally, the Lagrange Class stores quadrature weights for handling integration. Quadratures currently supported are Legendre-Gauss, Legendre-Gauss-Lobatto, Chebyshev-Gauss, and Chebyshev-Gauss-Lobatto.
SELF provides Fortran derived types for storing data for scalars, vectors, and tensors in one, two, and three dimensions. These derived types provide support for CPU and GPU data storage and include routines for host and device updates. Each object stores function data on element interiors and element boundaries to facilitate the calculation of numerical fluxes (DG methods) or the coupling of elements (CG methods).
Further, each SELF_Data and SELF_MappedData type binds a Lagrange_Class object and provides easy to use routines for interpolation and derivative (div, grad, curl) operations. The SELF_MappedData type extends SELF_Data to provide the same routines but with metric terms included in the calculations.