Meshes

A Mesh is a division of a domain into elements.

Mesh types

ClimaCore.Meshes.AbstractMeshType
AbstractMesh{dim}

A Mesh is an object which represents how we discretize a domain into elements.

It should be lightweight (i.e. exists on all MPI ranks), e.g for meshes stored in a file, it would contain the filename.

Face and vertex numbering

In 1D, faces and vertices are the same, and both are numbered [1,2].

In 2D, a face is a line segment between to vertices, and both are numbered [1,2,3,4], in a counter-clockwise direction.

 v4        f3        v3
   o-----------------o
   |                 |	    face    vertices
   |                 |	      f1 =>  v1 v2
f4 |                 | f2     f2 =>  v2 v3
   |                 |	      f3 =>  v3 v4
   |                 |        f4 =>  v4 v1
   |                 |
   o-----------------o
  v1       f1        v2

Interface

A subtype of AbstractMesh should define the following methods:

The following types/methods are provided by AbstractMesh:

source
ClimaCore.Meshes.IntervalMeshType
IntervalMesh <: AbstractMesh

A 1D mesh on an IntervalDomain.

Constuctors

IntervalMesh(domain::IntervalDomain, faces::AbstractVector)

Construct a 1D mesh with face locations at faces.

IntervalMesh(domain::IntervalDomain[, stretching=Uniform()]; nelems=)

Constuct a 1D mesh on domain with nelems elements, using stretching. Possible values of stretching are:

source
ClimaCore.Meshes.RectilinearMeshType
RectilinearMesh <: AbstractMesh2D

Constructors

RectilinearMesh(domain::RectangleDomain, n1, n2)

Construct a RectilinearMesh of equally-spaced n1 by n2 elements on domain.

RectilinearMesh(intervalmesh1::IntervalMesh1, intervalmesh2::IntervalMesh2)

Construct the product mesh of intervalmesh1 and intervalmesh2.

source
ClimaCore.Meshes.AbstractCubedSphereType
AbstractCubedSphere <: AbstractMesh2D

This is an abstract type of cubed-sphere meshes on SphereDomains. A cubed-sphere mesh has 6 panels, laid out as follows:

                                          :   Panel 1   :
                            +-------------+-------------+
                            |     +x1     |     +x1     |
                            |             |             |
                            |    Panel    |    Panel    |
                            |+x3   5   -x3|-x2   6   +x2|
                            |     -x2     |     -x3     |
                            |             |             |
                            |     -x1     |     -x1     |
              +-------------+-------------+-------------+
              |     -x2     |     -x2     |
              |             |             |
              |    Panel    |    Panel    |
              |+x1   3   -x1|+x3   4   -x3|
              |     +x3     |     -x1     |
              |             |             |
              |     +x2     |     +x2     |
+-------------+-------------+-------------+
|     +x3     |     +x3     |
|             |             |
|    Panel    |    Panel    |
|-x2   1   +x2|+x1   2   -x1|
|     +x1     |     +x2     |
|             |             |
|     -x3     |     -x3     |
+-------------+-------------+
:   Panel 6   :

This is the same panel ordering used by the S2 Geometry library (though we use 1-based instead of 0-based numering).

Elements are indexed by a CartesianIndex{3} object, where the components are:

  • horizontal element index (left to right) within each panel.
  • vertical element index (bottom to top) within each panel.
  • panel number

Subtypes should have the following fields:

  • domain: a SphereDomain
  • ne: number of elements across each panel

External links

source
ClimaCore.Meshes.EquiangularCubedSphereType
EquiangularCubedSphere <: AbstractCubedSphere

An equiangular gnomonic mesh proposed by [8]. Uses the element indexing convention of AbstractCubedSphere.

Constructors

EquiangularCubedSphere(
    domain::Domains.SphereDomain,
    ne::Integer,
    localelementmap=NormalizedBilinearMap()
    )

Constuct an EquiangularCubedSphere on domain with ne elements across each panel.

source
ClimaCore.Meshes.EquidistantCubedSphereType
EquidistantCubedSphere <: AbstractCubedSphere

An equidistant gnomonic mesh outlined in [9] and [10]. Uses the element indexing convention of AbstractCubedSphere.

Constructors

EquidistantCubedSphere(domain::Domains.SphereDomain, ne::Integer)

Constuct an EquidistantCubedSphere on domain with ne elements across each panel.

source
ClimaCore.Meshes.ConformalCubedSphereType
ConformalCubedSphere <: AbstractCubedSphere

A conformal mesh outlined in [9]. Uses the element indexing convention of AbstractCubedSphere.

Constructors

ConformalCubedSphere(domain::Domains.SphereDomain, ne::Integer)

Constuct a ConformalCubedSphere on domain with ne elements across each panel.

source

Local element map

Mesh stretching

ClimaCore.Meshes.ExponentialStretchingType
ExponentialStretching(H::FT)

Apply exponential stretching to the domain when constructing elements. H is the scale height (a typical atmospheric scale height H ≈ 7.5km).

For an interval $[z_0,z_1]$, this makes the elements uniformally spaced in $\zeta$, where

\[\zeta = \frac{1 - e^{-\eta/h}}{1-e^{-1/h}},\]

