API
DataLayouts
ClimaCore.DataLayouts
— ModuleClimaCore.DataLayouts
Notation:
i,j
are horizontal node indices within an elementk
is the vertical node index within an elementf
is the field index (1 if field is scalar, >1 if it is a vector field)v
is the vertical element index in a stackh
is the element stack index
Data layout is specified by the order in which they appear, e.g. IJKFVH
indexes the underlying array as [i,j,k,f,v,h]
ClimaCore.DataLayouts.DataF
— TypeDataF{S, A} <: Data0D{S}
Backing DataLayout
for 0D point data.
ClimaCore.DataLayouts.IF
— TypeIF{S, Ni, A} <: DataSlab1D{S, Ni}
Backing DataLayout
for 1D spectral element slab data.
Nodal element data (I) are contiguous for each S
datatype struct field (F) for a single element slab.
A DataSlab1D
view can be returned from other Data1D
objects by calling slab(data, idx...)
.
ClimaCore.DataLayouts.IJF
— TypeIJF{S, Nij, A} <: DataSlab2D{S, Nij}
Backing DataLayout
for 2D spectral element slab data.
Nodal element data (I,J) are contiguous for each S
datatype struct field (F) for a single element slab.
A DataSlab2D
view can be returned from other Data2D
objects by calling slab(data, idx...)
.
ClimaCore.DataLayouts.VF
— TypeVF{S, A} <: DataColumn{S}
Backing DataLayout
for 1D FV column data.
Column level data (V) are contiguous for each S
datatype struct field (F).
A DataColumn
view can be returned from other Data1DX
, Data2DX
objects by calling column(data, idx...)
.
ClimaCore.DataLayouts.IFH
— TypeIFH{S, Ni, A} <: Data1D{S, Ni}
Backing DataLayout
for 1D spectral element slabs.
Element nodal point (I) data is contiguous for each datatype S
struct field (F), for each 1D mesh element (H).
ClimaCore.DataLayouts.IJFH
— TypeIJFH{S, Nij, A} <: Data2D{S, Nij}
Backing DataLayout
for 2D spectral element slabs.
Element nodal point (I,J) data is contiguous for each datatype S
struct field (F), for each 2D mesh element slab (H).
ClimaCore.DataLayouts.VIFH
— TypeVIFH{S, Ni, A} <: Data1DX{S, Ni}
Backing DataLayout
for 1D spectral element slab + extruded 1D FV column data.
Column levels (V) are contiguous for every element nodal point (I) for each datatype S
struct field (F), for each 1D mesh element slab (H).
ClimaCore.DataLayouts.VIJFH
— TypeVIJFH{S, Nij, A} <: Data2DX{S, Nij}
Backing DataLayout
for 2D spectral element slab + extruded 1D FV column data.
Column levels (V) are contiguous for every element nodal point (I, J) for each S
datatype struct field (F), for each 2D mesh element slab (H).
Geometry
Coordinates
ClimaCore.Geometry.AbstractPoint
— TypeAbstractPoint
Represents a point in space.
The following types are supported:
XPoint(x)
YPoint(y)
ZPoint(z)
XYPoint(x, y)
XZPoint(x, z)
XYZPoint(x, y, z)
LatLongPoint(lat, long)
LatLongZPoint(lat, long, z)
Cartesian1Point(x1)
Cartesian2Point(x2)
Cartesian3Point(x3)
Cartesian12Point(x1, x2)
Cartesian13Point(x1, x3)
Cartesian123Point(x1, x2, x3)
ClimaCore.Geometry.float_type
— Functionfloat_type(T)
Return the floating point type backing T
: T
can either be an object or a type.
Points represent locations in space, specified by coordinates in a given coordinate system (Cartesian, spherical, etc), whereas vectors, on the other hand, represent displacements in space.
An analogy with time works well: times (also called instants or datetimes) are locations in time, while, durations are displacements in time.
Note 1: Latitude and longitude are specified via angles (and, therefore, trigonometric functions: cosd
, sind
, acosd
, asind
, tand
,...) in degrees, not in radians. Moreover, lat
(usually denoted by $\theta$) $\in [-90.0, 90.0]$, and long
(usually denoted by $\lambda$) $\in [-180.0, 180.0]$.
Note 2:: In a Geometry.LatLongZPoint(lat, long, z)
, z
represents the elevation above the surface of the sphere with radius R (implicitly accounted for in the geoemtry).
Note 3: There are also a set of specific Cartesian points (Cartesian1Point(x1)
, Cartesian2Point(x2)
, etc). These are occasionally useful for converting everything to a full Cartesian domain (e.g. for visualization purposes). These are distinct from XYZPoint
as ZPoint
can mean different things in different domains.
Domains
Types
ClimaCore.Domains.AbstractDomain
— TypeAbstractDomain
A domain represents a region of space.
ClimaCore.Domains.IntervalDomain
— TypeIntervalDomain(coord⁻, coord⁺; periodic=true)
IntervalDomain(coord⁻, coord⁺; boundary_names::Tuple{Symbol,Symbol})
Construct a IntervalDomain
, the closed interval is given by coord⁻
, coord⁺
coordinate arguments.
Either a periodic
or boundary_names
keyword argument is required.
ClimaCore.Domains.RectangleDomain
— TypeRectangleDomain(x1::ClosedInterval, x2::ClosedInterval;
x1boundary::Tuple{Symbol,Symbol},
x2boundary::Tuple{Symbol,Symbol},
x1periodic = false,
x2periodic = false,
)
Construct a RectangularDomain
in the horizontal. If a given x1 or x2 boundary is not periodic, then x1boundary
or x2boundary
boundary name keyword arguments must be supplied.
ClimaCore.Domains.SphereDomain
— TypeSphereDomain(radius)
A domain representing the surface of a sphere with radius radius
.
Interfaces
ClimaCore.Domains.boundary_names
— Functionboundary_names(obj::Union{AbstractDomain, AbstractMesh, AbstractTopology})
A tuple or vector of unique boundary names of a spatial domain.
Meshes
A Mesh
is a division of a domain into elements.
Mesh types
ClimaCore.Meshes.AbstractMesh
— TypeAbstractMesh{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:
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/methods are provided by AbstractMesh
:
ClimaCore.Meshes.IntervalMesh
— TypeIntervalMesh <: 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:
ClimaCore.Meshes.RectilinearMesh
— TypeRectilinearMesh <: 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
.
ClimaCore.Meshes.AbstractCubedSphere
— TypeAbstractCubedSphere <: AbstractMesh2D
This is an abstract type of cubed-sphere meshes on SphereDomain
s. 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
: aSphereDomain
ne
: number of elements across each panel
External links
ClimaCore.Meshes.EquiangularCubedSphere
— TypeEquiangularCubedSphere <: AbstractCubedSphere
An equiangular gnomonic mesh proposed by C. Ronchi, R. Iacono, P. S. Paolucci (1996). 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.
ClimaCore.Meshes.EquidistantCubedSphere
— TypeEquidistantCubedSphere <: AbstractCubedSphere
An equidistant gnomonic mesh outlined in M. Rančić, R. J. Purser, F. Mesinger (1996) and Ramachandran D Nair, Stephen J Thomas, Richard D Loft (2005). 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.
ClimaCore.Meshes.ConformalCubedSphere
— TypeConformalCubedSphere <: AbstractCubedSphere
A conformal mesh outlined in M. Rančić, R. J. Purser, F. Mesinger (1996). 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.
Local element map
ClimaCore.Meshes.LocalElementMap
— TypeLocalElementMap
An abstract type of mappings from the reference element to a physical domain.
ClimaCore.Meshes.IntrinsicMap
— TypeIntrinsicMap()
This LocalElementMap
uses the intrinsic mapping of the cubed sphere to map the reference element to the physical domain.
ClimaCore.Meshes.NormalizedBilinearMap
— TypeNormalizedBilinearMap()
The LocalElementMap
for meshes on spherical domains of O. Guba, M. A. Taylor, P. A. Ullrich, J. R. Overfelt, M. N. Levy (2014). It uses bilinear interpolation between the Cartesian coordinates of the element vertices, then normalizes the result to lie on the sphere.
Mesh stretching
ClimaCore.Meshes.Uniform
— TypeUniform()
Use uniformly-sized elements.
ClimaCore.Meshes.ExponentialStretching
— TypeExponentialStretching(H)
Apply exponential stretching to the domain when constructing elements. H
is the scale height (a typical atmospheric scale height H ≈ 7.5
km).
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.
Interfaces
ClimaCore.Meshes.domain
— Functiondomain(topology)
Returns the domain of the topology
from the underlying mesh
Meshes.domain(mesh::AbstractMesh)
The Domains.AbstractDomain
on which the mesh is defined.
ClimaCore.Meshes.elements
— FunctionMeshes.elements(mesh::AbstractMesh)
An iterator over the elements of a mesh. Elements of a mesh can be of any type.
ClimaCore.Meshes.nelements
— Functionnelements(mesh::AbstractMesh)
The number of elements in the mesh.
ClimaCore.Meshes.is_boundary_face
— FunctionMeshes.is_boundary_face(mesh::AbstractMesh, elem, face::Int)::Bool
Determine whether face face
of element elem
is on the boundary of mesh
.
elem
should be an element of elements(mesh)
.
ClimaCore.Meshes.boundary_face_name
— FunctionMeshes.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.
ClimaCore.Meshes.opposing_face
— Functionopelem, opface, reversed = Meshes.opposing_face(mesh::AbstractMesh, elem, face::Int)
The element and face (opelem
, opface
) that oppose face face
of element elem
.
ClimaCore.Meshes.coordinates
— FunctionMeshes.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.
ClimaCore.Meshes.containing_element
— Functionelem = 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.
ClimaCore.Meshes.reference_coordinates
— Functionξ = 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.
ClimaCore.Meshes.SharedVertices
— TypeMeshes.SharedVertices(mesh, elem, vert)
An iterator over (element, vertex) pairs that are shared with (elem,vert)
.
ClimaCore.Meshes.face_connectivity_matrix
— FunctionM = 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.
ClimaCore.Meshes.vertex_connectivity_matrix
— FunctionM = 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.
ClimaCore.Meshes.linearindices
— FunctionMeshes.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.
Topologies
A Topology
determines the ordering and connections between elements of a mesh.
Types
ClimaCore.Topologies.AbstractTopology
— TypeAbstractTopology
Subtypes of AbstractHorizontalTopology
define connectiveness of a mesh in the horizontal domain.
Interfaces
ClimaCore.Topologies.IntervalTopology
— TypeIntervalTopology(mesh::IntervalMesh)
A sequential topology on an Meshes.IntervalMesh
.
ClimaCore.Topologies.Topology2D
— TypeTopology2D(mesh::AbstractMesh2D, elemorder=Mesh.elements(mesh))
This is a generic non-distributed topology for 2D meshes. elemorder
is a vector or other linear ordering of the Mesh.elements(mesh)
.
ClimaCore.Topologies.spacefillingcurve
— Functionspacefillingcurve(mesh::Meshes.AbstractCubedSphere)
Generate element ordering, elemorder
, based on a space filling curve for a CubedSphere
mesh.
spacefillingcurve(mesh::Meshes.RectilinearMesh)
Generate element ordering, elemorder
, based on a space filling curve for a Rectilinear
mesh.
Interfaces
Missing docstring for Topologies.domain
. Check Documenter's build log for details.
ClimaCore.Topologies.mesh
— Functionmesh(topology)
Returns the mesh underlying the topology
ClimaCore.Topologies.nlocalelems
— Functionnlocalelems(topology)
The number of local elements in topology
.
ClimaCore.Topologies.vertex_coordinates
— Function(c1,c2,c3,c4) = vertex_coordinates(topology, elem)
The coordinates of the 4 vertices of element elem
.
ClimaCore.Topologies.opposing_face
— Function(opelem, opface, reversed) = opposing_face(topology, elem, face)
The opposing face of face number face
of element elem
in topology
.
opelem
is the opposing element number, 0 for a boundary, negative for a ghost elementopface
is the opposite face number, or boundary face number if a boundaryreversed
indicates whether the opposing face has the opposite orientation.
ClimaCore.Topologies.interior_faces
— Functioninterior_faces(topology::AbstractTopology)
An iterator over the interior faces of topology
. Each element of the iterator is a 5-tuple the form
(elem1, face1, elem2, face2, reversed)
where elemX, faceX
are the element and face numbers, and reversed
indicates whether they have opposing orientations.
ClimaCore.Topologies.boundary_tags
— Functionboundary_tags(topology)
A Tuple
or NamedTuple
of the boundary tags of the topology. A boundary tag is an integer that uniquely identifies a boundary.
ClimaCore.Topologies.boundary_tag
— Functionboundary_tag(topology, name::Symbol)
The boundary tag of the topology for boundary name name
. A boundary tag is an integer that uniquely identifies a boundary.
ClimaCore.Topologies.boundary_faces
— Functionboundary_faces(topology, boundarytag)
An iterator over the faces of topology
which face the boundary with tag boundarytag
. Each element of the iterator is an (elem, face)
pair.
Missing docstring for Topologies.vertices
. Check Documenter's build log for details.
ClimaCore.Topologies.local_neighboring_elements
— Functionlocal_neighboring_elements(topology::AbstractTopology, lidx::Integer)
An iterator of the local element indices (lidx) of the local elements which are neighbors of the local element lidx
in topology
(excluding lidx
itself).
ClimaCore.Topologies.ghost_neighboring_elements
— Functionghost_neighboring_elements(topology::AbstractTopology, ridx::Integer)
An iterator of the receive buffer indices (ridx) of the ghost elements which are neighbors of the local element lidx
in topology
.
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:
ClimaCore.Spaces
— ModuleFinite 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.
Quadratures
ClimaCore.Spaces.Quadratures.QuadratureStyle
— TypeClimaCore.Spaces.Quadratures.GLL
— TypeGLL{Nq}()
Gauss-Legendre-Lobatto quadrature using Nq
quadrature points.
ClimaCore.Spaces.Quadratures.GL
— TypeGL{Nq}()
Gauss-Legendre quadrature using Nq
quadrature points.
ClimaCore.Spaces.Quadratures.Uniform
— TypeUniform{Nq}()
Uniformly-spaced quadrature.
ClimaCore.Spaces.Quadratures.degrees_of_freedom
— Functiondegrees_of_freedom(QuadratureStyle) -> Int
Returns the degreesoffreedom of the QuadratureStyle
concrete type
ClimaCore.Spaces.Quadratures.polynomial_degree
— Functionpolynomial_degree(QuadratureStyle) -> Int
Returns the polynomial degree of the QuadratureStyle
concrete type
ClimaCore.Spaces.Quadratures.quadrature_points
— Functionpoints, weights = quadrature_points(::Type{FT}, quadrature_style)
The points and weights of the quadrature rule in floating point type FT
.
ClimaCore.Spaces.Quadratures.barycentric_weights
— Functionbarycentric_weights(x::SVector{Nq}) where {Nq}
The barycentric weights associated with the array of point locations x
:
\[w_j = \frac{1}{\prod_{k \ne j} (x_i - x_j)}\]
See Jean-Paul Berrut, Lloyd N Trefethen (2004), equation 3.2.
ClimaCore.Spaces.Quadratures.interpolation_matrix
— Functioninterpolation_matrix(x::SVector, r::SVector{Nq})
The matrix which interpolates the Lagrange polynomial of degree Nq-1
through the points r
, to points x
. The matrix coefficients are computed using the Barycentric formula of Jean-Paul Berrut, Lloyd N Trefethen (2004), section 4:
\[I_{ij} = \begin{cases} 1 & \text{if } x_i = r_j, \\ 0 & \text{if } x_i = r_k \text{ for } k \ne j, \\ \frac{\displaystyle \frac{w_j}{x_i - r_j}}{\displaystyle \sum_k \frac{w_k}{x_i - r_k}} & \text{otherwise,} \end{cases}\]
where $w_j$ are the barycentric weights, see barycentric_weights
.
ClimaCore.Spaces.Quadratures.differentiation_matrix
— Functiondifferentiation_matrix(r::SVector{Nq, T}) where {Nq, T}
The spectral differentiation matrix for the Lagrange polynomial of degree Nq-1
interpolating at points r
.
The matrix coefficients are computed using the Jean-Paul Berrut, Lloyd N Trefethen (2004), section 9.3:
\[D_{ij} = \begin{cases} \displaystyle \frac{w_j}{w_i (x_i - x_j)} &\text{ if } i \ne j \\ -\sum_{k \ne j} D_{kj} &\text{ if } i = j \end{cases}\]
where $w_j$ are the barycentric weights, see barycentric_weights
.
differentiation_matrix(FT, quadstyle::QuadratureStyle)
The spectral differentiation matrix at the quadrature points of quadstyle
, using floating point types FT
.
ClimaCore.Spaces.Quadratures.orthonormal_poly
— FunctionV = orthonormal_poly(points, quad)
V_{ij}
contains the j-1
th Legendre polynomial evaluated at points[i]
. i.e. it is the mapping from the modal to the nodal representation.
Internals
ClimaCore.Spaces.dss_transform
— Functiondss_transform(arg, local_geometry, weight, I...)
Transfrom arg[I...]
to a basis for direct stiffness summation (DSS). Transformations only apply to vector quantities.
local_geometry[I...]
is the relevantLocalGeometry
object. If it isnothing
, then no transformation is performedweight[I...]
is the relevant DSS weights. Ifweight
isnothing
, then the result is simply summation.
See Spaces.weighted_dss!
.
ClimaCore.Spaces.dss_untransform
— Functiondss_untransform(T, targ, local_geometry, I...)
Transform targ[I...]
back to a value of type T
after performing direct stiffness summation (DSS).
See Spaces.weighted_dss!
.
ClimaCore.Spaces.dss_interior_faces!
— Functiondss_interior_faces!(topology, data [, local_geometry_data=nothing, local_weights=nothing])
Perform DSS on the local interior faces of the topology.
If local_geometry
is nothing
, no transformations are applied to vectors. If local_weights
is nothing
, no weighting is applied (i.e. colocated values are summed).
Part of Spaces.weighted_dss!
.
ClimaCore.Spaces.dss_local_vertices!
— Functiondss_local_vertices!(topology, data[, local_geometry_data=nothing, local_weights=nothing])
Perform DSS on the local vertices of the topology.
Part of Spaces.weighted_dss!
.
ClimaCore.Spaces.dss_ghost_faces!
— Functiondss_ghost_faces!(topology, data, ghost_data,
local_geometry_data=nothing, ghost_geometry_data=nothing,
local_weights=nothing ghost_weights=nothing;
update_ghost=false)
Perform DSS on the ghost faces of the topology. ghost_data
should contain the ghost element data.
Part of Spaces.weighted_dss!
.
ClimaCore.Spaces.dss_ghost_vertices!
— Functiondss_ghost_vertices!(topology, data, ghost_data,
local_geometry_data=nothing, ghost_geometry_data=nothing,
local_weights=nothing ghost_weights=nothing;
update_ghost=false)
Perform DSS on the ghost faces of the topology. ghost_data
should contain the ghost element data.
Part of Spaces.weighted_dss!
.
RecursiveApply
ClimaCore.RecursiveApply
— ModuleRecursiveApply
This module contains operators to recurse over nested Tuple
s or NamedTuple
s.
To extend to another type T
, define RecursiveApply.rmap(fn, args::T...)
ClimaCore.RecursiveApply.tuplemap
— Functiontuplemap(fn::Function, tup)
A map
impl for mapping function fn
a tuple argument tup
Currently just calls Base.map
behind the scenes but is left stubbed out for potential specialization in the future.
tuplemap(fn::Function, tup1, tup2)
A map
impl for mapping function fn
over tup1
, tup2
tuple arguments.
Currently just calls Base.map
behind the scenes but is left stubbed out for potential specialization in the future.
Fields
ClimaCore.Fields.Field
— TypeField(values, space)
A set of values
defined at each point of a space
.
ClimaCore.Fields.coordinate_field
— Functioncoordinate_field(space::AbstractSpace)
Construct a Field
of the coordinates of the space.
ClimaCore.Fields.local_geometry_field
— Functionlocal_geometry_field(space::AbstractSpace)
Construct a Field
of the LocalGeometry
of the space.
Base.zeros
— Methodzeros(space::AbstractSpace)
Construct a field on space
that is zero everywhere.
Base.ones
— Methodones(space::AbstractSpace)
Construct a field on space
that is one everywhere.
Base.sum
— Methodsum([f=identity,]v::Field)
Approximate integration of v
or f.(v)
over the domain. In an AbstractSpectralElementSpace
, an integral over the entire space is computed by summation over the elements of the integrand multiplied by the Jacobian determinants and the quadrature weights at each node within an element. Hence, sum
is computed by summation of the field values multiplied by the Jacobian determinants and quadrature weights:
\[\sum_i f(v_i) W_i J_i \approx \int_\Omega f(v) \, d \Omega\]
where $v_i$ is the value at each node, and $f$ is the identity function if not specified.
Statistics.mean
— Methodmean([f=identity, ]v::Field)
The mean value of field
or f.(field)
over the domain, weighted by area. Similar to sum
, in an AbstractSpectralElementSpace
, this is computed by summation of the field values multiplied by the Jacobian determinants and quadrature weights:
\[\frac{\sum_i f(v_i) W_i J_i}{\sum_i W_i J_i} \approx \frac{\int_\Omega f(v) \, d \Omega}{\int_\Omega \, d \Omega}\]
where $v_i$ is the Field value at each node, and $f$ is the identity function if not specified.
LinearAlgebra.norm
— Methodnorm(v::Field, p=2; normalize=true)
The approximate $L^p$ norm of v
, where $L^p$ represents the space of measurable functions for which the p-th power of the absolute value is Lebesgue integrable, that is:
\[\| v \|_p = \left( \int_\Omega |v|^p d \Omega \right)^{1/p}\]
where $|v|$ is defined to be the absolute value if $v$ is a scalar-valued Field, or the 2-norm if it is a vector-valued Field or composite Field (see LinearAlgebra.norm). Similar to sum
and mean
, in an AbstractSpectralElementSpace
, this is computed by summation of the field values multiplied by the Jacobian determinants and quadrature weights. If normalize=true
(the default), then internally the discrete norm is divided by the sum of the Jacobian determinants and quadrature weights:
\[\left(\frac{\sum_i |v_i|^p W_i J_i}{\sum_i W_i J_i}\right)^{1/p} \approx \left(\frac{\int_\Omega |v|^p \, d \Omega}{\int_\Omega \, d \Omega}\right)^{1/p}\]
If p=Inf
, then the norm is the maximum of the absolute values
\[\max_i |v_i| \approx \sup_{\Omega} |v|\]
Consequently all norms should have the same units for all $p$ (being the same as calling norm
on a single value).
If normalize=false
, then the denominator term is omitted, and so the result will be the norm as described above multiplied by the length/area/volume of $\Omega$ to the power of $1/p$.
ClimaCore.Fields.set!
— Functionset!(f::Function, field::Field, args = ())
Apply function f
to populate values in field field
. f
must have a function signature with signature f(::LocalGeometry[, args...])
. Additional arguments may be passed to f
with args
.
ClimaCore.Fields.ColumnIndex
— TypeColumnIndex(ij,h)
An index into a column of a field. This can be used as an argument to getindex
of a Field
, to return a field on that column.
Example
colidx = ColumnIndex((1,1),1)
field[colidx]
ClimaCore.Fields.bycolumn
— FunctionFields.bycolumn(fn, space)
Call fn(colidx)
to every ColumnIndex
colidx
of space
. This can be used to apply multiple column-wise operations in a single pass, making use of multiple threads.
Example
∇ = GradientF2C()
div = DivergenceC2F()
bycolumn(axes(f)) do colidx
@. ∇f[colidx] = ∇(f[colidx])
@. df[colidx] = div(∇f[colidx])
end
Limiters
Interfaces
ClimaCore.Limiters.QuasiMonotoneLimiter
— TypeQuasiMonotoneLimiter
This limiter is inspired by the one presented in Guba et al Oksana Guba, Mark Taylor, Amik St-Cyr (2014). In the reference paper, it is denoted by OP1, and is outlined in eqs. (37)-(40). Quasimonotone here is meant to be monotone with respect to the spectral element nodal values. This limiter involves solving a constrained optimization problem (a weighted least square problem up to a fixed tolerance) that is completely local to each element.
As in HOMME, the implementation idea here is the following: we need to find a grid field which is closest to the initial field (in terms of weighted sum), but satisfies the min/max constraints. So, first we find values that do not satisfy constraints and bring these values to a closest constraint. This way we introduce some change in the tracer mass, which we then redistribute so that the l2 error is smallest. This redistribution might violate constraints; thus, we do a few iterations (until abs(Δtracer_mass) <= rtol * tracer_mass
).
ρq
: tracer density Field, whereq
denotes tracer concentration per unit mass. This can be a scalar field, or a struct-valued field.ρ
: fluid density Field (scalar).
Constructor
limiter = QuasiMonotoneLimiter(ρq::Field; rtol = eps(eltype(parent(ρq))))
Creates a limiter instance for the field ρq
with relative tolerance rtol
.
Usage
Call compute_bounds!
on the input fields:
compute_bounds!(limiter, ρq, ρ)
Then call apply_limiter!
on the output fields:
apply_limiter!(ρq, ρ, limiter)
ClimaCore.Limiters.compute_bounds!
— Functioncompute_bounds!(limiter::QuasiMonotoneLimiter, ρq::Field, ρ::Field)
Compute the desired bounds for the tracer concentration per unit mass q
, based on the tracer density, ρq
, and density, ρ
, fields.
This is computed by
compute_element_bounds!
- starts the ghost exchange (if distributed)
compute_neighbor_bounds_local!
- completes the ghost exchange (if distributed)
compute_neighbor_bounds_ghost!
(if distributed)
ClimaCore.Limiters.apply_limiter!
— Functionapply_limiter!(ρq, ρ, limiter::QuasiMonotoneLimiter)
Apply the limiter on the tracer density ρq
, using the computed desired bounds on the concentration q
and density ρ
as an optimal weight. This iterates over each element, calling apply_limit_slab!
. If the limiter fails to converge for any element, a warning is issued.
Internals
ClimaCore.Limiters.compute_element_bounds!
— Functioncompute_element_bounds!(limiter::QuasiMonotoneLimiter, ρq, ρ)
Given two fields ρq
and ρ
, computes the min and max of q
in each element, storing it in limiter.q_bounds
.
Part of compute_bounds!
.
ClimaCore.Limiters.compute_neighbor_bounds_local!
— Functioncompute_neighbor_bounds_local!(limiter::QuasiMonotoneLimiter, topology)
Update the field limiter.q_bounds_nbr
based on limiter.q_bounds
in the local neighbors.
Part of compute_bounds!
.
ClimaCore.Limiters.compute_neighbor_bounds_ghost!
— Functioncompute_neighbor_bounds_ghost!(limiter::QuasiMonotoneLimiter, topology)
Update the field limiter.q_bounds_nbr
based on limiter.q_bounds
in the ghost neighbors. This should be called after the ghost exchange has completed.
Part of compute_bounds!
.
ClimaCore.Limiters.apply_limit_slab!
— Functionapply_limit_slab!(slab_ρq, slab_ρ, slab_WJ, slab_q_bounds, rtol)
Apply the computed bounds of the tracer concentration (slab_q_bounds
) in the limiter to slab_ρq
, given the total mass slab_ρ
, metric terms slab_WJ
, and relative tolerance rtol
. Return whether the tolerance condition could be satisfied.
InputOutput
Writers
ClimaCore.InputOutput.HDF5Writer
— TypeHDF5Writer(filename::AbstractString[, context::ClimaComms.AbstractCommsContext])
An AbstractWriter
for writing to HDF5-formatted files using the ClimaCore storage conventions. An internal cache is used to avoid writing duplicate domains, meshes, topologies and spaces to the file. Use HDF5Reader
to load the data from the file.
The optional context
can be used for writing distributed fields: in this case, the MPICommsContext
used passed as an argument: this must match the context used for distributing the Field
.
The default Julia HDF5 binaries are not built with MPI support. To use the distributed functionality, you will need to configure HDF5.jl with an MPI-enabled HDF5 library, see the HDF5.jl documentation.
Interface
Usage
writer = InputOutput.HDF5Writer(filename)
InputOutput.write!(writer, Y, "Y")
close(writer)
ClimaCore.InputOutput.write!
— Functionwrite!(writer::AbstractWriter, obj[, preferredname])
Write the object obj
using writer
. An optional preferredname
can be provided, otherwise defaultname
will be used to generate a name. The name of the object will be returned.
A cache of domains, meshes, topologies and spaces is kept: if one of these objects has already been written, then the file will not be modified: instead the name under which the object was first written will be returned. Note that Field
s and FieldVector
s are not cached, and so can be written multiple times.
write!(filename, fvpair)
Write 'FieldVector' data, specified by a pair fieldvectorname => fieldvector
, to HDF5 file 'filename'.
Readers
ClimaCore.InputOutput.HDF5Reader
— TypeHDF5Reader(filename::AbstractString[, context::ClimaComms.AbstractCommsContext])
An AbstractReader
for reading from HDF5 files created by HDF5Writer
. The reader object contains an internal cache of domains, meshes, topologies and spaces that are read so that duplicate objects are not created.
The optional context
can be used for reading distributed fields: in this case, the MPICommsContext
used passed as an argument: resulting Field
s will be distributed using this context. As with HDF5Writer
, this requires a HDF5 library with MPI support.
Interface
Usage
reader = InputOutput.HDF5Reader(filename)
Y = read_field(reader, "Y")
Y.c |> propertynames
Y.f |> propertynames
ρ_field = read_field(reader, "Y.c.ρ")
w_field = read_field(reader, "Y.f.w")
close(reader)
To explore the contents of the reader
, use either
julia> reader |> propertynames
e.g, to explore the components of the space
,
julia> reader.space_cache
Dict{Any, Any} with 3 entries:
"center_extruded_finite_difference_space" => CenterExtrudedFiniteDifferenceSpace:…
"horizontal_space" => SpectralElementSpace2D:…
"face_extruded_finite_difference_space" => FaceExtrudedFiniteDifferenceSpace:…
Once "unpacked" as shown above, ClimaCorePlots
or ClimaCoreMakie
can be used to visualise fields. ClimaCoreTempestRemap
supports interpolation onto user-specified grids if necessary.
ClimaCore.InputOutput.read_domain
— Functionread_domain(reader::AbstractReader, name)
Reads a domain named name
from reader
. Domain objects are cached in the reader to avoid creating duplicate objects.
ClimaCore.InputOutput.read_mesh
— Functionread_mesh(reader::AbstractReader, name)
Reads a mesh named name
from reader
, or from the reader cache if it has already been read.
ClimaCore.InputOutput.read_topology
— Functionread_topology(reader::AbstractReader, name)
Reads a topology named name
from reader
, or from the reader cache if it has already been read.
ClimaCore.InputOutput.read_space
— Functionread_space(reader::AbstractReader, name)
Reads a space named name
from reader
, or from the reader cache if it has already been read.
ClimaCore.InputOutput.read_field
— Functionread_field(reader, name)
Reads a Field
or FieldVector
named name
from reader
. Fields are not cached, so that reading the same field multiple times will create multiple distinct objects.