terra::fv::hex::operators Namespace Reference

Namespaces
namespace	detail
namespace	fct_detail

Classes
struct	ComputeDtStableKernel
	Kokkos kernel that computes the local maximum stable explicit dt for each FV cell. More...
struct	FCTAntidiffKernel
	Kokkos kernel that computes only the pre-scaled antidiffusive fluxes from \(T^n\). More...
struct	FCTCorrectionKernel
	Kokkos kernel that applies the Zalesak-limited antidiffusive correction to \(T^L\). More...
struct	FCTLimiterKernel
	Kokkos kernel that computes the Zalesak nodal correction factors \(R_i^+\) and \(R_i^-\) from the low-order predictor \(T^L\) and the antidiffusive fluxes. More...
struct	FCTPredictorKernel
	Kokkos kernel for the explicit low-order predictor and antidiffusive flux assembly. More...
struct	FVFCTBuffers
	All working storage for a single FCT timestep, allocated once and reused every step. More...
class	Mass
class	UnsteadyAdvectionDiffusion

Functions
template<typename ScalarT >
ScalarT	compute_dt_stable (const grid::shell::DistributedDomain &domain, const linalg::VectorQ1Vec< ScalarT, 3 > &vel, const grid::Grid4DDataVec< ScalarT, 3 > &cell_centers, const grid::Grid3DDataVec< ScalarT, 3 > &grid, const grid::Grid2DDataScalar< ScalarT > &radii, const ScalarT diffusivity=ScalarT(0))
	Compute the largest explicit time step that keeps the FCT low-order predictor stable.
template<typename ScalarT >
void	fct_predictor (const grid::shell::DistributedDomain &domain, const linalg::VectorFVScalar< ScalarT > &T_old, const linalg::VectorQ1Vec< ScalarT, 3 > &vel, const grid::Grid4DDataVec< ScalarT, 3 > &cell_centers, const grid::Grid3DDataVec< ScalarT, 3 > &grid, const grid::Grid2DDataScalar< ScalarT > &radii, const ScalarT dt, FVFCTBuffers< ScalarT > &bufs, const ScalarT diffusivity=ScalarT(0), const grid::Grid4DDataScalar< ScalarT > &source={}, const bool subtract_divergence=true)
	Stage 1: low-order predictor + antidiffusive flux computation.
template<typename ScalarT >
void	fct_limiter (const grid::shell::DistributedDomain &domain, FVFCTBuffers< ScalarT > &bufs)
	Stage 2: compute Zalesak \(R^+\)/ \(R^-\) correction factors.
template<typename ScalarT >
void	fct_correction (const grid::shell::DistributedDomain &domain, linalg::VectorFVScalar< ScalarT > &T_new, FVFCTBuffers< ScalarT > &bufs)
	Stage 3: apply the Zalesak-limited antidiffusive correction.
template<typename ScalarT >
void	fct_explicit_step (const grid::shell::DistributedDomain &domain, linalg::VectorFVScalar< ScalarT > &T, const linalg::VectorQ1Vec< ScalarT, 3 > &vel, const grid::Grid4DDataVec< ScalarT, 3 > &cell_centers, const grid::Grid3DDataVec< ScalarT, 3 > &grid, const grid::Grid2DDataScalar< ScalarT > &radii, const ScalarT dt, FVFCTBuffers< ScalarT > &bufs, const ScalarT diffusivity=ScalarT(0), const grid::Grid4DDataScalar< ScalarT > &source={}, const bool subtract_divergence=true, const grid::Grid4DDataScalar< grid::shell::ShellBoundaryFlag > &boundary_mask={}, const DirichletBCs< ScalarT > &bcs={})
	One complete explicit FCT advection–diffusion timestep.
template<typename ScalarT >
void	upwind_explicit_step (const grid::shell::DistributedDomain &domain, linalg::VectorFVScalar< ScalarT > &T, const linalg::VectorQ1Vec< ScalarT, 3 > &vel, const grid::Grid4DDataVec< ScalarT, 3 > &cell_centers, const grid::Grid3DDataVec< ScalarT, 3 > &grid, const grid::Grid2DDataScalar< ScalarT > &radii, const ScalarT dt, FVFCTBuffers< ScalarT > &bufs, const ScalarT diffusivity=ScalarT(0), const grid::Grid4DDataScalar< ScalarT > &source={}, const bool subtract_divergence=true)
	One explicit first-order upwind advection–diffusion timestep (no FCT correction).
