Spaces

Coordinate axis labels used in documentation and error messages.

Corresponds to the x, y, and z spatial dimensions (1st, 2nd, and 3rd dimensions).

See also: _BRAMBLE_var2symbol

Subscript Unicode symbols for x, y, z coordinates used in operator notation.

These symbols are used to generate directional operator aliases via metaprogramming.

Examples

• D₊ₓ - forward difference in x-direction
• M₋ᵧ - backward average in y-direction
• jump₂ - jump in the z-direction

See also: _BRAMBLE_var2label

AbstractSpaceType{N}

Abstract supertype for all function spaces defined on a mesh.

This is the top-level abstraction for a grid-based function space. The parameter N represents the number of components of the field (e.g., N=1 for a scalar field, N=3 for a 3D vector field).

InnerProductType

Abstract type for selecting which discrete inner product formula to use.

Different inner product types correspond to different weight distributions on the grid, used in various finite difference schemes and stability analyses. The choice of inner product affects energy estimates and numerical stability properties.

Subtypes

• Innerh: Standard $L^2$ inner product using cell measures (volumes)
• Innerplus: Modified inner product using staggered grid spacings

Background

In finite difference methods, different inner products arise naturally from:

• Summation-by-parts (SBP) operators
• Energy method stability analysis
• Discrete integration formulas

The standard inner product (Innerh) uses cell volumes as weights, while the modified inner products (Innerplus) use combinations of forward/backward spacings, appearing in discrete energy estimates for difference operators.

Usage

# Compute standard L² inner product
result = innerₕ(uₕ, vₕ)  # Uses Innerh() internally

# Compute modified inner product in x-direction
result = inner₊ₓ(uₕ, vₕ)  # Uses Innerplus() internally

See also: Innerh, Innerplus, innerₕ, inner₊ₓ

Innerh <: InnerProductType

Selector for the standard discrete $L^2$ inner product weighted by cell measures.

The weights are the volumes (1D: lengths, 2D: areas, 3D: volumes) of grid cells, denoted $|\square_k|$. This is the most common inner product for finite difference methods and corresponds to the trapezoid rule for integration on non-uniform grids.

Mathematical form

For a 2D grid:

\[(u_h, v_h)_h = \sum_{i,j} |\square_{i,j}| u_h(x_i, y_j) v_h(x_i, y_j)\]

where $|\square_{i,j}|$ is the area of the cell centered at $(x_i, y_j)$.

Example

# Compute L² inner product
result = innerₕ(uₕ, vₕ)  # Uses Innerh() internally

# Compute L² norm
norm_value = normₕ(uₕ)  # Equivalent to sqrt(innerₕ(uₕ, uₕ))

See also: InnerProductType, Innerplus, innerₕ, normₕ

Innerplus <: InnerProductType

Selector for modified discrete $L^2$ inner products using staggered grid spacings.

These inner products use a combination of forward spacings $h_i$ and centered cell widths $h_{i+1/2}$, appearing naturally in energy estimates for finite difference operators. Different spatial directions may have different weight formulas.

The modified inner products are used for:

• Proving discrete energy stability
• Analyzing discrete conservation properties
• Constructing stable finite difference schemes

Mathematical form

For a 2D grid in the x-direction:

\[(u_h, v_h)_{+x} = \sum_{i,j} h_{x,i} h_{y,j+1/2} u_h(x_i, y_j) v_h(x_i, y_j)\]

Example

# These functions use Innerplus internally
result_x = inner₊ₓ(uₕ, vₕ)  # Modified inner product, x-direction
result_y = inner₊ᵧ(uₕ, vₕ)  # Modified inner product, y-direction

See also: InnerProductType, Innerh, inner₊ₓ, weights

SpaceWeights(innerh::VT, innerplus::NTuple{D, VT})
SpaceWeights{D, VT}(innerh::VT, innerplus::NTuple{D, VT})

