Meshes

AbstractMesh{dim}

Abstract supertype of meshes, which describe how a domain is discretized into elements. Subtypes: IntervalMesh, RectilinearMesh, and the AbstractCubedSphere meshes.

A mesh is lightweight and exists on all MPI ranks; a mesh stored in a file, for example, holds only the filename.

Notes

Face and vertex numbering. In 1D, faces and vertices coincide, and both are numbered [1, 2]. In 2D, a face is a line segment between two vertices, and both are numbered [1, 2, 3, 4] in counter-clockwise order.

 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 defines the following methods:

• domain(mesh)
• elements(mesh)
• is_boundary_face(mesh, elem, face)
• boundary_face_name(mesh, elem, face)
• opposing_face(mesh, elem, face)
• coordinates(mesh, elem, vert)
• containing_element (optional)

The following types and methods are provided for every AbstractMesh:

• SharedVertices(mesh, elem, vert)
• face_connectivity_matrix(mesh[,elemorder])
• vertex_connectivity_matrix(mesh[,elemorder])

AbstractMesh1D

Alias for AbstractMesh{1} (see AbstractMesh): the supertype of one-dimensional meshes, i.e. IntervalMesh, whose elements are intervals.

AbstractMesh2D

Alias for AbstractMesh{2} (see AbstractMesh): the supertype of two-dimensional meshes, i.e. RectilinearMesh and the AbstractCubedSphere meshes, whose elements are quadrilaterals.

IntervalMesh <: AbstractMesh

One-dimensional mesh on an IntervalDomain. Elements are numbered 1:nelems.

Fields

• stretch: The stretching rule used to place the faces.
• domain: The IntervalDomain.
• faces: Vector of face coordinates, of length nelems + 1.
• meta: Stretching-specific metadata (e.g. the solved stretching parameter), or nothing.
• reverse_mode: Whether the smallest element is at the top of the domain.

Constructor

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

Construct a 1D mesh with face locations at faces.

IntervalMesh(domain::IntervalDomain, stretching = Uniform(); nelems, reverse_mode = false)

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

• Uniform()
• ExponentialStretching(H)
• GeneralizedExponentialStretching(dz_bottom, dz_top)
• HyperbolicTangentStretching(dz_surface)

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

RectilinearMesh <: AbstractMesh2D

Rectilinear mesh of a RectangleDomain: the product of two IntervalMeshes. Elements are indexed by CartesianIndex{2}.

Constructor

RectilinearMesh(domain::RectangleDomain, n1, n2)

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

RectilinearMesh(intervalmesh1::IntervalMesh, intervalmesh2::IntervalMesh)

Construct the product mesh of intervalmesh1 and intervalmesh2.

AbstractCubedSphere <: AbstractMesh2D

Abstract supertype of cubed-sphere meshes on SphereDomains. Subtypes: EquiangularCubedSphere, EquidistantCubedSphere, and ConformalCubedSphere. 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 panel ordering of the S2 Geometry library, with 1-based instead of 0-based numbering.

Elements are indexed by a CartesianIndex{3}, whose components are:

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

Fields

Subtypes have the following fields:

• domain: The SphereDomain.
• ne: Number of elements across each panel.
• localelementmap: The LocalElementMap from the reference element to the sphere.

Notes

External references:

• S2Geometry library.
• MIT GCM exch2.

EquiangularCubedSphere <: AbstractCubedSphere

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

Constructor

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

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

EquidistantCubedSphere <: AbstractCubedSphere

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

Constructor

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

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

ConformalCubedSphere <: AbstractCubedSphere

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

Constructor

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

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

LocalElementMap

Abstract supertype of mappings from the reference element to the physical domain. Subtypes: IntrinsicMap and NormalizedBilinearMap.

IntrinsicMap()

LocalElementMap that uses the intrinsic mapping of the cubed sphere to map the reference element to the physical domain.

NormalizedBilinearMap()