template<typename ScalarT >
void	fct_antidiff (const grid::shell::DistributedDomain &domain, const linalg::VectorFVScalar< ScalarT > &T_old, const linalg::VectorQ1Vec< ScalarT, 3 > &vel, const grid::Grid4DDataVec< ScalarT, 3 > &cell_centers, const grid::Grid3DDataVec< ScalarT, 3 > &grid, const grid::Grid2DDataScalar< ScalarT > &radii, const ScalarT dt, FVFCTBuffers< ScalarT > &bufs)
	Computes pre-scaled antidiffusive fluxes \(\tilde{f}_{ij}\) from \(T^n\) (semi-implicit FCT stage 1).
template<typename ScalarT >
void	fct_semiimplicit_step (const grid::shell::DistributedDomain &domain, linalg::VectorFVScalar< ScalarT > &T, const linalg::VectorFVScalar< ScalarT > &T_L, const linalg::VectorQ1Vec< ScalarT, 3 > &vel, const grid::Grid4DDataVec< ScalarT, 3 > &cell_centers, const grid::Grid3DDataVec< ScalarT, 3 > &grid, const grid::Grid2DDataScalar< ScalarT > &radii, const ScalarT dt, FVFCTBuffers< ScalarT > &bufs)
	One semi-implicit FCT advection–diffusion timestep.

Function Documentation

compute_dt_stable()

template<typename ScalarT >

ScalarT terra::fv::hex::operators::compute_dt_stable	(	const grid::shell::DistributedDomain &	domain,
		const linalg::VectorQ1Vec< ScalarT, 3 > &	vel,
		const grid::Grid4DDataVec< ScalarT, 3 > &	cell_centers,
		const grid::Grid3DDataVec< ScalarT, 3 > &	grid,
		const grid::Grid2DDataScalar< ScalarT > &	radii,
		const ScalarT	diffusivity = ScalarT( 0 )
	)

Compute the largest explicit time step that keeps the FCT low-order predictor stable.

Performs a Kokkos parallel reduction over all non-ghost FV cells followed by an MPI_Allreduce to obtain the global minimum across all MPI ranks.

The result is exact (derived from the actual face-normal velocity fluxes and cell volumes), unlike the approximate estimate \(h_\text{min,radial} / u_\text{max}\) which ignores lateral cell sizes and diffusion stiffness on non-orthogonal cells.

Note

The resulting time step size can still be too large to ensure an "accurate" time discretization. It can be a good idea to scale it down further in practice.

Typical usage:

const ScalarType dt = dt_scaling * fv::hex::operators::compute_dt_stable(
    domain, u, cell_centers.grid_data(), coords_shell, coords_radii, diffusivity);

Parameters

domain	Distributed domain.
vel	Q1 nodal velocity (read-only; no ghost-layer update required).
cell_centers	Pre-computed cell centres with ghost layers (from initialize_cell_centers).
grid	Lateral node coordinates of the unit-sphere surface.
radii	Radial shell radii.
diffusivity	Physical diffusivity \(\kappa\) (default 0 = pure advection).

Returns

The minimum over all cells of \(M_{ii}/\lambda_i\).

fct_antidiff()

template<typename ScalarT >

void terra::fv::hex::operators::fct_antidiff	(	const grid::shell::DistributedDomain &	domain,
		const linalg::VectorFVScalar< ScalarT > &	T_old,
		const linalg::VectorQ1Vec< ScalarT, 3 > &	vel,
		const grid::Grid4DDataVec< ScalarT, 3 > &	cell_centers,
		const grid::Grid3DDataVec< ScalarT, 3 > &	grid,
		const grid::Grid2DDataScalar< ScalarT > &	radii,
		const ScalarT	dt,
		FVFCTBuffers< ScalarT > &	bufs
	)

Computes pre-scaled antidiffusive fluxes \(\tilde{f}_{ij}\) from \(T^n\) (semi-implicit FCT stage 1).