Holds the diagonal weight vectors for a grid space's discrete inner products, both the standard $L^2$ weights and the staggered ones, precomputed once rather than recomputed on every call.

Fields

• innerh::VT: weight vector for the standard discrete $L^2$ inner product (:innerₕ), based on cell measures ($|\square_k|$).
• innerplus::NTuple{D, VT}: tuple of weight vectors for modified, staggered inner products (:inner₊ₓ, :inner₊ᵧ, etc.), with one vector for each spatial dimension.

For a detailed explanation of the mathematical formulas corresponding to these weights, please refer to the documentation for ScalarGridSpace.

eltype(Wₕ::ScalarGridSpace) -> Type
eltype(::Type{<:ScalarGridSpace}) -> Type

Returns the element type of vectors in this space (e.g., Float64).

See also: backend

__innerplus_weights!(policy, v, innerplus_per_component)

Builds the weights for the modified discrete $L^2$ inner product on the space of grid functions ScalarGridSpace. The result is stored in vector v.

_innerh_weights!(u, Ωₕ::AbstractMeshType)

Builds the weights for the standard discrete $L^2$ inner product, $inner_h(\cdot, \cdot)$, on the space of grid functions, following the order of the points provided by indices(Ωₕ). The values are stored in vector u.

_innerplus_mean_weights!(u::VT, Ωₕ, component::Int = 1) where VT

Builds a set of weights based on the half spacings, associated with the component-th direction, for the modified discrete $L^2$ inner product on the space of grid functions, following the order of the points. The values are stored in vector u.

_innerplus_weights!(u::VT, Ωₕ, component = 1) where VT

Builds a set of weights based on the spacings, associated with the component-th direction, for the modified discrete $L^2$ inner product on the space of grid functions, following the order of the points provided by indices(Ωₕ). The values are stored in vector u.

^(Wₕ::ScalarGridSpace, ::Val{N}) where N -> CompositeGridSpace{N}
^(Wₕ::ScalarGridSpace, N::Int) -> CompositeGridSpace{N}

Constructs an N-component vector grid space from a scalar grid space Wₕ using mathematical exponentiation syntax: Vₕ = Wₕ^2 or Vₕ = Wₕ^dim(mesh).

Wₕ^1 is Wₕ, for both the Int and Val spellings.

The Int spelling only accepts 1 <= N <= 3, each branch-unswitched to a literal Wₕ^Val(N) call: Julia's return-type inference only union-splits up to 3 concrete types, so this is the largest range that stays inferrable even when N is a runtime value the compiler cannot constant-fold (e.g. threaded through a generic function argument). It also covers every mesh dimension this package supports (Wₕ^dim(mesh)). For N > 3, call Wₕ^Val(N) directly, which has no upper bound.

_shares_one_mesh(comps::Tuple) -> Bool

