Limiter

cell_mean(u::AbstractLazy, dΩ::Measure)
cell_mean(u::MultiFieldFEFunction, dΩ::Measure)
cell_mean(u::AbstractLazy, cache::CellMeanCache)
cell_mean(u::MultiFieldFEFunction, cache::Tuple{Vararg{CellMeanCache}})
cell_mean(u::MeshData{CellData}, dΩ::Measure)

Compute the mean value of u in each cell covered by the cell Measure dΩ, i.e. the L2-projection of u onto the piecewise-constant (degree 0, discontinuous) space:

\[ ūᵢ = \frac{1}{|Kᵢ|} ∫_{Kᵢ}(u) dΩ\]

The integrals are evaluated with the quadrature rule of dΩ, which must therefore be accurate enough — as for projection_l2!.

When the mean values are computed repeatedly — e.g. at each time step of an explicit scheme — the underlying mass matrix can be assembled and factorized only once with build_cell_mean_cache, the mean values being then obtained with cell_mean(u, cache).

Arguments

• u: the field whose mean values are sought: a SingleFieldFEFunction or, more generally, any lazy expression (AbstractLazy, e.g. a PhysicalFunction); it can also be a MultiFieldFEFunction (the mean values are then computed field by field), or even a MeshData{CellData}, in which case it is returned unchanged.
• dΩ: a Measure built on a CellDomain (e.g. Measure(CellDomain(mesh), 2d+1)).
• cache: a cache previously built with build_cell_mean_cache, allowing to reuse the assembled and factorized mass matrix instead of rebuilding it at each call.

Returns

A MeshCellData (i.e. a MeshData{<:CellData}) containing one mean value per cell — a scalar for a single-component field, an SVector for a multi-component field. For a MultiFieldFEFunction, a tuple of MeshCellData (one per field) is returned.

Examples

mesh = rectangle_mesh(20, 4)
degree = 2
fes = TrialFESpace(FunctionSpace(:Lagrange, degree), mesh, :discontinuous)
u = FEFunction(fes, mesh, PhysicalFunction(x -> x[1]))

dΩ = Measure(CellDomain(mesh), 2 * degree + 1)

u_mean = cell_mean(u, dΩ) # one mean value per cell

# repeated evaluations (e.g. at each time step): reuse a prebuilt cache
cache = build_cell_mean_cache(u, dΩ)
u_mean = cell_mean(u, cache)

Devs notes

• This function should be moved to a more appropriate file (TODO)

build_cell_mean_cache(u::AbstractLazy, dΩ::Measure)
build_cell_mean_cache(u::MultiFieldFEFunction, dΩ::Measure)

Build the cache required to compute the cell mean values of u on the cells covered by the Measure dΩ: the piecewise-constant discontinuous FESpace in which the mean values live, together with its factorized mass matrix.

This is useful when the mean values are computed repeatedly — e.g. at each time step of an explicit scheme — as the mass matrix is then assembled and factorized only once:

cache = build_cell_mean_cache(u, dΩ)
u_mean = cell_mean(u, cache) # cheap: reuses the cache

For a MultiFieldFEFunction, a tuple of caches (one per field) is returned, to be passed to the corresponding cell_mean method. Note that cell_mean(u, dΩ) simply builds the cache and calls cell_mean(u, cache).

linear_scaling_limiter(
    u::SingleFieldFEFunction,
    u_mean::MeshData{<:CellData},
    domainOrMeasure::Union{AbstractCellDomain, AbstractMeasure{<:AbstractCellDomain}};
    targetMeasures::NTuple{N, AbstractMeasure} = default_target_measures(u, domainOrMeasure),
    periodicBCs::Union{NTuple{Nbc, BoundaryFaceDomain{M, <:PeriodicBCType}}, Nothing} = nothing,
    bounds::Union{Tuple{<:Number, <:Number}, Nothing} = nothing,
    DMPrelax = zero(get_dof_type(u)),
    coefmax = one(get_dof_type(u)),
    checkvalues = false,
)

Apply the linear scaling limiter (see "Maximum-principle-satisfying and positivity-preserving high-order schemes for conservation laws: survey and new developments", Zhang & Shu) to a scalar discontinuous FEFunction, in order to enforce a local discrete maximum principle (DMP), i.e. to prevent the creation of spurious extrema — typically in the neighborhood of discontinuities.

In each cell, the limited solution is obtained by scaling the fluctuation around the cell mean with a scalar coefficient lim_u ∈ [0, coefmax]:

u_limited = u_mean + lim_u * (u - u_mean)