Ghost layers of T_old are exchanged via MPI before the kernel runs. On exit bufs.antidiff holds \(\tilde{f}_{ij} = (\Delta t / M_{ii})\,f_{ij}\) for every real cell and all 6 faces.

Parameters

domain	Distributed domain (used for MPI ghost-layer exchange).
T_old	Transported scalar \(T^n\) at time level \(n\).
vel	Q1 nodal velocity field \(\mathbf{u}\).
cell_centers	Pre-computed cell centres (ghost layers populated via initialize_cell_centers).
grid	Lateral node coordinates of the unit-sphere surface.
radii	Radial shell radii.
dt	Time step \(\Delta t\).
bufs	FCT scratch arrays; only antidiff and ghost_T are written.

fct_correction()

template<typename ScalarT >

void terra::fv::hex::operators::fct_correction	(	const grid::shell::DistributedDomain &	domain,
		linalg::VectorFVScalar< ScalarT > &	T_new,
		FVFCTBuffers< ScalarT > &	bufs
	)

Stage 3: apply the Zalesak-limited antidiffusive correction.

Reads bufs.T_L, bufs.antidiff, bufs.R_plus, bufs.R_minus (all read-only). Writes \(T^{n+1} = T^L + \sum_j \alpha_{ij}\,\tilde{f}_{ij}\) into T_new.

Precondition

fct_limiter must have been called first so that ghost layers of \(R^\pm\) are current.

Parameters

domain	Distributed domain.
T_new	Output: high-order corrected scalar at \(t_{n+1}\).
bufs	FCT scratch arrays (reads T_L, antidiff, R_plus, R_minus).

fct_explicit_step()

template<typename ScalarT >

void terra::fv::hex::operators::fct_explicit_step	(	const grid::shell::DistributedDomain &	domain,
		linalg::VectorFVScalar< ScalarT > &	T,
		const linalg::VectorQ1Vec< ScalarT, 3 > &	vel,
		const grid::Grid4DDataVec< ScalarT, 3 > &	cell_centers,
		const grid::Grid3DDataVec< ScalarT, 3 > &	grid,
		const grid::Grid2DDataScalar< ScalarT > &	radii,
		const ScalarT	dt,
		FVFCTBuffers< ScalarT > &	bufs,
		const ScalarT	diffusivity = ScalarT( 0 ),
		const grid::Grid4DDataScalar< ScalarT > &	source = {},
		const bool	subtract_divergence = true,
		const grid::Grid4DDataScalar< grid::shell::ShellBoundaryFlag > &	boundary_mask = {},
		const DirichletBCs< ScalarT > &	bcs = {}
	)

One complete explicit FCT advection–diffusion timestep.

Executes the three stages in sequence:

1. fct_predictor — low-order upwind + diffusion predictor \(T^L\) and antidiff fluxes.
2. fct_limiter — Zalesak \(R^\pm\) factors.
3. fct_correction — limited antidiffusive correction \(T^{n+1} = T^L + \sum_j \alpha_{ij}\,\tilde{f}_{ij}\).

The scheme is LED (local extremum diminishing) and conservation-consistent: the limited fluxes are antisymmetric, so cell integrals are preserved up to boundary terms.

For a non-linear iteration (e.g. defect-correction), call the three stages individually and repeat stages 2–3 with an updated \(T^L\).

Stability: requires \(\mathrm{CFL} = \Delta t\,\max_i(\sum_j \beta_{ij}^+)/M_{ii} < 1\).

Parameters

domain	Distributed domain.
T	Transported scalar — \(T^n\) on entry, \(T^{n+1}\) on exit.
vel	Q1 nodal velocity field \(\mathbf{u}\).
cell_centers	Pre-computed cell centres (ghost layers populated via initialize_cell_centers).
grid	Lateral node coordinates of the unit-sphere surface.
radii	Radial shell radii.
dt	Time step \(\Delta t\).
bufs	Pre-allocated FCT scratch arrays.
diffusivity	Physical diffusivity \(\kappa \geq 0\) (default 0 = pure advection).
source	Volumetric source term \(f\) [T/time]; null view = no source.
subtract_divergence	Subtract discrete divergence error (default true); see fct_predictor.
boundary_mask	Node-based boundary flag array (Q1 layout). When provided (non-null), Dirichlet BCs are enforced on \(T^L\) between the predictor and limiter so the Zalesak \(R^\pm\) factors see the correct boundary values. Default: null (no enforcement).
bcs	Prescribed boundary values. Ignored when boundary_mask is null.

