Spaces

A Space represents a discretized function space over some domain. Currently two main discretizations are supported: Spectral Element Discretization (both Continuous Galerkin and Discontinuous Galerkin types) and a staggered Finite Difference Discretization. Combination of these two in the horizontal/vertical directions, respectively, is what we call a hybrid space.

Sketch of a 2DX hybrid discretization:

3D hybrid discretization in a Cartesian domain

ClimaCore.SpacesModule
Meshes
  • domain
  • topology
  • coordinates
  • metric terms (inverse partial derivatives)
  • quadrature rules and weights

References / notes

source

Finite Difference Spaces

ClimaCore.jl supports staggered Finite Difference discretizations. Finite Differences discretize an interval domain by approximating the function by a value at either the center of each element (also referred to as cell) (CenterFiniteDifferenceSpace), or the interfaces (faces in 3D, edges in 2D or points in 1D) between elements (FaceFiniteDifferenceSpace).

Users should construct either the center or face space from the mesh, then construct the other space from the original one: this internally reuses the same data structures, and avoids allocating additional memory.

Spectral Element Spaces

ClimaCore.Spaces.SpectralElementSpace1DType
SpectralElementSpace1D(grid::SpectralElementGrid1D)
SpectralElementSpace1D(
    topology::Topologies.IntervalTopology,
    quadrature_style::Quadratures.QuadratureStyle;
    kwargs...
)
source

Extruded Finite Difference Spaces

ClimaCore.Spaces.ExtrudedFiniteDifferenceSpaceType
ExtrudedFiniteDifferenceSpace(grid, staggering)

ExtrudedFiniteDifferenceSpace(
    horizontal_space::AbstractSpace,
    vertical_space::FiniteDifferenceSpace,
    hypsography::Grids.HypsographyAdaption = Grids.Flat();
    deep::Bool = false,
)

An extruded finite-difference space, where the extruded direction is staggered, containing grid information at either

source

Utilities

ClimaCore.Spaces.areaFunction
Spaces.area(space::Spaces.AbstractSpace)

The length/area/volume of space. This is computed as the sum of the quadrature weights $W_i$ multiplied by the Jacobian determinants $J_i$:

\[\sum_i W_i J_i \approx \int_\Omega \, d \Omega\]

If space is distributed, this uses a ClimaComms.allreduce operation.

source