The cell mean is thus preserved by construction, and lim_u = 0 corresponds to a first-order (piecewise constant) solution.

Arguments

• u: the scalar FEFunction to limit; its FESpace must be discontinuous.
• u_mean: the cell mean values of u, as a MeshCellData (i.e. MeshData{<:CellData}). It is typically obtained with cell_mean(u, dΩ) where dΩ is a Measure on the cells of the mesh.
• domainOrMeasure: either the AbstractCellDomain — typically the whole CellDomain(mesh) — or the cell Measure dΩ used to compute u_mean. It is only used to retrieve the mesh and the neighboring cells, needed to compute the min/max of the neighboring cell mean values. Restricted CellDomains (a subset of cells) are not supported yet: the neighbor search currently runs over all mesh faces.
• targetMeasures: the list of Measures on which the min and max of u are evaluated (at the quadrature nodes of each measure). By default it is built by default_target_measures as (dΩ, dΓ, dΓbc): the cell measure, a measure on the interior faces and — when the mesh defines boundary names — one on the boundary faces. It can also be provided explicitly, as in (dΩ, dΓ). Note that a Measure defined on a periodic BoundaryFaceDomain is not supported yet (it currently raises an error).

Keyword arguments

• periodicBCs: not yet supported: any value other than nothing currently raises an error.
• bounds: optional (lower, upper) strong physical bounds imposed on the solution. The bounds are clipped to the range allowed by the local cell mean: if a cell mean does not satisfy the bounds, the admissible window of that cell reduces to its mean and the limiter degenerates to a first-order (piecewise constant) solution in that cell.
• DMPrelax: non-negative relaxation parameter (possibly broadcastable to the number of components) that enlarges the admissible range of neighboring cell mean values by - DMPrelax / + DMPrelax (default: 0, i.e. strict DMP).
• coefmax: upper bound of the limiter coefficient (default: 1, i.e. the fluctuation is never amplified).
• checkvalues: if true, check that each cell mean lies within the min/max values of u in the cell and within the admissible window, and error otherwise (default: false). This is a safety net — e.g. against NaNs — as both ranges are initialized from the cell means themselves and are thus consistent by construction.

Returns

A tuple (lim_u, u_lim) where:

• lim_u is a MeshCellData containing one limiter coefficient in [0, coefmax] per cell;
• u_lim is the limited field u_mean + lim_u * (u - u_mean) given as a lazy expression (AbstractLazy), and not as an FEFunction. It can be used directly inside a weak form (e.g. ∫(u_lim ⋅ v)dΩ); to obtain the dofs of the limited solution, project it for instance with projection_l2!.

Examples

mesh = rectangle_mesh(20, 4)
degree = 2
fes = TrialFESpace(FunctionSpace(:Lagrange, degree), mesh, :discontinuous)
u = FEFunction(fes, mesh, PhysicalFunction(x -> x[1])) # the DG field to limit

Ω = CellDomain(mesh)
dΩ = Measure(Ω, 2 * degree + 1)

u_mean = cell_mean(u, dΩ)
lim_u, u_lim = linear_scaling_limiter(u, u_mean, Ω; bounds = (0.0, 1.0))

u_limited = FEFunction(fes)
projection_l2!(u_limited, u_lim, dΩ)

default_target_measures(u::SingleFieldFEFunction, Ω::AbstractCellDomain)
default_target_measures(u::SingleFieldFEFunction, dΩ::AbstractMeasure{<:AbstractCellDomain})

Build the default list of target Measures used by linear_scaling_limiter to evaluate the min and max of the solution.

When a cell domain Ω is provided, a cell measure dΩ is first built on that domain and the method below applies. The quadrature degree of dΩ is the degree of the FunctionSpace of u; for a :Lagrange space the quadrature type follows the space (:Legendre, :Lobatto or :Uniform), and for any other space a uniform quadrature is used. Beware that uniform quadratures are currently only defined on hypercubes (Line, Square, Cube, any degree) and on Triangle (degree 1): any other combination — e.g. a Triangle of degree ≥ 2 or a Tetra — is not supported yet and raises an error when the measures are used.

When a cell measure dΩ is provided, the returned measures are:

• dΩ itself, on which the min/max of the solution are evaluated inside the cells;
• a measure on the interior faces of the mesh, on which they are evaluated on both sides of each face;
• a measure on the boundary faces of the mesh (all boundary names), if the mesh defines boundary names.

The face measures are built with the same quadrature as dΩ. They are required to "see" the values reached on the cell boundaries: for a Gauss-type quadrature (e.g. :Legendre) the nodes of dΩ lie strictly inside the cells.