fct_limiter()

template<typename ScalarT >

void terra::fv::hex::operators::fct_limiter	(	const grid::shell::DistributedDomain &	domain,
		FVFCTBuffers< ScalarT > &	bufs
	)

Stage 2: compute Zalesak \(R^+\)/ \(R^-\) correction factors.

Ghost layers of bufs.T_L are exchanged via MPI so that the neighbourhood min/max stencil is correct for subdomain-boundary cells. Then FCTLimiterKernel is launched, and finally ghost layers of bufs.R_plus and bufs.R_minus are exchanged so that Stage 3 can read neighbour factors.

This function can be called multiple times in a non-linear iteration loop with a fixed bufs.antidiff computed once from \(T^n\) but an updated bufs.T_L iterate.

Parameters

domain	Distributed domain (used for ghost-layer exchange).
bufs	FCT scratch arrays; reads T_L and antidiff, writes R_plus and R_minus.

fct_predictor()

template<typename ScalarT >

void terra::fv::hex::operators::fct_predictor	(	const grid::shell::DistributedDomain &	domain,
		const linalg::VectorFVScalar< ScalarT > &	T_old,
		const linalg::VectorQ1Vec< ScalarT, 3 > &	vel,
		const grid::Grid4DDataVec< ScalarT, 3 > &	cell_centers,
		const grid::Grid3DDataVec< ScalarT, 3 > &	grid,
		const grid::Grid2DDataScalar< ScalarT > &	radii,
		const ScalarT	dt,
		FVFCTBuffers< ScalarT > &	bufs,
		const ScalarT	diffusivity = ScalarT( 0 ),
		const grid::Grid4DDataScalar< ScalarT > &	source = {},
		const bool	subtract_divergence = true
	)

Stage 1: low-order predictor + antidiffusive flux computation.

Ghost layers of T_old are exchanged via MPI before the kernel runs. On exit:

• bufs.T_L holds \(T^L\) (first-order upwind + explicit diffusion at \(t+\Delta t\)).
• bufs.antidiff holds \(\tilde{f}_{ij} = (\Delta t / M_{ii})\,f_{ij}\) for every face.

Parameters

domain	Distributed domain (used for MPI ghost-layer exchange).
T_old	Transported scalar \(T^n\) at time level \(n\).
vel	Q1 velocity field \(\mathbf{u}\) (nodal, read-only).
cell_centers	Pre-computed cell centres (with ghost layers filled via initialize_cell_centers).
grid	Lateral node coordinates of the unit-sphere surface.
radii	Radial shell radii.
dt	Time step \(\Delta t\).
bufs	Pre-allocated FCT scratch arrays (updated in-place).
diffusivity	Physical diffusivity \(\kappa\) (default 0 = pure advection).
source	Volumetric source term \(f\) [T/time] per cell. A default-constructed (null) view means no source. Added as \(\Delta t \cdot f_i\) to the low-order predictor.
subtract_divergence	When true (default), subtract \(T_i\,(\sum_j \beta_{ij})/M_{ii}\) from the predictor, converting \(\nabla\cdot(\mathbf{u}T)\) to \(\mathbf{u}\cdot\nabla T\). This is always correct for temperature: when \(\nabla\cdot\mathbf{u}=0\) the correction vanishes; when \(\nabla\cdot\mathbf{u}\neq 0\) it removes an unphysical source term.

fct_semiimplicit_step()

template<typename ScalarT >

void terra::fv::hex::operators::fct_semiimplicit_step	(	const grid::shell::DistributedDomain &	domain,
		linalg::VectorFVScalar< ScalarT > &	T,
		const linalg::VectorFVScalar< ScalarT > &	T_L,
		const linalg::VectorQ1Vec< ScalarT, 3 > &	vel,
		const grid::Grid4DDataVec< ScalarT, 3 > &	cell_centers,
		const grid::Grid3DDataVec< ScalarT, 3 > &	grid,
		const grid::Grid2DDataScalar< ScalarT > &	radii,
		const ScalarT	dt,
		FVFCTBuffers< ScalarT > &	bufs
	)

