Software Design


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.

Data Types

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

Operation Types

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

  1. Collocation

  2. 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

Formulation Assumptions

  • 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)

Data & Object layout

SELF Memory Management

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.

Array indexing and data layout

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

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)





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.

Derived Types for Spectral Element Data

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.