LocalElementMap for meshes on spherical domains of [50]. It uses bilinear interpolation between the Cartesian coordinates of the element vertices, then normalizes the result to lie on the sphere.

Uniform()

Use uniformly-sized elements.

ExponentialStretching(H::FT)

Apply exponential stretching to the domain when constructing elements. H is the scale height [m]; a typical atmospheric scale height is H ≈ 7.5e3.

For an interval $[z_0, z_1]$, the elements are uniformly spaced in $ζ$, 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).

Construct a stretched mesh via

IntervalMesh(interval_domain, ExponentialStretching(H); nelems, reverse_mode = false)

The resulting faces are reference heights without terrain warping.

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 the target element spacings at the bottom and at the top of the vertical column domain [m]. In typical atmosphere configurations, dz_bottom is the smallest spacing and dz_top the largest; in typical land configurations, dz_bottom is the largest and dz_top the smallest. For land configurations, use reverse_mode = true (default false).

Construct a generalized stretched mesh via

IntervalMesh(
    interval_domain,
    GeneralizedExponentialStretching(dz_bottom, dz_top);
    nelems,
    FT_solve = Float64,
    tol = 1e-3,
    reverse_mode = false,
)

where FT_solve is the float type and tol the residual tolerance of the root solve for the stretching parameters. nelems must be at least 2. The resulting faces are reference heights without terrain warping.

HyperbolicTangentStretching(dz_surface::FT)

Apply hyperbolic tangent stretching to the domain when constructing elements. dz_surface is the target element spacing at the surface [m]: in typical atmosphere configurations, the spacing at the bottom of the vertical column domain; in typical land configurations, the spacing at the top.

For an interval $[z_0, z_1]$, the elements are uniformly spaced in $ζ$, 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.

Construct a stretched mesh via

IntervalMesh(
    interval_domain,
    HyperbolicTangentStretching(dz_surface);
    nelems,
    FT_solve = Float64,
    tol = nothing,
    reverse_mode = false,
)

where FT_solve is the float type and tol the residual tolerance of the root solve for $γ$ (default 1e-6 * dz_surface). nelems must be at least 2. For land configurations, use reverse_mode = true. The resulting faces are reference heights without terrain warping.

truncate_mesh(parent_mesh::IntervalMesh, trunc_domain::IntervalDomain)

Construct an IntervalMesh on trunc_domain by truncating parent_mesh. The result has as many elements as parent_mesh has faces below the top of trunc_domain, and uses a GeneralizedExponentialStretching with the bottom and top spacings of the truncated parent faces.

Meshes.domain(mesh::AbstractMesh)

Return the domain (a subtype of Domains.AbstractDomain) on which mesh is defined.

Meshes.elements(mesh::AbstractMesh)

Return an iterator over the elements of mesh. Elements can be of any type.

nelements(mesh::AbstractMesh)

Return the number of elements in mesh.

n_elements_per_panel_direction(::AbstractCubedSphere)

Return the number of elements per panel along each direction.

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.

Elements of a 2D mesh are quadrilaterals, so vert must be in 1:4; elements of a 1D mesh have two vertices, so vert must be in 1:2.

Meshes.containing_element(mesh::AbstractMesh, coord)

Return the element of mesh containing the coordinate coord. If the coordinate falls on the boundary between two or more elements, one of them is chosen arbitrarily.

Meshes.reference_coordinates(mesh::AbstractMesh, elem, coord)

Return an SVector ξ of coordinates in the reference element such that Meshes.coordinates(mesh, elem, ξ) == coord. Used for interpolation to a specific point.

Meshes.linearindices(elemorder)

Given a data structure elemorder with elemorder[i] = elem that orders elements, return the inverse map orderindex such that orderindex[elem] = i.

The result is a LinearIndices for CartesianIndices, a dense Int array for a vector of CartesianIndex, and a Dict otherwise.

Meshes.element_horizontal_length_scale(mesh::AbstractMesh)

Return the approximate length scale of the elements of mesh, in the units of the domain coordinates.