One semi-implicit FCT advection–diffusion timestep.

The semi-implicit scheme decouples the implicit solve from the flux-correction step, removing the CFL restriction while retaining the LED property. The caller is responsible for providing the implicit low-order predictor \(T^L\), obtained by solving the backward-Euler system

\[ \bigl(M + \Delta t\,A_{\mathrm{upw}} + \Delta t\,D\bigr)\,T^L = M\,T^n, \]

where \(A_{\mathrm{upw}}\) is the upwind advection matrix and \(D\) is the two-point diffusion matrix. See UnsteadyAdvectionDiffusion in advection_diffusion.hpp for the concrete operator and FGMRES for the solver.

Steps executed internally:

1. fct_antidiff: ghost-update \(T^n\), compute \(\tilde{f}_{ij}\) from \(T^n\).
2. Copy the externally provided \(T^L\) into bufs.T_L.
3. fct_limiter: ghost-update \(T^L\), compute \(R^\pm\), ghost-update \(R^\pm\).
4. fct_correction: \(T^{n+1} = T^L + \sum_j \alpha_{ij}\,\tilde{f}_{ij}\) → written into \(T\).

Stability: unconditionally stable for the low-order predictor; the Zalesak correction cannot introduce new extrema beyond those already present in \(T^L\).

Parameters

domain	Distributed domain.
T	Transported scalar — \(T^n\) on entry, \(T^{n+1}\) on exit.
T_L	Implicit low-order predictor satisfying the backward-Euler system above.
vel	Q1 nodal velocity field \(\mathbf{u}\) (needed only for antidiff geometry).
cell_centers	Pre-computed cell centres (ghost layers populated via initialize_cell_centers).
grid	Lateral node coordinates of the unit-sphere surface.
radii	Radial shell radii.
dt	Time step \(\Delta t\).
bufs	Pre-allocated FCT scratch arrays.

upwind_explicit_step()

template<typename ScalarT >

void terra::fv::hex::operators::upwind_explicit_step	(	const grid::shell::DistributedDomain &	domain,
		linalg::VectorFVScalar< ScalarT > &	T,
		const linalg::VectorQ1Vec< ScalarT, 3 > &	vel,
		const grid::Grid4DDataVec< ScalarT, 3 > &	cell_centers,
		const grid::Grid3DDataVec< ScalarT, 3 > &	grid,
		const grid::Grid2DDataScalar< ScalarT > &	radii,
		const ScalarT	dt,
		FVFCTBuffers< ScalarT > &	bufs,
		const ScalarT	diffusivity = ScalarT( 0 ),
		const grid::Grid4DDataScalar< ScalarT > &	source = {},
		const bool	subtract_divergence = true
	)

One explicit first-order upwind advection–diffusion timestep (no FCT correction).

Equivalent to calling fct_predictor and then copying \(T^L \to T\). The result is the maximally diffusive (but stable) low-order solution. Use this as a baseline or whenever sharp-gradient preservation is not required.

Physical diffusion (two-point flux) is included when diffusivity > 0.

Stability: requires \(\mathrm{CFL} = \Delta t\,\max_i(\sum_j \beta_{ij}^+)/M_{ii} < 1\).

Parameters

domain	Distributed domain.
T	Transported scalar — \(T^n\) on entry, \(T^{n+1}\) on exit.
vel	Q1 nodal velocity field \(\mathbf{u}\).
cell_centers	Pre-computed cell centres (ghost layers populated via initialize_cell_centers).
grid	Lateral node coordinates of the unit-sphere surface.
radii	Radial shell radii.
dt	Time step \(\Delta t\).
bufs	Pre-allocated FCT scratch arrays (uses T_L and ghost_T).
diffusivity	Physical diffusivity \(\kappa \geq 0\) (default 0 = pure advection).
source	Optional volumetric source term \(f\) [T/time] (default: none).
subtract_divergence	Subtract discrete divergence error (default true); see fct_predictor.