CFDomains

Parent type of SpectralDomain and FDDomain

FDDomain <: AbstractDomain`

struct HVLayout{rank} end
layout = HVLayout(rank)

Singleton type describing a multi-layer data layout where horizontal layers are contiguous in memory. rank is the number of horizontal indices (1 or 2).

filter = HyperDiffusion(fieldtype::Symbol, domain, niter::Int, nu)                  # user-friendly
filter = HyperDiffusion{fieldtype, D, F, X}(domain::D, niter::Int, nu::F, extra::X) # internal

Return a filter that applies Laplacian diffusion iterated niter times with hyperdiffusive coefficient nu ≥ 0 on fields of type fieldtype. Supported field types are :scalar, :vector, :vector_curl, :vector_div.

The filter is to be used as:

# out-of-place, allocates `scratch`, allocates and returns `filtered_data`
filtered_data = filter(void, in_data, void)

# mutating, non-allocating
out_data = filter(out_data, in_data, scratch)

# in-place, non-allocating
inout_data = filter(inout_data, inout_data, scratch)

filter(out_data, in_data, scratch) returns hyperdiffusion!(filter.domain, filter, out_data, in_data, scratch). If scratch::Void, then

scratch = scratch_space(filter, in_data)

The user-friendly HyperDiffusion may be specialized for specific types of domain. The specialized method should call the internal constructor which takes an extra argument extra, stored as filter.extra for later use by hyperdiffusion!. The default user-friendly constructor sets extra=nothing.

Similarly, scratch_space(filter, in_data) calls scratch_hyperdiff(f.domain, Val(fieldtype), in_data) which by default returns nothing. scratch_hyperdiff may be specialized for specific domain types and fieldtype.

abstract type MassCoordinate{N} <: VerticalCoordinate{N} end

Parent type for a mass-based vertical coordinate. Children types should specialize mass_level. See also VerticalCoordinate and mass_coordinate.

abstract type PressureCoordinate{N} <: VerticalCoordinate{N} end

Parent type for a pressure-based vertical coordinate. Children types should specialize pressure_level and mass_level. See also VerticalCoordinate.

multi_layer_domain = Shell(nz::Int, layer::AbstractDomain, layout)

Return a multi-layer domain made of nz layers with data layout specified by layout. Unless you know what you are doing, it is recommended to use rather:

multi_layer_domain = shell(nz::Int, layer::AbstractDomain)

which gets the data layout from data_layout(layer). Otherwise, multi_layer_domain may be non-optimal or non-usable.

sigma = SigmaCoordinate(N, ptop) <: PressureCoordinate{N}

Pressure based sigma-coordinate for N levels with top pressure ptop. Pressure levels are linear in vertical coordinate k : k/N = (ps-p)/(ps-ptop) where k ranges from 0 (ground) to N (model top).

SpectralDomain <: AbstractDomain

parent type of SpectralSphere

Parent type for spherical domains using spherical harmonics.

struct VHLayout{rank} end
layout = VHLayout(rank)

Singleton type describing a multi-layer data layout where vertical columns are contiguous. rank is the number of horizontal indices (1 or 2).

abstract type VerticalCoordinate{N} end

Parent type for generalized vertical coordinates ranging from 0 to N. See also PressureCoordinate.

field = allocate_field(kind::Symbol, domain::AbstractDomain, precision::Type)

Allocate a field of the given kind and precision over the given domain. Typical values for precision are Float32, Float64 or ForwardDiff.Dual. Depending on the domain, valid values for kind may include :scalar, :vector, :scalar_spec, scalar_spat (spectral/spatial representation of a scalar field), :vector_spec, vector_spat (spectral/spatial representation of a vector field).

Internally, allocate_field(kind::Symbol, domain, F) returns allocate_field(Val(kind), domain, F). To specialize allocate_field for MyDomain <: AbstractDomain, one must provide methods for:

allocate_field(::Val{kind}, domain::MyDomain, F)

where symbol kind::Symbol is one of the valid field kinds for that domain.

fields = allocate_fields(kinds::Tuple, domain::AbstractDomain, F::Type)
fields = allocate_fields(kinds::NamedTuple, domain::AbstractDomain, F::Type)

Allocate a (named) tuple of fields according to the provided kinds. For instance:

fields = allocate_fields((:vector, :scalar), domain, F)
fields = allocate_fields((a=:vector, b=:scalar), domain, F)

are equivalent to, respectively:

fields = (allocate_field(:vector, domain, F), allocate_field(:scalar, domain, F))
fields = (a=allocate_field(:vector, domain, F), b=allocate_field(:scalar, domain, F))

layout = data_layout(shell::Shell)

Return layout describing the data layout of multi-layer domain shell.

layout = data_layout(domain::Domain)

Return layout describing the preferred data layout for a shell made of layers of type Domain.

Typical values for layout are the singletons HVLayout() (layers are contiguous in memory) and VHLayout (columns are contiguous in memory).

new_data = hyperdiff_shell!(layer, layout, filter, out_data, in_data, scratch)

Apply hyperdiffusive filtertoindata. The result is written intooutdataand returned. Ifdomain::SpectralDomain, the filter applies to spectral coefficients. Ifout_data::Void, it is adequately allocated. SeeHyperDiffusion.

new_data = hyperdiffusion!(domain, filter, out_data, in_data, scratch)

Apply hyperdiffusive filtertoindata. The result is written intooutdataand returned. Ifdomain::SpectralDomain, the filter applies to spectral coefficients. Ifout_data::Void, it is adequately allocated. SeeHyperDiffusion.

Estimates the largest eigenvalue -lambda=dx^-2 of the scalar Laplace operator and returns dx which is a (non-dimensional) length on the unit sphere characterizing the mesh resolution. By design, the Courant number for the wave equation with unit wave speed solved with time step dt is dt/dx.

mcoord = mass_coordinate(pcoord::PressureCoordinate, metric_cov=1)

Return the mass-based coordinate deduced from pcoord and the covariant metric factor metric_cov. This object mcoord can then be used with:

m = mass_level(k, masstot, vcoord)

With metric_cov==1, masstot should be in Pa ( kg/m²⋅(m/s²) ). With metric_cov the covariant metric factor, masstot should be in kg⋅(m/s²). See also mass_level.

mcoord = mass_coordinate(mcoord::MassCoordinate)

Return mcoord itself, unchanged. Interim function for backwards compatibility.

m = mass_level(k, masstot, mcoord::MassCoordinate{N})

Return mass m in level k/2 as prescribed by vertical coordinate mcoord and total mass masstot.

So-called full levels correspond to odd values k=1,2...2N-1 while interfaces between full levels (so-called half-levels) correspond to even values k=0,2...2N

masstot may be:

• per unit surface, with unit Pa
• per unit non-dimensional surface (e.g. on the unit sphere), with unit kg.(m/s²))

Which convention is appropriate depends on the metric_factor provided when constructing mcoord. See also mass_coordinate.

periodize!(data, box::AbstractBox, mgr)

Enforce horizontally-periodic boundary conditions on array data representing grid point values in box. data may also be a collection, in which case periodize! is applied to each element of the collection. Call periodize! on data obtained by computations involving horizontal averaging/differencing.

p = pressure_level(k, ps, vcoord::PressureCoordinate{N})

Returns pressure p corresponding to level k/2 as prescribed by vertical coordinate vcoord and surface pressure ps.

So-called full levels correspond to odd values k=1,2...2N-1 while interfaces between full levels (so-called half-levels) correspond to even values k=0,2...2N

scratch = scratch_hyperdiff(domain, Val(fieldtype), field)

Return scratch space used to apply hyperdiffusion on domain for a field of a certain fieldtype. See HyperDiffusion and hyperdiffusion!.

scratch = scratch_space(filter::AbstractFilter, field, [scratch])

If scratch is omitted or ::Void, return scratch space for applying filter to field, allocated by scratch_hyperdiff. Otherwise just return scratch. See also HyperDiffusion.

multi_layer_domain = shell(nz::Int, layer::AbstractDomain)

Return a multi-layer domain made of nz layers with data layout specified by data_layout(layer).

interpolated = interpolate!(mgr::LoopManager, domain::Shell, field, coord, refs, increasing)

Interpolate field defined on 3D domain to reference values refs of coord. If coord increases with level number, increasing==true and vice-versa. refs must be sorted according to increasing. out may be void. Example:

T_ref = interpolate!(mgr, void, domain, [850., 500.], temperature, pressure, false)

With these arguments, interpolate! is not meant to be specialized for a specific domain type. Instead, specialize: interpolated = interpolate!(mgr, out, layout::Layout, refs, field, coord, increasing::Val)

g = Fix(f, args)

Return callable g such that g(x,...) calls f by prepending args... before x...:

g(x...) == f(args..., x...)

This is similar to Base.Fix1, with several arguments.

trisk = TRiSK(vsphere, edge, Val(N))
U_perp[edge]    = trisk(U)        # linear, single-layer, U::AbstractVector
U_perp[k, edge] = trisk(U, k)     # linear, multi-layer, U::AbstractMatrix
qU[edge]        = trisk(U, q)     # nonlinear, single-layer
qU[k, edge]     = trisk(U, q, k)  # nonlinear, multi-layer

Compute TRiSK operator U⟂ or q×U at edge::Int of vsphere::VoronoiSphere.

U is a contravariant vector field known at edges. U_perp is a covariant vector field known at edges.

N=sphere.trisk_deg[edge] is the number of edges involved in the TRiSK stencil and must be provided as a compile-time constant for performance. This may be done via the macro @unroll from ManagedLoops.

@inbounds propagates into TRiSK and trisk.

avg = average_ie(vsphere, edge)
qe[edge] = avg(qi)         # single-layer, qi::AbstractVector
qe[k, edge] = avg(qi, k)   # multi-layer, qi::AbstractMatrix

Interpolate scalar field at edge::Int of vsphere::VoronoiSphere by a centered average (second-order accurate).

qi is a scalar field known at primal cells. qe is a scalar field known at edges.

@inbounds propagates into average_ie and avg.

avg = average_iv(vsphere, dual_cell)
qv[dual_cell] = avg(qi)         # single-layer, qi::AbstractVector
qv[k, dual_cell] = avg(qi, k)   # multi-layer, qi::AbstractMatrix

Interpolate scalar field at dual_cell::Int of vsphere::VoronoiSphere by an area-weighted average (first-order accurate).

qi is a scalar field known at primal cells. qv is a scalar field known at dual cells.

@inbounds propagates into average_iv and avg.

avg = average_ve(vsphere, edge)
qe[edge] = avg(qv)         # single-layer, qv::AbstractVector
qe[k, edge] = avg(qv, k)   # multi-layer, qv::AbstractMatrix

Interpolate scalar field at edge::Int of vsphere::VoronoiSphere by a centered average (first-order accurate).

qv is a scalar field known at dual cells. qe is a scalar field known at edges.

@inbounds propagates into average_ve and avg.

cflux = centered_flux(vsphere, edge)
flux[edge] = cflux(mass, ucov)         # single-layer, ucov::AbstractVector
flux[k, edge] = cflux(mass, ucov, k)   # multi-layer, ucov::AbstractMatrix

Compute centered flux at edge::Int of vsphere::VoronoiSphere, with respect to the unit sphere.

mass is a scalar field known at primal cells. ucov is a covariant vector field known at edges. flux is a contravariant vector field known at edges.

If ucov is defined with respect to a physical metric (e.g. in m²⋅s⁻¹) which is conformal, multiply cflux by the contravariant physical metric factor (in m⁻²). mass being e.g. in kg, on gets a flux in kg⋅s⁻¹ which can be fed into divergence.

@inbounds propagates into centered_flux and cflux.

op = curl(vsphere, dual_cell)
curlu[dual_cell] = op(ucov)         # single-layer, ucov::AbstractVector
curlu[k, dual_cell] = op(ucov, k)   # multi-layer, ucov::AbstractMatrix

Compute curl of ucov at dual_cell::Int of vsphere::VoronoiSphere. ucov is a covariant vector field known at edges. curlu is a density (two-form) over dual cells. To obtain a scalar, divide curlu by the dual cell area vsphere.Av

@inbounds propagates into curl and op.

div = divergence(vsphere, cell, Val(N))
dvg[cell] = div(flux) # single-layer, flux::AbstractVector
dvg[k, cell] = div(flux, k)  # multi-layer, flux::AbstractMatrix

Compute divergence with respect to the unit sphere of flux at cell::Int of vsphere::VoronoiSphere. flux is a contravariant vector field known at edges. dvg is a scalar field known at primal cells.

N=sphere.primal_deg[cell] is the number of cell edges and must be provided as a compile-time constant for performance. This may be done via the macro @unroll from ManagedLoops.

@inbounds propagates into divergence and div.

vsphere = dot_product(vsphere::VoronoiSphere) # optional, returns only relevant fields as a named tuple
dot_prod = dot_product(vsphere, cell::Int, v::Val{N})

# single-layer, `ucov` and `vcov` are ::AbstractVector
dp[cell] = dot_prod(ucov, vcov) 

# multi-layer, `ucov` and `vcov` are ::AbstractMatrix
dp[k, cell] = dot_prod(ucov, vcov, k)

Compute dot product with respect to the unit sphere of ucov, vcov at cell::Int of vsphere::VoronoiSphere. ucov and vcov are covariant vector fields known at edges.

N=sphere.primal_deg[cell] is the number of cell edges and must be provided as a compile-time constant for performance. This may be done via the macro @unroll from ManagedLoops.

@inbounds propagates into dot_product and dot_prod.

grad = gradient(vsphere, edge)
gradcov[edge] = grad(q)         # single-layer, q::AbstractVector
gradcov[k, edge] = grad(q, k)   # multi-layer, q::AbstractMatrix

Compute gradient of q at edge::Int of vsphere::VoronoiSphere.

q is a scalar field known at primal cells. gradcov is a covariant vector field known at edges. gradcov is numerically zero-curl.

@inbounds propagates into gradient and div.

grad = gradient3d(vsphere, cell, Val(N))
gradq[ij] = grad(q)         # single-layer, q::AbstractVector
gradq[k, ij] = grad(q, k)   # multi-layer, q::AbstractMatrix

Compute 3D gradient of q at cell::Int of vsphere::VoronoiSphere. q is a scalar field known at primal cells. gradq is a 3D vector field yielding a 3-uple at each primal cell.

N=sphere.primal_deg[cell] is the number of cell edges and must be provided as a compile-time constant for performance. This may be done via the macro @unroll from ManagedLoops.

@inbounds propagates into gradient3d and grad.

grad = gradperp(vsphere, edge)
flux[edge] = grad(psi)         # single-layer, q::AbstractVector
flux[k, edge] = grad(psi, k)   # multi-layer, q::AbstractMatrix

Compute grad⟂ of streamfunction psi at edge::Int of vsphere::VoronoiSphere.

psi is a scalar field known at dual cells. flux is a contravariant vector field known at edges. flux is numerically non-divergent.

@inbounds propagates into gradperp and grad.

op = perp(vsphere, ij)
U_perp[ij] = op(U)              # single-layer, U::AbstractVector
U_perp[k, ij] = op(U, k::Int)   # multi-layer, U::AbstractMatrix

Compute the perp operator U⟂ at edge::Int of vsphere::VoronoiSphere. Unlike the TRiSK operator, this operator is not antisymmetric but it has a smaller stencil and is numerically consistent.

Array U represents a vector field U by its components normal to edges of primal cells. U_perp represents similarly U⟂. Equivalently, it represents U by its components normal to edges of dual cells.

@inbounds propagates into perp and op.