where $\eta = \frac{z - z_0}{z_1-z_0}$, and $h = \frac{H}{z_1-z_0}$ is the non-dimensional scale height. If reverse_mode is true, the smallest element is at the top, and the largest at the bottom (this is typical for land model configurations).

Then, the user can define a stretched mesh via

ClimaCore.Meshes.IntervalMesh(interval_domain, ExponentialStretching(H); nelems::Int, reverse_mode = false)

faces contain reference z without any warping.

source
ClimaCore.Meshes.GeneralizedExponentialStretchingType
GeneralizedExponentialStretching(dz_bottom::FT, dz_top::FT)

Apply a generalized form of exponential stretching to the domain when constructing elements. dz_bottom and dz_top are target element grid spacings at the bottom and at the top of the vertical column domain (m). In typical atmosphere configurations, dz_bottom is the smallest grid spacing and dz_top the largest one. On the other hand, for typical land configurations, dz_bottom is the largest grid spacing and dz_top the smallest one.

For land configurations, use reverse_mode = true (default value false).

Then, the user can define a generalized stretched mesh via

ClimaCore.Meshes.IntervalMesh(interval_domain, GeneralizedExponentialStretching(dz_bottom, dz_top); nelems::Int, reverse_mode = false)

faces contain reference z without any warping.

source
ClimaCore.Meshes.HyperbolicTangentStretchingType
HyperbolicTangentStretching(dz_surface::FT)

Apply a hyperbolic tangent stretching to the domain when constructing elements. dz_surface is the target element grid spacing at the surface. In typical atmosphere configuration, it is the grid spacing at the bottom of the vertical column domain (m). On the other hand, for typical land configurations, it is the grid spacing at the top of the vertical column domain.

For an interval $[z_0,z_1]$, this makes the elements uniformally spaced in $\zeta$, where

\[\eta = 1 - \frac{tanh[\gamma(1-\zeta)]}{tanh(\gamma)},\]

where $\eta = \frac{z - z_0}{z_1-z_0}$. The stretching parameter $\gamma$ is chosen to achieve a given resolution dz_surface at the surface.

Then, the user can define a stretched mesh via

ClimaCore.Meshes.IntervalMesh(interval_domain, HyperbolicTangentStretching(dz_surface); nelems::Int, reverse_mode)

reverse_mode is default to false for atmosphere configurations. For land configurations, use reverse_mode = true.

faces contain reference z without any warping.

source

Mesh utilities

ClimaCore.Meshes.truncate_meshFunction
truncate_mesh(
    parent_mesh::AbstractMesh,
    trunc_domain::IntervalDomain{CT},
)

Constructs an IntervalMesh, truncating the given parent_mesh defined on a truncated trunc_domain. The truncation preserves the number of degrees of freedom covering the space from the trunc_domain's z_bottom to z_top, adjusting the stretching.

source

Interfaces

ClimaCore.Meshes.elementsFunction
Meshes.elements(mesh::AbstractMesh)

An iterator over the elements of a mesh. Elements of a mesh can be of any type.

source
ClimaCore.Meshes.boundary_face_nameFunction
Meshes.boundary_face_name(mesh::AbstractMesh, elem, face::Int)::Union{Symbol,Nothing}

The name of the boundary facing face of element elem, or nothing if it is not on the boundary.

source
ClimaCore.Meshes.opposing_faceFunction
opelem, opface, reversed = Meshes.opposing_face(mesh::AbstractMesh, elem, face::Int)

The element and face (opelem, opface) that oppose face face of element elem.

source
ClimaCore.Meshes.coordinatesFunction
Meshes.coordinates(mesh, elem, vert::Int)
Meshes.coordinates(mesh, elem, ξ::SVector)

Return the physical coordinates of a point in an element elem of mesh. The position of the point can either be a vertex number vert or the coordinates ξ in the reference element.

source
ClimaCore.Meshes.containing_elementFunction
elem = Meshes.containing_element(mesh::AbstractMesh, coord)

The element elem in mesh containing the coordinate coord. If the coordinate falls on the boundary between two or more elements, an arbitrary element is chosen.

source
ClimaCore.Meshes.reference_coordinatesFunction
ξ = Meshes.reference_coordinates(mesh::AbstractMesh, elem, coord)

An SVector of coordinates in the reference element such that

Meshes.coordinates(mesh, elem, ξ) == coord

This can be used for interpolation to a specific point.

source
ClimaCore.Meshes.face_connectivity_matrixFunction
M = Meshes.face_connectivity_matrix(mesh, elemorder = elements(mesh))

Construct a Bool-valued SparseCSCMatrix containing the face connections of mesh. Elements are indexed according to elemorder.

Note that M[i,i] == true only if two distinct faces of element i are connected.

source
ClimaCore.Meshes.vertex_connectivity_matrixFunction
M = Meshes.vertex_connectivity_matrix(mesh, elemorder = elements(mesh))

Construct a Bool-valued SparseCSCMatrix containing the vertex connections of mesh. Elements are indexed according to elemorder.

Note that M[i,i] == true only if two distinct vertices of element i are connected.

source
ClimaCore.Meshes.linearindicesFunction
Meshes.linearindices(elemorder)

Given a data structure elemorder[i] = elem that orders elements, construct the inverse map from orderindex = linearindices(elemorder) such that orderindex[elem] = i.

This will try to use the most efficient structure available.

source