Operators
Operators can compute spatial derivative operations.
- for performance reasons, we need to be able to "fuse" multiple operators and
function applications
- Julia provides a tool for this: broadcasting, with a very flexible API
Can think of operators are "pseudo-functions": can't be called directly, but act similar to functions in the context of broadcasting. They are matrix-free, in the sense that we define the action of the operator directly on a field, without explicitly assembling the matrix representing the discretized operator.
Spectral element operators
Differential Operators
ClimaCore.Operators.Gradient
— Typegrad = Gradient()
grad.(f)
Compute the (strong) gradient of f
on each element, returning a CovariantVector
-field.
The $i$th covariant component of the gradient is the partial derivative with respect to the reference element:
\[(\nabla f)_i = \frac{\partial f}{\partial \xi^i}\]
Discretely, this can be written in matrix form as
\[D_i f\]
where $D_i$ is the derivative matrix along the $i$th dimension.
References
- Mark A Taylor, Aimé Fournier (2010), equation 16
ClimaCore.Operators.Divergence
— Typediv = Divergence()
div.(u)
Computes the per-element spectral (strong) divergence of a vector field $u$.
The divergence of a vector field $u$ is defined as
\[\nabla \cdot u = \sum_i \frac{1}{J} \frac{\partial (J u^i)}{\partial \xi^i}\]
where $J$ is the Jacobian determinant, $u^i$ is the $i$th contravariant component of $u$.
This is discretized by
\[\sum_i I \left\{\frac{1}{J} \frac{\partial (I\{J u^i\})}{\partial \xi^i} \right\}\]
where $I\{x\}$ is the interpolation operator that projects to the unique polynomial interpolating $x$ at the quadrature points. In matrix form, this can be written as
\[J^{-1} \sum_i D_i J u^i\]
where $D_i$ is the derivative matrix along the $i$th dimension
References
- Mark A Taylor, Aimé Fournier (2010), equation 15
ClimaCore.Operators.WeakDivergence
— Typewdiv = WeakDivergence()
wdiv.(u)
Computes the "weak divergence" of a vector field u
.
This is defined as the scalar field $\theta \in \mathcal{V}_0$ such that for all $\phi\in \mathcal{V}_0$
\[\int_\Omega \phi \theta \, d \Omega = - \int_\Omega (\nabla \phi) \cdot u \,d \Omega\]
where $\mathcal{V}_0$ is the space of $u$.
This arises as the contribution of the volume integral after by applying integration by parts to the weak form expression of the divergence
\[\int_\Omega \phi (\nabla \cdot u) \, d \Omega = - \int_\Omega (\nabla \phi) \cdot u \,d \Omega + \oint_{\partial \Omega} \phi (u \cdot n) \,d \sigma\]
It can be written in matrix form as
\[ϕ^\top WJ θ = - \sum_i (D_i ϕ)^\top WJ u^i\]
which reduces to
\[θ = -(WJ)^{-1} \sum_i D_i^\top WJ u^i\]
where
- $J$ is the diagonal Jacobian matrix
- $W$ is the diagonal matrix of quadrature weights
- $D_i$ is the derivative matrix along the $i$th dimension
ClimaCore.Operators.WeakGradient
— Typewgrad = WeakGradient()
wgrad.(f)
Compute the "weak gradient" of f
on each element.
This is defined as the the vector field $\theta \in \mathcal{V}_0$ such that for all $\phi \in \mathcal{V}_0$
\[\int_\Omega \phi \cdot \theta \, d \Omega = - \int_\Omega (\nabla \cdot \phi) f \, d\Omega\]
where $\mathcal{V}_0$ is the space of $f$.
This arises from the contribution of the volume integral after by applying integration by parts to the weak form expression of the gradient
\[\int_\Omega \phi \cdot (\nabla f) \, d \Omega = - \int_\Omega f (\nabla \cdot \phi) \, d\Omega + \oint_{\partial \Omega} f (\phi \cdot n) \, d \sigma\]
In matrix form, this becomes
\[{\phi^i}^\top W J \theta_i = - ( J^{-1} D_i J \phi^i )^\top W J f\]
which reduces to
\[\theta_i = -W^{-1} D_i^\top W f\]
where $D_i$ is the derivative matrix along the $i$th dimension.
ClimaCore.Operators.Curl
— Typecurl = Curl()
curl.(u)
Computes the per-element spectral (strong) curl of a covariant vector field $u$.
Note: The vector field $u$ needs to be excliclty converted to a CovaraintVector
, as then the Curl
is independent of the local metric tensor.
The curl of a vector field $u$ is a vector field with contravariant components
\[(\nabla \times u)^i = \frac{1}{J} \sum_{jk} \epsilon^{ijk} \frac{\partial u_k}{\partial \xi^j}\]
where $J$ is the Jacobian determinant, $u_k$ is the $k$th covariant component of $u$, and $\epsilon^{ijk}$ are the Levi-Civita symbols. In other words
\[\begin{bmatrix} (\nabla \times u)^1 \\ (\nabla \times u)^2 \\ (\nabla \times u)^3 \end{bmatrix} = \frac{1}{J} \begin{bmatrix} \frac{\partial u_3}{\partial \xi^2} - \frac{\partial u_2}{\partial \xi^3} \\ \frac{\partial u_1}{\partial \xi^3} - \frac{\partial u_3}{\partial \xi^1} \\ \frac{\partial u_2}{\partial \xi^1} - \frac{\partial u_1}{\partial \xi^2} \end{bmatrix}\]
In matrix form, this becomes
\[\epsilon^{ijk} J^{-1} D_j u_k\]
Note that unused dimensions will be dropped: e.g. the 2D curl of a Covariant12Vector
-field will return a Contravariant3Vector
.
References
- Mark A Taylor, Aimé Fournier (2010), equation 17
ClimaCore.Operators.WeakCurl
— Typewcurl = WeakCurl()
wcurl.(u)
Computes the "weak curl" on each element of a covariant vector field u
.
Note: The vector field $u$ needs to be excliclty converted to a CovaraintVector
, as then the WeakCurl
is independent of the local metric tensor.
This is defined as the vector field $\theta \in \mathcal{V}_0$ such that for all $\phi \in \mathcal{V}_0$
\[\int_\Omega \phi \cdot \theta \, d \Omega = \int_\Omega (\nabla \times \phi) \cdot u \,d \Omega\]
where $\mathcal{V}_0$ is the space of $f$.
This arises from the contribution of the volume integral after by applying integration by parts to the weak form expression of the curl
\[\int_\Omega \phi \cdot (\nabla \times u) \,d\Omega = \int_\Omega (\nabla \times \phi) \cdot u \,d \Omega - \oint_{\partial \Omega} (\phi \times u) \cdot n \,d\sigma\]
In matrix form, this becomes
\[{\phi_i}^\top W J \theta^i = (J^{-1} \epsilon^{kji} D_j \phi_i)^\top W J u_k\]
which, by using the anti-symmetry of the Levi-Civita symbol, reduces to
\[\theta^i = - \epsilon^{ijk} (WJ)^{-1} D_j^\top W u_k\]
Interpolation Operators
ClimaCore.Operators.Interpolate
— Typei = Interpolate(space)
i.(f)
Interpolates f
to the space
. If space
has equal or higher polynomial degree as the space of f
, this is exact, otherwise it will be lossy.
In matrix form, it is the linear operator
\[I = \bigotimes_i I_i\]
where $I_i$ is the barycentric interpolation matrix in the $i$th dimension.
See also Restrict
.
ClimaCore.Operators.Restrict
— Typer = Restrict(space)
r.(f)
Computes the projection of a field f
on $\mathcal{V}_0$ to a lower degree polynomial space space
($\mathcal{V}_0^*$). space
must be on the same topology as the space of f
, but have a lower polynomial degree.
It is defined as the field $\theta \in \mathcal{V}_0^*$ such that for all $\phi \in \mathcal{V}_0^*$
\[\int_\Omega \phi \theta \,d\Omega = \int_\Omega \phi f \,d\Omega\]
In matrix form, this is
\[\phi^\top W^* J^* \theta = (I \phi)^\top WJ f\]
where $W^*$ and $J^*$ are the quadrature weights and Jacobian determinant of $\mathcal{V}_0^*$, and $I$ is the interpolation operator (see Interpolate
) from $\mathcal{V}_0^*$ to $\mathcal{V}_0$. This reduces to
\[\theta = (W^* J^*)^{-1} I^\top WJ f\]
DSS
ClimaCore.Spaces.weighted_dss!
— FunctionSpaces.weighted_dss!(f::Field[, ghost_buffer = Spaces.create_ghost_buffer(field)])
Apply weighted direct stiffness summation (DSS) to f
. This operates in-place (i.e. it modifies the f
). ghost_buffer
contains the necessary information for communication in a distributed setting, see Spaces.create_ghost_buffer
.
This is a projection operation from the piecewise polynomial space $\mathcal{V}_0$ to the continuous space $\mathcal{V}_1 = \mathcal{V}_0 \cap \mathcal{C}_0$, defined as the field $\theta \in \mathcal{V}_1$ such that for all $\phi \in \mathcal{V}_1$
\[\int_\Omega \phi \theta \,d\Omega = \int_\Omega \phi f \,d\Omega\]
In matrix form, we define $\bar \theta$ to be the unique global node representation, and $Q$ to be the "scatter" operator which maps to the redundant node representation $\theta$
\[\theta = Q \bar \theta\]
Then the problem can be written as
\[(Q \bar\phi)^\top W J Q \bar\theta = (Q \bar\phi)^\top W J f\]
which reduces to
\[\theta = Q \bar\theta = Q (Q^\top W J Q)^{-1} Q^\top W J f\]
weighted_dss!(data, space, ghost_buffer = nothing)
Computes weighted dss of data
. This function consists of
1). weighted_dss_start!
which loads the send buffer and starts communications,
2). weighted_dss_internal!
which progresses the communication and performs dss operations on interior vertices and faces, and
3). weighted_dss_ghost!
which completes the communication and performs dss on ghost vertices and faces
These constituent functions can also be called separately for increased overlap between computation and communication when merging multiple dss operations.
ClimaCore.Spaces.create_ghost_buffer
— FunctionSpaces.create_ghost_buffer(field::Field)
Create a buffer for communicating neighbour information of field
.
ClimaCore.Spaces.weighted_dss_start!
— Functionweighted_dss_start!(data, space, ghost_buffer)
Create the ghost buffer, if necessary, load the send buffer and start communication.
Part of Spaces.weighted_dss!
.
ClimaCore.Spaces.weighted_dss_internal!
— Functionweighted_dss_internal!(data, space, ghost_buffer)
Perform weighted dss for internal faces and vertices. Progress the communication.
Part of Spaces.weighted_dss!
.
ClimaCore.Spaces.weighted_dss_ghost!
— Functionweighted_dss_ghost!(data, space, ghost_buffer)
Complete communication. Perform weighted dss for ghost faces and vertices.
Part of Spaces.weighted_dss!
.
Finite difference operators
Finite difference operators are similar with some subtle differences:
- they can change staggering (center to face, or vice versa)
- they can span multiple elements
- no DSS is required
- boundary handling may be required
We use the following convention:
- centers are indexed by integers
1, 2, ..., n
- faces are indexed by half integers
half, 1+half, ..., n+half
ClimaCore.Operators.FiniteDifferenceOperator
— TypeFiniteDifferenceOperator
An abstract type for finite difference operators. Instances of this should define:
See also BoundaryCondition
for how to define the boundaries.
Interpolation operators
ClimaCore.Operators.InterpolateC2F
— TypeI = InterpolateC2F(;boundaries..)
I.(x)
Interpolate a center-valued field x
to faces, using the stencil
\[I(x)[i] = \frac{1}{2} (x[i+\tfrac{1}{2}] + x[i-\tfrac{1}{2}])\]
Supported boundary conditions are:
SetValue(x₀)
: set the value at the boundary face to bex₀
. On the left boundary the stencil is
\[I(x)[\tfrac{1}{2}] = x₀\]
SetGradient(v)
: set the value at the boundary such that the gradient isv
. At the left boundary the stencil is
\[I(x)[\tfrac{1}{2}] = x[1] - \frac{1}{2} v³\]
Extrapolate
: use the closest interior point as the boundary value. At the left boundary the stencil is
\[I(x)[\tfrac{1}{2}] = x[1]\]
ClimaCore.Operators.InterpolateF2C
— TypeInterpolateF2C()
Interpolate from face to center mesh. No boundary conditions are required (or supported).
ClimaCore.Operators.WeightedInterpolateC2F
— TypeWI = WeightedInterpolateC2F(; boundaries)
WI.(w, x)
Interpolate a center-valued field x
to faces, weighted by a center-valued field w
, using the stencil
\[WI(w, x)[i] = \frac{ w[i+\tfrac{1}{2}] x[i+\tfrac{1}{2}] + w[i-\tfrac{1}{2}] x[i-\tfrac{1}{2}]) }{ 2 (w[i+\tfrac{1}{2}] + w[i-\tfrac{1}{2}]) }\]
Supported boundary conditions are:
SetValue(val)
: set the value at the boundary face to beval
.SetGradient
: set the value at the boundary such that the gradient isval
.Extrapolate
: use the closest interior point as the boundary value
ClimaCore.Operators.WeightedInterpolateF2C
— TypeWI = WeightedInterpolateF2C(; boundaries)
WI.(w, x)
Interpolate a face-valued field x
to centers, weighted by a face-valued field w
, using the stencil
\[WI(w, x)[i] = \frac{ w[i+\tfrac{1}{2}] x[i+\tfrac{1}{2}] + w[i-\tfrac{1}{2}] x[i-\tfrac{1}{2}]) }{ 2 (w[i+\tfrac{1}{2}] + w[i-\tfrac{1}{2}]) }\]
No boundary conditions are required (or supported)
ClimaCore.Operators.UpwindBiasedProductC2F
— TypeU = UpwindBiasedProductC2F(;boundaries)
U.(v, x)
Compute the product of the face-valued vector field v
and a center-valued field x
at cell faces by upwinding x
according to the direction of v
.
More precisely, it is computed based on the sign of the 3rd contravariant component, and it returns a Contravariant3Vector
:
\[U(\boldsymbol{v},x)[i] = \begin{cases} v^3[i] x[i-\tfrac{1}{2}]\boldsymbol{e}_3 \textrm{, if } v^3[i] > 0 \\ v^3[i] x[i+\tfrac{1}{2}]\boldsymbol{e}_3 \textrm{, if } v^3[i] < 0 \end{cases}\]
where $\boldsymbol{e}_3$ is the 3rd covariant basis vector.
Supported boundary conditions are:
SetValue(x₀)
: set the value ofx
to bex₀
in a hypothetical ghost cell on the other side of the boundary. On the left boundary the stencil is\[U(\boldsymbol{v},x)[\tfrac{1}{2}] = \begin{cases} v^3[\tfrac{1}{2}] x_0 \boldsymbol{e}_3 \textrm{, if } v^3[\tfrac{1}{2}] > 0 \\ v^3[\tfrac{1}{2}] x[1] \boldsymbol{e}_3 \textrm{, if } v^3[\tfrac{1}{2}] < 0 \end{cases}\]
Extrapolate()
: set the value ofx
to be the same as the closest interior point. On the left boundary, the stencil is\[U(\boldsymbol{v},x)[\tfrac{1}{2}] = U(\boldsymbol{v},x)[1 + \tfrac{1}{2}]\]
ClimaCore.Operators.Upwind3rdOrderBiasedProductC2F
— TypeU = Upwind3rdOrderBiasedProductC2F(;boundaries)
U.(v, x)
Compute the product of a face-valued vector field v
and a center-valued field x
at cell faces by upwinding x
, to third-order of accuracy, according to v
\[U(v,x)[i] = \begin{cases} v[i] \left(-2 x[i-\tfrac{3}{2}] + 10 x[i-\tfrac{1}{2}] + 4 x[i+\tfrac{1}{2}] \right) / 12 \textrm{, if } v[i] > 0 \\ v[i] \left(4 x[i-\tfrac{1}{2}] + 10 x[i+\tfrac{1}{2}] -2 x[i+\tfrac{3}{2}] \right) / 12 \textrm{, if } v[i] < 0 \end{cases}\]
This stencil is based on Louis J. Wicker, William C. Skamarock (2002), eq. 4(a).
Supported boundary conditions are:
FirstOrderOneSided(x₀)
: uses the first-order downwind scheme to computex
on the left boundary, and the first-order upwind scheme to computex
on the right boundary.ThirdOrderOneSided(x₀)
: uses the third-order downwind reconstruction to computex
on the left boundary,
and the third-order upwind reconstruction to compute x
on the right boundary.
These boundary conditions do not define the value at the actual boundary faces, and so this operator cannot be materialized directly: it needs to be composed with another operator that does not make use of this value, e.g. a DivergenceF2C
operator, with a SetValue
boundary.
ClimaCore.Operators.FCTBorisBook
— TypeU = FCTBorisBook(;boundaries)
U.(v, x)
Correct the flux using the flux-corrected transport formulation by Boris and Book Jay P Boris, David L Book (1973).
Input arguments:
- a face-valued vector field
v
- a center-valued field
x
\[Ac(v,x)[i] = s[i] \max \left\{0, \min \left[ |v[i] |, s[i] \left( x[i+\tfrac{3}{2}] - x[i+\tfrac{1}{2}] \right) , s[i] \left( x[i-\tfrac{1}{2}] - x[i-\tfrac{3}{2}] \right) \right] \right\},\]
where $s[i] = +1$ if $v[i] \geq 0$ and $s[i] = -1$ if $v[i] \leq 0$, and $Ac$ represents the resulting corrected antidiffusive flux. This formulation is based on Jay P Boris, David L Book (1973), as reported in Dale R Durran (2010) section 5.4.1.
Supported boundary conditions are:
FirstOrderOneSided(x₀)
: uses the first-order downwind reconstruction to computex
on the left boundary, and the first-order upwind reconstruction to computex
on the right boundary.
Similar to the Upwind3rdOrderBiasedProductC2F
operator, these boundary conditions do not define the value at the actual boundary faces, and so this operator cannot be materialized directly: it needs to be composed with another operator that does not make use of this value, e.g. a DivergenceF2C
operator, with a SetValue
boundary.
ClimaCore.Operators.FCTZalesak
— TypeU = FCTZalesak(;boundaries)
U.(A, Φ, Φᵗᵈ)
Correct the flux using the flux-corrected transport formulation by Zalesak Steven T Zalesak (1979).
Input arguments:
- a face-valued vector field
A
- a center-valued field
Φ
- a center-valued field
Φᵗᵈ
\[Φ_j^{n+1} = Φ_j^{td} - (C_{j+\frac{1}{2}}A_{j+\frac{1}{2}} - C_{j-\frac{1}{2}}A_{j-\frac{1}{2}})\]
This stencil is based on Steven T Zalesak (1979), as reported in Dale R Durran (2010) section 5.4.2, where $C$ denotes the corrected antidiffusive flux.
Supported boundary conditions are:
FirstOrderOneSided(x₀)
: uses the first-order downwind reconstruction to computex
on the left boundary, and the first-order upwind reconstruction to computex
on the right boundary.
Similar to the Upwind3rdOrderBiasedProductC2F
operator, these boundary conditions do not define the value at the actual boundary faces, and so this operator cannot be materialized directly: it needs to be composed with another operator that does not make use of this value, e.g. a DivergenceF2C
operator, with a SetValue
boundary.
ClimaCore.Operators.LeftBiasedC2F
— TypeL = LeftBiasedC2F(;boundaries)
L.(x)
Interpolate a center-value field to a face-valued field from the left.
\[L(x)[i] = x[i-\tfrac{1}{2}]\]
Only the left boundary condition should be set. Currently supported is:
SetValue(x₀)
: set the value to bex₀
on the boundary.
\[L(x)[\tfrac{1}{2}] = x_0\]
ClimaCore.Operators.RightBiasedC2F
— TypeR = RightBiasedC2F(;boundaries)
R.(x)
Interpolate a center-valued field to a face-valued field from the right.
\[R(x)[i] = x[i+\tfrac{1}{2}]\]
Only the right boundary condition should be set. Currently supported is:
SetValue(x₀)
: set the value to bex₀
on the boundary.
\[R(x)[n+\tfrac{1}{2}] = x_0\]
ClimaCore.Operators.LeftBiasedF2C
— TypeL = LeftBiasedF2C(;boundaries)
L.(x)
Interpolate a face-value field to a center-valued field from the left.
\[L(x)[i+\tfrac{1}{2}] = x[i]\]
Only the left boundary condition should be set. Currently supported is:
SetValue(x₀)
: set the value to bex₀
on the boundary.
\[L(x)[1] = x_0\]
ClimaCore.Operators.RightBiasedF2C
— TypeR = RightBiasedF2C(;boundaries)
R.(x)
Interpolate a face-valued field to a center-valued field from the right.
\[R(x)[i] = x[i+\tfrac{1}{2}]\]
Only the right boundary condition should be set. Currently supported is:
SetValue(x₀)
: set the value to bex₀
on the boundary.
\[R(x)[n+\tfrac{1}{2}] = x_0\]
Derivative operators
ClimaCore.Operators.GradientF2C
— TypeG = GradientF2C(;boundaryname=boundarycondition...)
G.(x)
Compute the gradient of a face-valued field x
, returning a center-valued Covariant3
vector field, using the stencil:
\[G(x)[i]^3 = x[i+\tfrac{1}{2}] - x[i-\tfrac{1}{2}]\]
We note that the usual division factor $1 / \Delta z$ that appears in a first-order finite difference operator is accounted for in the LocalVector
basis. Hence, users need to cast the output of the GradientF2C
to a UVector
, VVector
or WVector
, according to the type of domain on which the operator is defined.
The following boundary conditions are supported:
- by default, the value of
x
at the boundary face will be used. SetValue(x₀)
: calculate the gradient assuming the value at the boundary isx₀
. For the left boundary, this becomes:
\[G(x)[1]³ = x[1+\tfrac{1}{2}] - x₀\]
Extrapolate()
: set the value at the center closest to the boundary
to be the same as the neighbouring interior value. For the left boundary, this becomes:
\[G(x)[1]³ = G(x)[2]³\]
ClimaCore.Operators.GradientC2F
— TypeG = GradientC2F(;boundaryname=boundarycondition...)
G.(x)
Compute the gradient of a center-valued field x
, returning a face-valued Covariant3
vector field, using the stencil:
\[G(x)[i]^3 = x[i+\tfrac{1}{2}] - x[i-\tfrac{1}{2}]\]
The following boundary conditions are supported:
SetValue(x₀)
: calculate the gradient assuming the value at the boundary isx₀
. For the left boundary, this becomes:\[G(x)[\tfrac{1}{2}]³ = 2 (x[1] - x₀)\]
SetGradient(v₀)
: set the value of the gradient at the boundary to bev₀
. For the left boundary, this becomes:\[G(x)[\tfrac{1}{2}] = v₀\]
ClimaCore.Operators.AdvectionF2F
— TypeA = AdvectionF2F(;boundaries)
A.(v, θ)
Vertical advection operator at cell faces, for a face-valued velocity field v
and face-valued variables θ
, approximating $v^3 \partial_3 \theta$.
It uses the following stencil
\[A(v,θ)[i] = \frac{1}{2} (θ[i+1] - θ[i-1]) v³[i]\]
No boundary conditions are currently supported.
ClimaCore.Operators.AdvectionC2C
— TypeA = AdvectionC2C(;boundaries)
A.(v, θ)
Vertical advection operator at cell centers, for cell face velocity field v
cell center variables θ
, approximating $v^3 \partial_3 \theta$.
It uses the following stencil
\[A(v,θ)[i] = \frac{1}{2} \{ (θ[i+1] - θ[i]) v³[i+\tfrac{1}{2}] + (θ[i] - θ[i-1])v³[i-\tfrac{1}{2}]\}\]
Supported boundary conditions:
SetValue(θ₀)
: set the value ofθ
at the boundary face to beθ₀
. At the lower boundary, this is:
\[A(v,θ)[1] = \frac{1}{2} \{ (θ[2] - θ[1]) v³[1 + \tfrac{1}{2}] + (θ[1] - θ₀)v³[\tfrac{1}{2}]\}\]
Extrapolate
: use the closest interior point as the boundary value. At the lower boundary, this is:
\[A(v,θ)[1] = (θ[2] - θ[1]) v³[1 + \tfrac{1}{2}] \}\]
ClimaCore.Operators.DivergenceF2C
— TypeD = DivergenceF2C(;boundaryname=boundarycondition...)
D.(v)
Compute the vertical contribution to the divergence of a face-valued field vector v
, returning a center-valued scalar field, using the stencil
\[D(v)[i] = (Jv³[i+\tfrac{1}{2}] - Jv³[i-\tfrac{1}{2}]) / J[i]\]
where Jv³
is the Jacobian multiplied by the third contravariant component of v
.
The following boundary conditions are supported:
- by default, the value of
v
at the boundary face will be used. SetValue(v₀)
: calculate the divergence assuming the value at the boundary isv₀
. For the left boundary, this becomes:
\[D(v)[1] = (Jv³[1+\tfrac{1}{2}] - Jv³₀) / J[i]\]
Extrapolate()
: set the value at the center closest to the boundary to be the same as the neighbouring interior value. For the left boundary, this becomes:
\[D(v)[1]³ = D(v)[2]³\]
ClimaCore.Operators.DivergenceC2F
— TypeD = DivergenceC2F(;boundaryname=boundarycondition...)
D.(v)
Compute the vertical contribution to the divergence of a center-valued field vector v
, returning a face-valued scalar field, using the stencil
\[D(v)[i] = (Jv³[i+\tfrac{1}{2}] - Jv³[i-\tfrac{1}{2}]) / J[i]\]
where Jv³
is the Jacobian multiplied by the third contravariant component of v
.
The following boundary conditions are supported:
SetValue(v₀)
: calculate the divergence assuming the value at the boundary isv₀
. For the left boundary, this becomes:\[D(v)[\tfrac{1}{2}] = \frac{1}{2} (Jv³[1] - Jv³₀) / J[i]\]
SetDivergence(x)
: set the value of the divergence at the boundary to bex
.\[D(v)[\tfrac{1}{2}] = x\]
ClimaCore.Operators.CurlC2F
— TypeC = CurlC2F(;boundaryname=boundarycondition...)
C.(v)
Compute the vertical-derivative contribution to the curl of a center-valued covariant vector field v
. It acts on the horizontal covariant components of v
(that is it only depends on $v₁$ and $v₂$), and will return a face-valued horizontal contravariant vector field (that is $C(v)³ = 0$).
Specifically it approximates:
\[\begin{align*} C(v)^1 &= -\frac{1}{J} \frac{\partial v_2}{\partial \xi^3} \\ C(v)^2 &= \frac{1}{J} \frac{\partial v_1}{\partial \xi^3} \\ \end{align*}\]
using the stencils
\[\begin{align*} C(v)[i]^1 &= - \frac{1}{J[i]} (v₂[i+\tfrac{1}{2}] - v₂[i-\tfrac{1}{2}]) \\ C(v)[i]^2 &= \frac{1}{J[i]} (v₁[i+\tfrac{1}{2}] - v₁[i-\tfrac{1}{2}]) \end{align*}\]
where $v₁` and$v₂$are the 1st and 2nd covariant components of$v$, and$J`` is the Jacobian determinant.
The following boundary conditions are supported:
SetValue(v₀)
: calculate the curl assuming the value of $v$ at the boundary isv₀
. For the left boundary, this becomes:\[C(v)[\tfrac{1}{2}]^1 = -\frac{2}{J[i]} (v_2[1] - (v₀)_2) C(v)[\tfrac{1}{2}]^2 = \frac{2}{J[i]} (v_1[1] - (v₀)_1)\]
SetCurl(v⁰)
: enforce the curl operator output at the boundary to be the contravariant vectorv⁰
.
Other
ClimaCore.Operators.SetBoundaryOperator
— TypeSetBoundaryOperator(;boundaries...)
This operator only modifies the values at the boundary:
SetValue(val)
: set the value to beval
on the boundary.
ClimaCore.Operators.FirstOrderOneSided
— TypeFirstOrderOneSided()
Use a first-order up/down-wind scheme to compute the value at the boundary.
ClimaCore.Operators.ThirdOrderOneSided
— TypeThirdOrderOneSided()
Use a third-order up/down-wind scheme to compute the value at the boundary.
Finite difference boundary conditions
ClimaCore.Operators.BoundaryCondition
— TypeBoundaryCondition
An abstract type for boundary conditions for FiniteDifferenceOperator
s.
Subtypes should define:
ClimaCore.Operators.SetValue
— TypeSetValue(val)
Set the value at the boundary to be val
. In the case of gradient operators, this will set the input value from which the gradient is computed.
ClimaCore.Operators.SetGradient
— TypeSetGradient(val)
Set the gradient at the boundary to be val
. In the case of gradient operators this will set the output value of the gradient.
ClimaCore.Operators.SetDivergence
— TypeSetDivergence(val)
Set the divergence at the boundary to be val
.
ClimaCore.Operators.Extrapolate
— TypeExtrapolate()
Set the value at the boundary to be the same as the closest interior point.
Internal APIs
ClimaCore.Operators.return_eltype
— Functionreturn_eltype(::Op, fields...)
Defines the element type of the result of operator Op
ClimaCore.Operators.return_space
— Functionreturn_space(::Op, spaces...)
Defines the space the operator Op
returns values on.
ClimaCore.Operators.stencil_interior_width
— Functionstencil_interior_width(::Op, args...)
Defines the width of the interior stencil for the operator Op
with the given arguments. Returns a tuple of 2-tuples: each 2-tuple should be the lower and upper bounds of the index offsets of the stencil for each argument in the stencil.
Example
stencil(::Op, arg1, arg2) = ((-half, 1+half), (0,0))
implies that at index i
, the stencil accesses arg1
at i-half
, i+half
and i+1+half
, and arg2
at index i
.
ClimaCore.Operators.stencil_interior
— Functionstencil_interior(::Op, loc, idx, args...)
Defines the stencil of the operator Op
in the interior of the domain at idx
; args
are the input arguments.
ClimaCore.Operators.boundary_width
— Functionboundary_width(::Op, ::BC, args...)
Defines the width of a boundary condition BC
on an operator Op
. This is the number of locations that are used in a modified stencil.
ClimaCore.Operators.stencil_left_boundary
— Functionstencil_left_boundary(::Op, ::BC, loc, idx, args...)
Defines the stencil of operator Op
at idx
near the left boundary, with boundary condition BC
.
ClimaCore.Operators.stencil_right_boundary
— Functionstencil_right_boundary(::Op, ::BC, loc, idx, args...)
Defines the stencil of operator Op
at idx
near the right boundary, with boundary condition BC
.