Whether every leaf element in comps (a composite's components(uₕ)) sits on the exact same mesh object.

A composite's leaves need not share a mesh — heterogeneous composites, whose leaves are built over differently-sized meshes, are a supported pattern (see the interpolation tutorial). When they do share one, a single evaluation of a vector-valued function at one leaf's grid points is valid for every leaf, which is what the scatter paths of Rₕ!/avgₕ! use this to decide; when they do not, each leaf needs its own evaluation, since there is no shared "grid point i" across differently-sized meshes.

leaf_count(Wₕ::AbstractSpaceType) -> Int

Returns the number of leaf scalar spaces in Wₕ. For a scalar space, returns 1. For a composite space, returns the total flattened number of scalar leaves.

leaf_spaces_offsets(Wₕ) -> Tuple

The scalar spaces underneath Wₕ paired with their offsets into the global degree of freedom vector, depth first and left to right, as a tuple of (space, offset).

A scalar space is its own only leaf, at offset zero.

getindex(uₕ::VectorElement{<:ScalarGridSpace{2}}, i::Integer, j::Integer)
getindex(uₕ::VectorElement{<:ScalarGridSpace{3}}, i::Integer, j::Integer, k::Integer)
getindex(uₕ::VectorElement{<:ScalarGridSpace{D}}, I::CartesianIndex{D}) where {D}
getindex(uₕ::VectorElement{<:ScalarGridSpace}, I::CartesianIndex)

Access field degrees of freedom by spatial grid coordinates or CartesianIndex.

Translates spatial grid coordinates directly into flat linear coefficient offsets using the mesh's LinearIndices with zero heap allocations and full @inbounds transparency.

Examples

Ωₕ = mesh(domain(interval(0.0, 1.0) × interval(0.0, 1.0)), (10, 10))
Wₕ = gridspace(Ωₕ)
uₕ = element(Wₕ, 0.0)

# Set and get via 2D coordinates
uₕ[2, 3] = 42.0
uₕ[2, 3] == 42.0

# Access via CartesianIndex
I = CartesianIndex(2, 3)
uₕ[I] == 42.0

See also: VectorElement, ScalarGridSpace, reshape

setindex!(uₕ::VectorElement{<:ScalarGridSpace{2}}, val, i::Integer, j::Integer) -> VectorElement
setindex!(uₕ::VectorElement{<:ScalarGridSpace{3}}, val, i::Integer, j::Integer, k::Integer) -> VectorElement
setindex!(uₕ::VectorElement{<:ScalarGridSpace{D}}, val, I::CartesianIndex{D}) where {D} -> VectorElement
setindex!(uₕ::VectorElement{<:ScalarGridSpace}, val, I::CartesianIndex) -> VectorElement

Mutate field degrees of freedom by spatial grid coordinates or CartesianIndex in-place.

Translates spatial grid coordinates directly into flat linear coefficient offsets using the mesh's LinearIndices with zero heap allocations and full @inbounds transparency.

Returns uₕ matching Base collection conventions.

_find_vec_in_broadcast(bc)

Internal helper to extract a VectorElement from a broadcast expression.

Recursively searches through the arguments of a broadcast expression tree to find a VectorElement instance. This is used by the broadcasting machinery to determine which function space should be used for the result.

Arguments

• bc: A broadcast expression, tuple of arguments, or individual value

Returns

• The first VectorElement found in the expression tree
• nothing if no VectorElement is found

Implementation Notes

Uses multiple dispatch to handle:

• Broadcasted objects: Extract and search arguments
• Tuples: Recursively search each element
• VectorElement: Return immediately (found!)
• Other types: Return nothing and continue searching

This enables broadcasts like uₕ .+ vₕ .* 2 to automatically preserve the space information.

CellAverage(f, nq::Val)

Project f by averaging it over each cell with nq quadrature points per direction. The rule behind avgₕ!.

PointValue(f)

Project f by evaluating it at each mesh point. The rule behind Rₕ!.

ProjectionRule

Supertype of the rules project! can apply: what a grid function's value at a point should be, given a continuous function. PointValue evaluates at the mesh point, CellAverage averages over the cell.

_rule_component(rule, k) -> ProjectionRule

rule restricted to leaf k of a composite space whose leaves do not share one mesh, so there is no shared grid point to evaluate once and scatter.

_rule_kernel(rule, space) -> callable

A concretely typed callable mapping a linear grid index to the scalar value rule gives at that point. Kept a named struct per rule rather than a closure: measured, an anonymous closure over the captures here takes a miscompiled path that allocates per grid point (gpena/Bramble.jl#64).

_rule_scatter_kernel(rule, space, ::Val{NC}) -> callable

As _rule_kernel, for a rule whose function returns all NC leaf values at once, so it is evaluated once per point and scattered across the leaves.

project!(uₕ::VectorElement, rule, markers = ()) -> uₕ

Project onto uₕ in place according to rule, optionally restricted to the union of the labelled marker regions, leaving every other entry zero.

Handles the space's shape and the execution policy; rule supplies only the per-point value. See PointValue and CellAverage, and Rₕ! / avgₕ! for the public spellings.

AVG_QUAD_POINTS

Default number of Gauss-Legendre points per direction, per cell, used by avgₕ. Six points are exact for polynomials up to degree eleven.

Unlike an adaptive rule, a fixed one does not tighten itself on coarse cells, so the default is chosen to be accurate on cells far coarser than any practical grid. Measured on 4 points spanning [-1, 4] with a function varying by a factor of e^5 across the domain (deliberately harsher than a real mesh), the worst error over 30 random grids was

points   1D        2D        3D        evaluations per cell (3D)
3        6.1e-5    1.6e-4    3.1e-4     27
4        2.8e-7    8.1e-7    3.4e-6     64
5        1.9e-9    4.1e-9    7.8e-9    125
6        5.5e-12   1.4e-11   1.6e-11   216

Cost is quad_points^D evaluations per cell. On a fine grid three points are usually ample; lower it with the quad_points keyword when the integrand is cheap to resolve and the cells are small.

_gauss_rule(::Val{N}, ::Type{T})

Returns (nodes, weights) for the N-point Gauss-Legendre rule on [0, 1] as NTuple{N,T}, so the per-cell loop that consumes them does not allocate.

The rule is built by QuadGK.gauss in the requested element type, so Float32 and BigFloat grids get a rule at their own precision rather than a rounded Float64 one. Weights sum to one, which makes the weighted sum over a cell the cell average directly.

⊗(A, B)

Kronecker product operator (alias for kron).

Computes the Kronecker product (tensor product) of matrices A and B. This operator is used extensively in constructing multidimensional shift operators.

Examples

I₂ = I(2)
I₃ = I(3)
result = I₂ ⊗ I₃  # 6×6 identity matrix

See also: shift

_Eye(be::Backend, npts, ::Val{i})

Internal helper to create identity or shifted diagonal matrices, in the matrix type the backend be chose, routed through backend_eye/matrix_type(be) rather than a package-specific lazy type, so the result always matches whatever matrix_type the caller's backend picked.

Arguments

• be::Backend: the backend whose matrix_type the result is built in.
• npts::Int: Size of the square matrix.
• ::Val{i}: Diagonal offset (0 = main diagonal, 1 = superdiagonal, -1 = subdiagonal).

Returns

• For i=0: Identity matrix of size npts × npts.
• For i≠0: Matrix with ones on the i-th diagonal, zeros elsewhere.

See also: shift

shift(Ωₕ::AbstractMeshType, ::Val{SHIFT_DIM}, ::Val{i})

Returns the matrix that shifts a grid function by i points along direction SHIFT_DIM, as a sparse operator over the flattened degrees of freedom of Ωₕ.

Arguments

• SHIFT_DIM: the direction to shift along, 1 for $x$, 2 for $y$, 3 for $z$.
• i: how far to shift. 1 is the superdiagonal, -1 the subdiagonal, 0 the identity. The stencil is truncated at the boundary rather than wrapped, so the matrix has n - |i| nonzeros per direction rather than n.

Tensor-Product Structure

A mesh is a tensor product of its one-dimensional meshes, and its degrees of freedom are flattened in column-major order, so a shift along one direction is the identity in every other direction. Writing $E_k$ for the identity of size $n_k$ and $S_k(i)$ for the one-dimensional shift by i, the operator is a Kronecker product with $S$ in one slot:

\[\begin{aligned} \text{along } x: &\quad E_z \otimes E_y \otimes S_x(i) \\ \text{along } y: &\quad E_z \otimes S_y(i) \otimes E_x \\ \text{along } z: &\quad S_z(i) \otimes E_y \otimes E_x \end{aligned}\]

_recursive_shift builds exactly this, recursing from the outermost dimension inwards and placing $S$ when it reaches SHIFT_DIM. The per-direction forms it generalises, each of which the test suite checks against shift:

# 1D
shift(Ωₕ, Val(1), Val(i))  ==  _Eye(be, nₓ, Val(i))

# 2D, on an nₓ × n_y grid
shift(Ωₕ, Val(1), Val(i))  ==  backend_eye(be, n_y) ⊗ _Eye(be, nₓ, Val(i))
shift(Ωₕ, Val(2), Val(i))  ==  _Eye(be, n_y, Val(i)) ⊗ backend_eye(be, nₓ)

# 3D, on an nₓ × n_y × n_z grid
shift(Ωₕ, Val(3), Val(i))  ==  _Eye(be, n_z, Val(i)) ⊗ backend_eye(be, nₓ * n_y)

i == 0 short-circuits to the identity of the whole grid without building any Kronecker product.

This is the building block of the difference, jump and average matrices: a backward difference is shift(Ωₕ, dim, Val(0)) - shift(Ωₕ, dim, Val(-1)), and the other families differ only in which pair of shifts they subtract or average.

See also: ⊗, diff₋ₓ.

_define_directional_alias!(base_op_name, alias_name, dir_string, suffix,
                           direction_index, what, formula; opening_sentence = "")

Defines alias_name(vₕ, uₕ) as base_op_name(vₕ, uₕ, Val(direction_index)) and attaches a docstring to it.

The in-place sibling of _define_directional_alias. Two generators rather than one because the shapes differ: the allocating alias takes a single argument that may be a mesh, a space or a grid function, while this one takes a destination and a source and is only ever about grid functions.

opening_sentence, given non-empty, replaces the generic "The $dir_string $what of uₕ along the $suffix direction, $$formula$, written into vₕ." the same way it does for _define_directional_alias – Dc!/Dₕ! have no backward/forward adjective to put in dir_string either.

_define_directional_alias(base_op_name, alias_name, dir_string, suffix,
                          direction_index, what, formula;
                          opening_sentence = "", formula_note = "", alias_note = "",
                          trailing_note = "")

Defines alias_name(arg) as base_op_name(arg, Val(direction_index)) and attaches a docstring to it.

what names the quantity, such as "finite difference", and formula is the LaTeX for it. Both are needed because the four operator families share this generator: describing every alias as a "difference" would be wrong for the averages, and would not separate the unscaled difference from the finite difference.

opening_sentence, given non-empty, replaces the generic "The $dir_string $what along the $suffix direction, $$formula$." with the caller's own wording: Dc and Dₕ have no backward/forward adjective to put in dir_string at all.

The three remaining keyword notes are each a sentence the docstring includes only when given (non-empty), one per insertion point a family may need: formula_note follows the opening sentence (the diff/finite-difference families use this to contrast the two, which does not apply to an average); alias_note follows the Alias for ... sentence, before arg is described (Dₕ uses this to compare itself with Dc); trailing_note follows the description of arg, before the closing "Accepts a grid function..." paragraph (Dstar₊, Dc and Dₕ use this for their boundary-behaviour and precondition caveats, which differ both in what happens at the ends – Dstar₊/Dc truncate, Dₕ falls back to a one-sided difference (gpena/Bramble.jl#183) – and in whether a mesh needs at least three points along the direction).

_define_grid_function_forms(base_name, apply_fn, extra_args, dir_instance;
                            docstring = "")

Defines the three methods every directional operator family applies a grid function through: base_name! on a scalar grid function, base_name! on a composite one, and the allocating base_name built on top of them.

The three are byte-identical across the families apart from which applicator they call and what it takes before the direction (gpena/Bramble.jl#101) — difference.jl and average.jl each generated them from their own @eval loop before this existed:

extra_args are spliced as bare identifiers, so each generated method names an ordinary top-level function rather than closing over one: that is what keeps every call site specialising to its own zero-allocation method, as the hand-written versions did.

The composite method recurses into the scalar one through apply_fn's own composite method rather than repeating the walk, so a leaf's own submesh is what each leaf is measured against (gpena/Bramble.jl#79).

docstring, given non-empty, is attached to the scalar base_name! method — the families whose prose lives here rather than on a separately hand-written matrix form use it.

_define_operator_aliases(base_name, alias_stem, dir_string, what, formula;
                         vectorial_alias = nothing, vectorial_note = "",
                         alias_kwargs = _no_alias_kwargs,
                         bang_alias_kwargs = _no_alias_kwargs)

Defines one family's whole alias surface: the per-coordinate alias_stem pair for every direction (Dcₓ/Dcₓ!, M₋ᵧ/M₋ᵧ!, …) and, when vectorial_alias is given, the ₕ alias over every coordinate at once.

_define_directional_alias/_define_directional_alias! and _define_vectorial_alias were already shared; the loop calling them was written out once per family in difference.jl and once more in average.jl (gpena/Bramble.jl#101). This is that loop.

alias_kwargs/bang_alias_kwargs are called as f(direction, suffix) and return the keyword arguments for that one alias — the notes each family needs are bespoke prose ("over the averaged spacing", "second order on a non-uniform grid, where Dcₓ is first"), so they are supplied per family rather than templated from what/formula.

vectorial_dir_string/vectorial_what default to dir_string/what and exist for the families that describe the tuple-valued alias differently from the per-coordinate ones: Dstar₊/Dc/Dₕ override every directional opening_sentence and so pass dir_string and what empty, while Dstar₊ₕ/Dcₕ/∇ₕ still want "the centered difference of arg along every coordinate".

_define_vectorial_alias(base_op_name, alias_name, dir_string, what; note = "")

Defines the ₕ alias that applies base_op_name along every coordinate and returns a tuple, one entry per spatial dimension. On a one-dimensional mesh it returns that single entry rather than a one-tuple.

The counterpart of _define_directional_alias for the tuple-valued aliases (∇₋ₕ, diff₋ₕ, M₋ₕ). The operator families generated the same three methods independently before this existed.

note, given non-empty, is an extra sentence appended after the worked 2D example – ∇ₕ uses this to place itself relative to ∇₋ₕ/∇₊ₕ, a comparison none of the other vectorial aliases need.

StarSpacings(h)

Lazy view of the averaged spacings $(h_i + h_{i+1})/2$ over a mesh's cached backward spacings h, which is what Dstar₊ₓ divides by.

Entry i reads h[i] and h[i+1], so it is defined for i < length(h). That is exactly the range the forward stencil's interior covers; the last point has no forward neighbour and the engine truncates it to zero without consulting this.

backward_derivative_weights!(v::AbstractVector, Ωₕ::AbstractMeshType, diff_dim::Val)

Computes the geometric weights for the backward finite difference operator and stores them in-place in vector v.

backward_difference(arg, dim_val::Val)

Constructs the unscaled backward difference operator, representing the operation $u_{i} - u_{i-1}$.

backward_difference_dim!(out, in, [h], dims, diff_dim)

Low-level, in-place function to compute the unscaled backward difference of vector in along dimension diff_dim, storing the result in out. This function computes $u_{i} - u_{i-1}$.

backward_finite_difference(arg, dim_val::Val)

Constructs the backward finite difference operator, which approximates the first derivative using the formula $\frac{u_{i} - u_{i-1}}{h_i}$.

centered_difference(Ωₕ::AbstractMeshType, dim_val::Val)

The centered difference along dim_val, as a sparse matrix.

Reaches one point either side, so both end rows are empty and the mesh needs at least three points along the direction.

cross_weighted_difference(Ωₕ::AbstractMeshType, dim_val::Val)

The cross-weighted centered difference along dim_val, as a sparse matrix.

A three-point stencil in the interior, but neither end row is empty: with no neighbour on the far side, row 1 agrees with D₊ₓ(Ωₕ, dim_val) and row n with D₋ₓ(Ωₕ, dim_val), each under its own diagonal weight alongside the interior cross-weighting. The mesh still needs at least three points along the direction.

forward_derivative_weights!(v::AbstractVector, Ωₕ::AbstractMeshType, diff_dim::Val)

Computes the geometric weights for the forward finite difference operator and stores them in-place in vector v.

forward_difference(arg, dim_val::Val)

Constructs the unscaled forward difference operator, representing the operation $u_{i+1} - u_i$.

forward_difference_dim!(out, in, [h], dims, diff_dim)

Low-level, in-place function to compute the unscaled forward difference of vector in along dimension diff_dim, storing the result in out. This function computes $u_{i+1} - u_i$.

forward_finite_difference(arg, dim_val::Val)

Constructs the forward finite difference operator, which approximates the first derivative using the formula $\frac{u_{i+1} - u_i}{h_i}$.

forward_star_difference(Ωₕ::AbstractMeshType, dim_val::Val)

The starred forward difference along dim_val, as a sparse matrix.

The forward difference scaled by the averaged spacing instead of the forward one. The last point along the direction has no forward neighbour, so its row is empty.

star_spacings(Ωₕ::Mesh1D)

Returns the averaged spacings $(h_i + h_{i+1})/2$ of Ωₕ as a StarSpacings view over its cached backward spacings. Allocates nothing.

Away from the first point this equals half_spacing(Ωₕ, i). At i = 1 it does not: the cached h₁ repeats the first interval, so this gives $x_2 - x_1$ where the cell width gives half of it, the boundary cell being a half cell.

jump!(vₕ, uₕ, dim_val::Val)

The jump across the interfaces along dim_val, written into vₕ, which is returned.

Forwards to forward_difference!, as the allocating form forwards to forward_difference: a jump across an interface and an unscaled forward difference are the same quantity.

jump(arg, dim_val::Val)

The jump across the interfaces along the direction dim_val, $u_{i+1} - u_i$.

arg is a mesh, a grid space or a VectorElement: the first two give the operator as a sparse matrix, the third applies it and returns a VectorElement.

Forwards to forward_difference, since a jump across an interface and an unscaled forward difference are the same quantity. The name records which of the two is meant.

jump_dim!(out, in, dims, jump_dim::Val)

In-place jump of in along jump_dim, written into out, computing $u_{i+1} - u_i$.

Forwards to forward_difference_dim!; the jump and the unscaled forward difference are the same arithmetic.

backward_average(arg, dim_val::Val)

Constructs or applies the backward averaging operator, representing the operation $\frac{u_{i-1} + u_{i}}{2}$.

backward_average_dim!(out, in, dims, average_dim)

Low-level, in-place function to compute the backward average of vector in along dimension average_dim, storing the result in out. This function computes $\frac{u_{i-1} + u_{i}}{2}$.

forward_average(arg, dim_val::Val)

Constructs or applies the forward averaging operator, representing the operation $\frac{u_{i} + u_{i+1}}{2}$.

forward_average_dim!(out, in, dims, average_dim)

Low-level, in-place function to compute the forward average of vector in along dimension average_dim, storing the result in out. This function computes $\frac{u_{i} + u_{i+1}}{2}$.

inner_Γ(uₕ::VectorElement, vₕ::VectorElement, labels::Symbol...)

Placeholder for the true, $(D-1)$-dimensional boundary integral $\int_\Gamma u\,v\,ds$ over the mesh regions labels name.

Not the same quantity as innerₕ(uₕ, vₕ; markers = labels). That is a masked sum of the existing cell measures (a $D$-dimensional quantity restricted to a set of points), and it scales like h, vanishing under refinement: 0.125 on a 5×5 mesh of the unit square restricted to :bottom. This function is meant for the mesh-independent surface integral a Neumann or Robin term needs (1.0 on that same mesh and region), which is a genuinely different quantity, not a bug in the other one. See point 11 of docs/form-unlock-plan.md.

Not yet implemented. No $(D-1)$-dimensional surface quadrature weight exists anywhere in the package today: only Innerh (cell measures) and Innerplus (directional, still $D$-dimensional). This always throws until that weight is built; kept unexported until it does something.