Procedures

Resolve the communicator a new decomposition will run on and register it with the process-wide lifecycle bookkeeping.

Perform one adaptation epoch on the model (see the module documentation). On return, adapted reports whether the mesh changed; when it did, the model is already rebound to the new mesh with the solution transferred (conservatively), and the caller should re-evaluate its time step (RecommendedTimeStep) before the next ForwardStep. When the leaf set is unchanged the model is untouched.

Perform one adaptation epoch on the model (see the module documentation). On return, adapted reports whether the mesh changed; when it did, the model is already rebound to the new mesh with the solution transferred (conservatively), and the caller should re-evaluate its time step (RecommendedTimeStep) before the next ForwardStep. When the leaf set is unchanged the model is untouched.

Mutate the forest from a per-leaf flag array (indexed over the current leaves in this%leaf order, e.g. the flag array produced by the refinement indicator):

Mutate the forest from a per-leaf flag array (indexed over the current leaves in this%leaf order, e.g. the flag array produced by SELF_RefinementIndicator_2D):

Register boundary conditions for the advection-diffusion model. Hyperbolic: mirror (no-normal-flow) wall condition on the solution. Parabolic: no-stress condition that zeros the normal gradient at the boundary, giving zero diffusive flux through walls. This is stable with the Bassi-Rebay (BR1) method.

Register boundary conditions for the advection-diffusion model. Hyperbolic: mirror (no-normal-flow) wall condition on the solution. Parabolic: no-stress condition that zeros the normal gradient at the boundary, giving zero diffusive flux through walls. This is stable with the Bassi-Rebay (BR1) method.

Register the (CPU) no-normal-flow and radiation boundary conditions. GPU builds call this parent and then overwrite both registrations with the device kernels in AdditionalInit_LinearEuler2D.

Register the (CPU) radiation boundary condition. GPU builds overwrite this registration with the device kernel in AdditionalInit_LinearEuler3D.

Adds a pressure-balanced warm bubble perturbation.

Adds a pressure-balanced warm bubble perturbation.

Allgather an integer array with perElem entries per element from the decomposition's contiguous rank-local element ranges into the global element ordering.

Allgather an integer array with perElem entries per element from the decomposition's contiguous rank-local element ranges into the global element ordering.

Allgather a real(prec) array with perElem entries per element from the decomposition's contiguous rank-local element ranges into the global element ordering.

Allgather a real(prec) array with perElem entries per element from the decomposition's contiguous rank-local element ranges into the global element ordering.

\addtogroup SELF_SupportRoutines @{ \fn AlmostEqual Compares two floating point numbers and determines if they are equal (to machine precision).

Push the controller-owned amplitude-gate settings onto the indicator. Called after every indicator Init, because Init is intent(out) and therefore resets them to the indicator's defaults. The setters carry the validation.

Push the controller-owned amplitude-gate settings onto the indicator. Called after every indicator Init, because Init is intent(out) and therefore resets them to the indicator's defaults. The setters carry the validation.

Execute a transfer plan on nodal element data: uNew(:,:,li,:) receives old element data copied, prolonged (exact interpolation), or restricted (conservative L2 projection) according to plan entry li. interp must be the solution interpolant the data lives on (its mortar operators drive the transfer). Runs once per adaptation epoch - not a per-time-step hot path - so clarity is preferred over fused loops here.

Execute a transfer plan on nodal element data: uNew(:,:,:,li,:) receives old element data copied, prolonged (exact interpolation), or restricted (conservative L2 projection) according to plan entry li. interp must be the solution interpolant the data lives on (its mortar operators drive the transfer). Runs once per adaptation epoch - not a per-time-step hot path - so clarity is preferred over fused loops here.

Apply the transfer plan on the device, writing solution%interior_gpu directly and moving no solution data across the PCIe/xGMI link - on any number of ranks.

Transfer the staged pre-regrid solution onto the regridded mesh through plan, filling the rank-local element range [eFirst,eLast] of the new solution.

Apply the transfer plan on the device, writing solution%interior_gpu directly and moving no solution data across the PCIe/xGMI link - on any number of ranks.

Transfer the staged pre-regrid solution onto the regridded mesh through plan, filling the rank-local element range [eFirst,eLast] of the new solution.

Execute the contiguous sub-range eFirst..eLast of a transfer plan: uNew(:,:,k,:) receives the data of new element eFirst+k-1. uOld is the full (global) old field; the output is only the requested slice. This is the gather-then-slice (AMR Stage-5 v1) decomposed-mesh entry point: each rank passes its own contiguous range of the new element ordering and a gathered global old solution, and fills exactly its rank-local storage. ApplyTransferPlanWindow is the point-to-point (v2) entry point, which needs only the sub-range of the old field the requested new range actually references; this routine is that one over the whole old field, so the two are bit-identical by construction rather than by inspection.

Execute the contiguous sub-range eFirst..eLast of a transfer plan: uNew(:,:,:,k,:) receives the data of new element eFirst+k-1. uOld is the full (global) old field; the output is only the requested slice. This is the gather-then-slice (AMR Stage-5 v1) decomposed-mesh entry point: each rank passes its own contiguous range of the new element ordering and a gathered global old solution, and fills exactly its rank-local storage. ApplyTransferPlanWindow is the point-to-point (v2) entry point, which needs only the sub-range of the old field the requested new range actually references; this routine is that one over the whole old field, so the two are bit-identical by construction rather than by inspection.

Execute the contiguous sub-range eFirst..eLast of a transfer plan against a WINDOW of the old field: uOld holds old elements oldFirst..oldLast only, indexed in the global old element numbering that plan%sourceElem and plan%family use, rather than the whole global field. uNew(:,:,k,:) receives the data of new element eFirst+k-1.

Execute the contiguous sub-range eFirst..eLast of a transfer plan against a WINDOW of the old field: uOld holds old elements oldFirst..oldLast only, indexed in the global old element numbering that plan%sourceElem and plan%family use, rather than the whole global field. uNew(:,:,:,k,:) receives the data of new element eFirst+k-1.

Enforce the 2:1 face balance condition: no leaf face may separate elements differing by more than one refinement level. Iterates to a fixed point - in each sweep, any leaf whose equal-or-larger neighbour is a leaf two or more levels coarser triggers refinement of that coarser neighbour; refinement can ripple, so sweeps repeat until nothing changes. The leaf set is rebuilt on return.

Enforce the 2:1 balance condition: no leaf face may separate elements differing by more than one refinement level. Iterates to a fixed point - in each sweep, any leaf whose equal-or-larger neighbour is a leaf two or more levels coarser triggers refinement of that coarser neighbour; refinement can ripple, so sweeps repeat until nothing changes. The leaf set is rebuilt on return.

Evaluates the bilinear map of face s of a trilinear (nGeo=1) element at the face coordinates (u,v) in [-1,1]^2, using the face trace-coordinate convention of BoundaryInterp: faces 1,6 -> (xi1,xi2); faces 2,4 -> (xi1,xi3); faces 3,5 -> (xi2,xi3).

LLF Riemann flux on GPU — fully device-resident.

LLF Riemann flux on GPU — fully device-resident.

LMARS interface flux on GPU.

LMARS (Low-Mach Approximate Riemann Solver, Chen et al. 2013) interface flux. No hydrostatic pressure split: gravity is folded into SourceMethod via the Souza non-conservative form using the geopotential carried in the state vector (variable index 5).

LMARS interface flux on GPU. No hydrostatic pressure split: gravity is handled by the Souza non-conservative source term using the geopotential carried as state variable index 6.

LMARS (Low-Mach Approximate Riemann Solver, Chen et al. 2013) interface flux. No hydrostatic pressure split: gravity is folded into SourceMethod via the Souza non-conservative form using the geopotential carried in the state vector (variable index 6).

Fill mesh%nodeCoords(:,:,:,e) for one element. Place the four corner nodes from the file's node table, then perform linear (nGeo=1) placement or transfinite interpolation (nGeo>1) using the (possibly curved) edge curves. For straight sides the edge curve is filled in here by linear interpolation between corners at Chebyshev-Gauss-Lobatto parametric coordinates.

Fill mesh%nodeCoords(:,:,:,:,e) for one hexahedral element by transfinite (Coons) interpolation of its six face grids. Faces flagged in the mesh file use the file's face-point grids; the remaining faces are bilinear patches of their four corner nodes evaluated at Chebyshev-Gauss-Lobatto parametric coordinates. Edge curves are extracted from the face grids (preferring a flagged face when an edge borders one flagged and one bilinear face) and the standard Boolean-sum formula x = Pxi + Peta + Pzeta - PxiPeta - PxiPzeta - PetaPzeta + PxiPeta*Pzeta combines face, edge, and corner contributions. For an element with all-straight faces this reduces to trilinear interpolation of the eight corners.

Axis-aligned bounding box of each element. The box is taken over the element's boundary nodes (geometry%x%boundary), which are evaluated at the element edges (reference coordinate = +/-1) and therefore reach the true element boundary and corners. Using the interior nodes alone is not sufficient: for Gauss quadrature the interior nodes are strictly inside [-1,1], so a box built from them under-covers the element by a fixed fraction of its size and rejects points lying on element edges/corners. The boundary nodes give the correct extent for both Gauss and Gauss-Lobatto node sets.

Axis-aligned bounding box of each element, taken over the boundary nodes (geometry%x%boundary) which are evaluated on the element faces (reference coordinate = +/-1) and so reach the true element boundary, edges and corners. See BuildElementBBoxes_2D for why the interior nodes alone are insufficient (Gauss interior nodes under-cover the element).

Fill newGeom for the emitted mesh, reusing the previous epoch's geometry for every element that did not change and generating only the rest (AMR Stage 6c). nReused reports how many elements were copied rather than computed.

Fill newGeom for the emitted mesh, reusing the previous epoch's geometry for every element that did not change and generating only the rest. nReused reports how many elements were copied rather than computed.

Build the aggregated halo exchange tables for the mesh this domain decomposition belongs to.

Build the nodal->modal transform Pmodal for a 1-D degree-N interpolant whose nodes are controlPoints(1:N+1). The transform is the exact inverse of the L2-normalized Legendre Vandermonde V(i,p) = Ltilde_{p-1}(x_i); Pmodal(ii,p) stores V^{-1}(p,ii) so that uhat(p) = sum_ii Pmodal(ii,p) * u(ii). The Vandermonde inverse is formed in double precision for conditioning and cast back to the working precision.

Build the nodal->modal transform Pmodal for a 1-D degree-N interpolant whose nodes are controlPoints(1:N+1). The transform is the exact inverse of the L2-normalized Legendre Vandermonde V(i,p) = Ltilde_{p-1}(x_i); Pmodal(ii,p) stores V^{-1}(p,ii) so that uhat(p) = sum_ii Pmodal(ii,p) * u(ii). The Vandermonde inverse is formed in double precision for conditioning and cast back to the working precision.

Build the old->new transfer plan for the adaptation epoch that took the forest from the leaf configuration (nOld, oldLeaf) - a snapshot of (forest%nLeaves, forest%leaf) taken before mutating - to its current leaf configuration. See the module documentation for the allowed mutations within one epoch. The plan is valid until the forest is mutated again.

Build the old->new transfer plan for the adaptation epoch that took the forest from the leaf configuration (nOld, oldLeaf) - a snapshot of (forest%nLeaves, forest%leaf) taken before mutating - to its current leaf configuration. See the module documentation for the allowed mutations within one epoch. The plan is valid until the forest is mutated again.

Base method for calculating entropy of a model When this method is not overridden, the entropy is simply set to 0.0. When you develop a model built on top of this abstract class or one of its children, it is recommended that you define a convex mathematical entropy function that is used as a measure of the model stability.

GPU implementation of the EC-DG 2-D tendency.

Computes du/dt = source - EC-DG flux divergence.

GPU-resident tendency for ESAtmo2D.

ESAtmo2D tendency = EC inviscid pipeline (parent) + optional constant-coefficient Laplacian diffusion (BR1 weak-form DG).

GPU-resident tendency for ESAtmo3D. The inviscid pipeline is identical to ECDGModel3D's GPU CalculateTendency; if either nu or kappa is positive, the constant-coefficient Laplacian divergence (BR1 weak-form) is then accumulated into fluxDivergence before forming dSdt.

ESAtmo3D tendency = EC inviscid pipeline (parent) + optional constant-coefficient Laplacian diffusion (BR1 weak-form DG).

Argument validation shared by the portable and backend Estimate implementations.

Argument validation shared by the portable and backend Estimate implementations.

Convert a (signed, possibly out-of-range) real cell index to an integer cell index clamped to [0, nCells-1].

Return to the automatic energy scale (the largest gate energy over the elements).

Return to the automatic energy scale (the largest gate energy over the elements).

\addtogroup SELF_SupportRoutines @{ \fn CompareArray Compares to INTEGER arrays and determines if they are identical.

Copy whole-element geometry blocks from src into myGeom: element srcIdx(k) of src becomes element dstIdx(k) of myGeom, for k = 1..n (AMR Stage 6c).

Copy whole-element geometry blocks from src into myGeom: element srcIdx(k) of src becomes element dstIdx(k) of myGeom, for k = 1..n (AMR Stage 6c).

Returns the two corner-node IDs delimiting a given local side using SELF's 2D side convention: Side 1 South = [CN1, CN2] Side 2 East = [CN2, CN3] Side 3 North = [CN4, CN3] Side 4 West = [CN1, CN4]

Count the rank-local elements currently carrying the requested flag value (SELF_AMR_REFINE / SELF_AMR_KEEP / SELF_AMR_COARSEN). Provided as a convenience for drivers and tests; a global count across MPI ranks is the caller's responsibility.

Count the rank-local elements currently carrying the requested flag value (SELF_AMR_REFINE / SELF_AMR_KEEP / SELF_AMR_COARSEN). Provided as a convenience for drivers and tests; a global count across MPI ranks is the caller's responsibility.

Report whether a link exists at datasetName in an open file.

GPU-resident fill of diffFlux%boundaryNormal with the SIPG-stabilised BR1 diffusive flux: f = -coeff(avg_grad . n)nmag + tau(uL - uR)nmag

Fill diffFlux%boundaryNormal with the SIPG-stabilised BR1 flux:

GPU-resident fill of diffFlux%boundaryNormal with the SIPG-stabilised BR1 diffusive flux: f = -coeff(avg_grad . n)nmag + tau(uL - uR)nmag tau = eta_penaltycoeff(N+1)^2/length_scale, computed on the host.

Fill diffFlux%boundaryNormal with the SIPG-stabilised BR1 flux:

GPU-resident fill of diffFlux%interior with the constant-coefficient Laplacian flux F_d(iVar) = -coeff(iVar) * d(s_iVar)/dx_d.

Fill diffFlux%interior with the constant-coefficient Laplacian flux at every interior node:

GPU-resident fill of diffFlux%interior with the constant-coefficient Laplacian flux F_d(iVar) = -coeff(iVar) * d(s_iVar)/dx_d, where coeff is 0 for rho, nu for momentum, kappa for rho*theta.

Fill diffFlux%interior with the constant-coefficient Laplacian flux at every interior node:

GPU override of DiracDelta_2D_Points_t. When the scalar is a device-resident MappedScalar2D, the kernel writes scalar%interior_gpu directly (host array left stale; caller may UpdateHost). For non-GPU scalar subtypes, defers to the inherited host implementation.

Scatter a discrete Dirac delta of unit strength (S = 1) onto a 2D MappedScalar, one variable per stored point. Variable p of scalar receives the delta associated with point p.

GPU override of DiracDelta_3D_Points_t. See DiracDelta_2D_Points.

Scatter a discrete Dirac delta of unit strength (S = 1) onto a 3D MappedScalar, one variable per stored point. See DiracDelta_2D_Points_t for the full specification; this is the direct 3D analogue.

GPU implementation of the reference-element split-form divergence. df must be a device pointer to a scalar field of size (N+1)(N+1)nElem*nVar.

Computes the split-form (two-point) divergence of a 2-D vector field in the reference element (computational coordinates).

GPU implementation of the reference-element split-form divergence. df must be a device pointer to a scalar field of size (N+1)^3 * nElem * nVar.

Computes the split-form (two-point) divergence of a 3-D vector field in the reference element (computational coordinates).

Identifies which of the six domain boundary planes of the SimpleMortarMesh contains face s of the element, and returns the corresponding boundary condition id. The domain is [0,3dx] x [0,2dx] x [0,2dx] with the region x > 2dx, above/below the small elements, outside the mesh; the x = 3dx plane is the domain east.

Create a six-element mesh with two 2:1 mortar interfaces, used to validate the mortar machinery in configurations the SimpleMortarMesh cannot reach:

Copy the device-resident migrated window to the host, for the SELF_AMR_MIGRATE_VERIFY diagnostic only. This is exactly the host traffic the device path exists to avoid, so it must never appear on the default path.

Copy the migrated window into a host array, for the SELF_AMR_MIGRATE_VERIFY diagnostic. Diagnostic-only and off the default path: on a GPU build the override is a device-to-host transfer of the whole window, which is exactly the traffic the device path exists to avoid.

Copy the device-resident migrated window to the host, for the SELF_AMR_MIGRATE_VERIFY diagnostic only. This is exactly the host traffic the device path exists to avoid, so it must never appear on the default path.

Copy the migrated window into a host array, for the SELF_AMR_MIGRATE_VERIFY diagnostic. Diagnostic-only and off the default path: on a GPU build the override is a device-to-host transfer of the whole window, which is exactly the traffic the device path exists to avoid.

Build outMesh (a conforming-or-mortar Mesh3D_t) from a 2:1-balanced forest. baseMesh is the mesh the forest was initialised from (supplies BC metadata and the communicator; on nRanks > 1 the forest must be rank-replicated so every rank emits identical global tables). The forest must already be balanced (MaxLevelJump <= 1); EmitMesh does not mutate it.

Build outMesh (a conforming-or-mortar Mesh2D_t) from a 2:1-balanced forest. baseMesh is the mesh the forest was initialised from (supplies BC metadata and the communicator; on nRanks > 1 the forest must be rank-replicated so every rank emits identical global tables). The forest must already be balanced (MaxLevelJump <= 1); EmitMesh does not mutate it.

Grow the node arrays (amortized doubling) so at least need nodes fit.

Grow the node arrays (amortized doubling) so at least need nodes fit.

High-water-mark device allocation (AMR Stage 6b). Grows ptr to hold neededBytes, reusing the existing allocation when it already does, and records the capacity in allocBytes. Contents are not preserved across a growth - every caller rewrites the buffer before reading it.

Grow the device buffers to hold the current logical arrays, reusing existing allocations when the bytes already fit (AMR Stage 6b).

Grow the device buffers to hold the current logical arrays, reusing existing allocations when the bytes already fit (AMR Stage 6b).

Grow each device buffer to hold the current logical array, reusing the existing allocation when the bytes already fit (AMR Stage 6b). A device pointer carries no shape, so unlike the host side no remapping is needed - only a byte-capacity check.

Grow each device buffer to hold the current logical array, reusing the existing allocation when the bytes already fit (AMR Stage 6b). A device pointer carries no shape, so unlike the host side no remapping is needed - only a byte-capacity check.

Grow each device buffer to hold the current logical array, reusing the existing allocation when the bytes already fit (AMR Stage 6b/6c). A device pointer carries no shape, so only byte capacity has to be tracked.

Grow each device buffer to hold the current logical array, reusing the existing allocation when the bytes already fit (AMR Stage 6b/6c). A device pointer carries no shape, so only byte capacity has to be tracked.

Grow each device buffer to hold the current logical array, reusing the existing allocation when the bytes already fit (AMR Stage 6b). A device pointer carries no shape, so only byte capacity has to be tracked.

Grow each device buffer to hold the current logical array, reusing the existing allocation when the bytes already fit (AMR Stage 6b). A device pointer carries no shape, so only byte capacity has to be tracked.

Guarantee that pool has at least needed elements, growing it by poolGrowth when it must reallocate. Existing contents are NOT preserved: every caller re-establishes the array contents after resizing (the solution through the AMR transfer, the geometry through GenerateFromMesh, everything else by being written before it is read), so copying the old data would be wasted bandwidth.

Prepare the cached GenerateFromMesh scratch for a mesh with nGeo/quadrature and nElem elements (AMR Stage 6c). The nGeo -> N interpolant is built once and reused; the node coordinate staging is resized, so an adapting run stops rebuilding either one per epoch.

Prepare the cached GenerateFromMesh scratch for a mesh with nGeo/quadrature and nElem elements (AMR Stage 6c). The nGeo -> N interpolant is built once and reused; the node coordinate staging is resized, so an adapting run stops rebuilding either one per epoch.

Quadratic entropy: eta(u) = u^2 / 2

Quadratic entropy: eta(u) = u^2 / 2

Mathematical entropy: total energy density (kinetic + internal).

Mathematical entropy: total energy density (kinetic + internal).

Acoustic energy for the PML model. The background density is the scalar this%rho0 (the PML variant does not carry a per-node background density; its variable 5 is the auxiliary phi_u, not rho0), and the sound speed is taken from s(4). This overrides the parent entropy_func, which reads rho0 from s(5).

The entropy function is the sum of kinetic and internal energy For the linear model, this is

The entropy function is the sum of kinetic and internal energy For the linear model, this is

GPU path: one device thread per element computes the modal transform, the raw smoothness ratio and the gate energy; those are copied back to the host, where the shared second phase applies the amplitude gate, the log10 and the thresholds.

Compute the per-element modal-energy indicator sigma_e and refine/keep/coarsen flag from the nodal solution field. ivar selects the driving variable in [1,solution%nVar]; passing SELF_AMR_ALLVARS (=0) reduces the indicator over all variables by taking, per element, the largest (least smooth) smoothness ratio S_e before the log10.

GPU path: one device block per element computes the modal transform, the raw smoothness ratio and the gate energy; those are copied back to the host, where the shared second phase applies the amplitude gate, the log10 and the thresholds.

Compute the per-element modal-energy indicator sigma_e and refine/keep/coarsen flag from the nodal solution field. ivar selects the driving variable in [1,solution%nVar]; passing SELF_AMR_ALLVARS (=0) reduces the indicator over all variables by taking, per element, the largest (least smooth) smoothness ratio S_e before the log10.

Evaluate a 2D MappedScalar at the cached points on the device. The caller is responsible for hipMalloc-ing values_dev with capacity nPointsnVarprec bytes. Layout is column-major (nPoints, nVar).

Evaluate a 2D MappedScalar at all located points by tensor-product Lagrange interpolation at the stored reference coordinates. Points with elements(p) == 0 receive a value of zero. When LocatePoints has cached the per-point basis at the matching polynomial degree, this routine reuses it and skips Lagrange-polynomial evaluation altogether.

Evaluate a 3D MappedScalar at all located points by tensor-product Lagrange interpolation. Points with elements(p) == 0 receive zero. When LocatePoints has cached the basis at the matching degree, this routine reuses it.

Migrate the pre-regrid solution into this rank's old-element window with point-to-point messages: the Stage-5 v2 replacement for allgathering the whole old field onto every rank.

Migrate the pre-regrid solution into this rank's old-element window with point-to-point messages: the Stage-5 v2 replacement for allgathering the whole old field onto every rank.

The complete host-memory migration: post the messages, copy the part of the window this rank already owns while they are in flight, then wait. This is the portable backend's whole implementation; the GPU backend calls PostOldWindowExchange / FinishOldWindowExchange itself so that it can serve the local part with a device-to-device copy instead.

Find the equal-or-larger face neighbour of node across its local face s using the classic octree ascend/descend search. Returns: nbr - neighbour node id (0 = physical domain boundary). It is either a LEAF at any level <= level(node), or an INTERNAL node at exactly level(node) (meaning the shared face is subdivided on the neighbour side, i.e. finer neighbours exist). ns - the neighbour's local face that faces node. nf - the flip between the two shared faces (0..7, sideInfo(4) convention), inherited from the base-mesh face where the search crosses a root boundary. With this, a 2:1 hanging face is exactly "nbr is a leaf with level(nbr) = level(node)-1", and finer neighbours are exactly "nbr is internal".

Find the equal-or-larger face neighbour of node across its local side s using the classic quadtree ascend/descend search. Returns: nbr - neighbour node id (0 = physical domain boundary). It is either a LEAF at any level <= level(node), or an INTERNAL node at exactly level(node) (meaning the shared face is subdivided on the neighbour side, i.e. finer neighbours exist). ns - the neighbour's local side that faces node. nf - the flip between the two shared edges (0 same direction, 1 reversed), inherited from the base-mesh face where the search crosses a root boundary. With this, a 2:1 hanging face is exactly "nbr is a leaf with level(nbr) = level(node)-1", and finer neighbours are exactly "nbr is internal".

Second phase of an estimate, shared by every backend so that the flags are identical whichever computed the spectra: apply the amplitude gate, take the log10 and set the refine/keep/coarsen flags.

Second phase of an estimate, shared by every backend so that the flags are identical whichever computed the spectra: apply the amplitude gate, take the log10 and set the refine/keep/coarsen flags.

Wait on the messages PostOldWindowExchange posted. A no-op when nothing was posted, which is the common case for a symmetric refinement: under a contiguous space-filling-curve repartition each rank's new range is sourced from its own old range, so the window is entirely local and not a single byte moves.

Coordinate form of MortarFaceMap: maps receiver-face coordinates (u,v) to donor-face coordinates (u2,v2) for each of the eight SELF face flips.

Interior flux. Variables 1-4 use the parent linear Euler 2D flux; auxiliary variables 5-7 carry zero flux in both directions and are evolved purely by the PML source term.

\addtogroup SELF_SupportRoutines @{ \fn ForwardShift Shift an array integers by one index forward, moving the last index to the first.

Forward steps the model using the associated tendency procedure and time integrator

Release the forest, the indicator, and any controller-emitted mesh/geometry. The caller-owned base mesh/geometry are untouched. A model still pointing at a controller-emitted mesh must not be used after this call.

Release the forest, the indicator, and any controller-emitted mesh/geometry. The caller-owned base mesh/geometry are untouched. A model still pointing at a controller-emitted mesh must not be used after this call.

Free the 1D DG model, including GPU BC arrays.

Free the 2D DG model, including GPU BC arrays.

Free the 3D DG model, including GPU BC arrays.

Free EC Advection 2D, including GPU BC arrays.

Free EC Advection 3D, including GPU BC arrays.

Frees all memory (host and device) associated with an instance of the Lagrange class

Frees all memory (host and device) associated with an instance of the Lagrange class

Frees all memory (host and device) associated with an instance of the Lagrange_t class

Generate geometry for nElem elements directly from their mesh node coordinates (AMR Stage 6c). GenerateFromMesh is a thin wrapper over this.

Generate geometry for nElem elements directly from their mesh node coordinates (AMR Stage 6c). GenerateFromMesh is a thin wrapper over this.

Returns the node associated with the given bcid. If the bcid is not found, a null pointer is returned.

Returns the current simulation time stored in the model % t attribute

Map a global element position (gx,gy,gz), with gx in [1,nX] etc., to the tile-ordered element id used throughout the structured mesh. This is the inverse of the (i,ti) decomposition: gx = i + nxpertile*(ti-1).

Heap sort of the halo side list by ascending 64-bit key, permuting the companion arrays alongside the keys. O(n log n), in place; the relative order of equal keys is irrelevant because global side ids are unique.

GPU-accelerated mirror BC for 2D EC Advection.

Mirror boundary condition on the solution: sets the exterior state equal to the interior state. With LLF this zeroes the upwind dissipation at the wall.

Prescribed-zero boundary state for the tracer.

No-normal-flow (wall) boundary condition.

GPU-accelerated no-normal-flow BC for 2D Entropy-Stable Atmosphere. Reflects normal momentum, mirrors density, rho*theta, and Phi.

No-normal-flow boundary condition for 2D linear Euler equations. Reflects the velocity vector about the boundary normal while preserving pressure, sound speed, background density, and the relaxation rate.

GPU-accelerated no-normal-flow BC for 2D Linear Euler.

No-normal-flow BC for the PML-augmented linear Euler model. Variables 1-4 are treated identically to the parent LinearEuler2D no-normal-flow BC. Auxiliary variables 5-7 are given a zero exterior state; since they carry zero Riemann flux the exterior value is mathematically irrelevant, but zeroing it keeps the boundary state clean for diagnostics.

No-normal-flow boundary condition for 2D linear shallow water equations. Reflects the velocity vector about the boundary normal while preserving the free surface elevation.

Radiation BC: zero acoustic perturbation in the exterior state; the sound speed (variable 4), background density (variable 5) and relaxation rate (variable 6) are copied from the interior side so the Riemann solver sees a consistent c and rho0 (impedance-matched, non-reflecting outflow) and the exterior state stays a faithful mirror of the interior one for diagnostics. The relaxation rate carries no flux, so its exterior value never enters the tendency.

GPU-accelerated radiation BC for 2D Linear Euler.

Radiation BC for the PML-augmented linear Euler model: zero acoustic perturbation in the exterior state, sound speed copied from interior so the Riemann solver sees a consistent c, and auxiliary variables set to zero.

GPU-accelerated mirror BC for 3D EC Advection.

Mirror boundary condition on the solution: sets the exterior state equal to the interior state. With LLF this zeroes the upwind dissipation at the wall.

Prescribed-zero boundary state for the tracer.

No-normal-flow (wall) boundary condition.

GPU-accelerated no-normal-flow BC for 3D Entropy-Stable Atmosphere. Reflects normal momentum, mirrors density and rho*theta.

Radiation BC: zero acoustic perturbation in the exterior state; the sound speed (variable 5), background density (variable 6) and relaxation rate (variable 7) are copied from the interior side so the Riemann solver sees a consistent c and rho0 and the exterior state stays a faithful mirror of the interior one for diagnostics. The relaxation rate carries no flux, so its exterior value never enters the tendency.

Attach the controller to an initialized model. The model's current mesh becomes the forest's base mesh (level 0); its geometry interpolant drives the indicator and the solution transfer. Thresholds are the sigma = log10 modal-energy-ratio cut-offs of the refinement indicator (refineThreshold > coarsenThreshold; see SELF_RefinementIndicator_2D). ivar is the driving solution variable (or SELF_AMR_ALLVARS). maxLevel >= 0 caps the refinement depth; nHalo >= 0 sets the refine-flag halo width in elements.

Attach the controller to an initialized model. The model's current mesh becomes the forest's base mesh (level 0); its geometry interpolant drives the indicator and the solution transfer. Thresholds are the sigma = log10 modal-energy-ratio cut-offs of the refinement indicator (refineThreshold > coarsenThreshold; see SELF_RefinementIndicator_3D). ivar is the driving solution variable (or SELF_AMR_ALLVARS). maxLevel >= 0 caps the refinement depth; nHalo >= 0 sets the refine-flag halo width in elements.

Initialize the 1D DG model, then upload BC element/side arrays to GPU.

Initialize the 2D DG model, then upload BC element/side arrays to GPU.

Initialize the 3D DG model, then upload BC element/side arrays to GPU.

Initialize EC Advection 2D, then upload BC element/side arrays to GPU.

Initialize EC Advection 3D, then upload BC element/side arrays to GPU.

Initialize an instance of the Lagrange class On output, all of the attributes for the Lagrange class are allocated and values are initialized according to the number of control points, number of target points, and the types for the control and target nodes. If a GPU is available, device pointers for the Lagrange attributes are allocated and initialized.

Initialize an instance of the Lagrange class On output, all of the attributes for the Lagrange class are allocated and values are initialized according to the number of control points, number of target points, and the types for the control and target nodes. If a GPU is available, device pointers for the Lagrange attributes are allocated and initialized.

Initialize an instance of the Lagrange_t class On output, all of the attributes for the Lagrange_t class are allocated and values are initialized according to the number of control points, number of target points, and the types for the control and target nodes. If a GPU is available, device pointers for the Lagrange_t attributes are allocated and initialized.

Initialize the forest with one root per base-mesh element (all leaves at level 0). The base geometry (node coordinates) is copied so leaf geometry can be regenerated after any amount of refinement without holding a reference to the mesh. Requires a single-rank mesh (a decomposed mesh only stores its local elements); a rank-replicated forest over a decomposed base is built by gathering the global tables and calling InitGlobal.

Allocate host + device storage for nPoints in nDim dimensions. The basis-cache device buffers are allocated lazily in UpdateDevice once the polynomial degree N is known (it is set by LocatePoints).

Allocate storage for a point cloud of size nPoints in nDim dimensions. nDim must be 2 or 3.

Initialize the forest with one root per base-mesh element (all leaves at level 0). The base geometry (node coordinates) is copied so leaf geometry can be regenerated after any amount of refinement without holding a reference to the mesh. Requires a single-rank mesh (a decomposed mesh only stores its local elements); a rank-replicated forest over a decomposed base is built by gathering the global tables and calling InitGlobal.

Allocate the indicator for an interpolant of degree interp%N and nElem elements and precompute the nodal->modal transform matrix from the interpolant control points.

Allocate the indicator for an interpolant of degree interp%N and nElem elements and precompute the nodal->modal transform matrix from the interpolant control points.

Allocate the interior array for a 2-D two-point vector field. The interior array has rank 6 with layout (n,i,j,nEl,nVar,idir).

Allocate the interior array for a 3-D two-point vector field. The interior array has rank 7 with layout (n,i,j,k,nEl,nVar,idir).

Build a rank-replicated forest over a decomposed base mesh: allgather the global node coordinates, side table, and material ids (by the decomposition's contiguous element ownership) and initialize the forest from the global tables. Collective over the mesh communicator; runs once, at controller initialization.

Build a rank-replicated forest over a decomposed base mesh: allgather the global node coordinates, side table, and material ids (by the decomposition's contiguous element ownership) and initialize the forest from the global tables. Collective over the mesh communicator; runs once, at controller initialization.

Initialize the forest directly from GLOBAL base-mesh tables (one root per global base element, all leaves at level 0). This is the initialization path for a rank-replicated forest over a decomposed base mesh: every rank passes the same gathered tables and holds an identical forest. rootNbr carries global element ids (0 = physical boundary), rootNbrSide/rootFlip decode the base sideInfo(4) pairing, rootBC the base boundary-condition id per face, and rootMaterial the base material id per element.

Initialize the forest directly from GLOBAL base-mesh tables (one root per global base element, all leaves at level 0). This is the initialization path for a rank-replicated forest over a decomposed base mesh (AMR Stage 5): every rank passes the same gathered tables and holds an identical forest. rootNbr carries global element ids (0 = physical boundary), rootNbrSide/rootFlip decode the base sideInfo(4) pairing, rootBC the base boundary-condition id per side, and rootMaterial the base material id per element.

Invert the n-by-n matrix A by Gauss-Jordan elimination with partial pivoting, in double precision. Used once at initialization for the small (N+1) Legendre Vandermonde; A is a well-conditioned Vandermonde in an orthonormal basis, so a direct solve is appropriate.

Invert the n-by-n matrix A by Gauss-Jordan elimination with partial pivoting, in double precision. Used once at initialization for the small (N+1) Legendre Vandermonde; A is a well-conditioned Vandermonde in an orthonormal basis, so a direct solve is appropriate.

Inverse of a SELF face flip: the flip that maps donor-face indices back to receiver-face indices. Flips 0..4 and 6 are involutions; 5 and 7 invert each other.

Physical geometry-node coordinates of leaf leafIndex, produced by repeated exact isoparametric subdivision of its root element along the octree path. geomInterp must be a degree-nGeo Lagrange interpolant on the mesh's geometry (quadrature) nodes - the same interpolant SELF_MeshRefinement_3D builds. Level 0 leaves return the root geometry directly. Pure function of (root coords, level, octant path), so regenerated leaf geometry is bit-identical across epochs.

Physical geometry-node coordinates of leaf leafIndex, produced by repeated exact isoparametric subdivision of its root element along the quadtree path. geomInterp must be a degree-nGeo Lagrange interpolant on the mesh's geometry (quadrature) nodes - the same interpolant SELF_MeshRefinement_2D builds. Level 0 leaves return the root geometry directly.

Run the host spatial-hash + Newton search (inherited from Points_t), then mirror per-point state to the device.

Locate each stored physical point inside the 2D SEMQuad geometry. On exit, elements(p) holds the (rank-local) element id and coordinates(p,:) holds the (s,t) reference coordinates, both for points that resolve. Points that resolve to no element are marked with elements(p) = 0.

Locate each stored physical point inside the 3D SEMHex geometry. On exit, elements(p) and coordinates(p,1:3) hold the located element id and reference (s,t,u) for points that resolve; elements(p) = 0 otherwise.

Numerically stable logarithmic mean (Ismail-Roe 2009 / Ranocha 2018).

Numerically stable logarithmic mean (Ismail-Roe 2009 / Ranocha 2018).

Size the backing pools for (Np,nVar,nElem) and remap the public arrays onto them at the exact logical shape. Bounds-remapping a rank-1 contiguous target to a rank-4 pointer is what keeps every extent and stride identical to a plain allocate while allowing the pool underneath to be larger and reused. See SELF_DataPool.

Size the backing pools for (Np,nVar,nElem) and remap the public arrays onto them at the exact logical shape. Bounds-remapping a rank-1 contiguous target to a rank-5 pointer is what keeps every extent and stride identical to a plain allocate while allowing the pool underneath to be larger and reused. See SELF_DataPool.

Size the backing pools for (Np,nVar,nElem) and remap the public arrays onto them at the exact logical shape. See SELF_DataPool for why this is done with pointer remapping rather than by over-allocating the element dimension.

Size the backing pools for (Np,nVar,nElem) and remap the public arrays onto them at the exact logical shape. See SELF_DataPool for why this is done with pointer remapping rather than by over-allocating the element dimension.

Size the backing pools for (Np,nVar,nElem) and remap the public arrays onto them at the exact logical shape. See SELF_DataPool for why this is done with pointer remapping rather than by over-allocating the element dimension.

Size the backing pools for (Np,nVar,nElem) and remap the public arrays onto them at the exact logical shape. See SELF_DataPool for why this is done with pointer remapping rather than by over-allocating the element dimension.

Scan the mesh boundary condition IDs and populate the elements/sides arrays for each registered boundary condition.

Scan the mesh sideInfo and populate the elements/sides arrays for each registered boundary condition.

Scan the mesh sideInfo and populate the elements/sides arrays for each registered boundary condition.

Computes the divergence of a 2-D vector using the weak form On input, the attribute of the vector is assigned and the attribute is set to the physical directions of the vector. This method will project the vector onto the contravariant basis vectors.

Computes the divergence of a 3-D vector using the weak form On input, the attribute of the vector is assigned and the attribute is set to the physical directions of the vector. This method will project the vector onto the contravariant basis vectors.

Computes the divergence of a 3-D vector using the weak form On input, the attribute of the vector is assigned and the attribute is set to the physical directions of the vector. This method will project the vector onto the contravariant basis vectors.

Calculates the gradient of a function using the weak form of the gradient and the average boundary state. This method will compute the average boundary state from the and attributes of

Calculates the gradient of a function using the weak form of the gradient and the average boundary state. This method will compute the average boundary state from the and attributes of

GPU implementation of the physical-space split-form divergence for a 2-D two-point vector on a curvilinear mesh.

Computes the physical-space divergence of a 2-D split-form vector field on a curvilinear mesh.

GPU implementation of the physical-space split-form divergence for a 3-D two-point vector on a curvilinear mesh.

Computes the physical-space divergence of a 3-D split-form vector field on a curvilinear mesh.

Calculates the gradient of a function using the strong form of the gradient in mapped coordinates.

Calculates the gradient of a function using the strong form of the gradient in mapped coordinates.

Calculates the gradient of a function using the strong form of the gradient in mapped coordinates.

Calculates the gradient of a function using the strong form of the gradient in mapped coordinates.

Determines whether face sA of element A coincides with face sB of element B, and if so with which flip: face_B(F_flip(u,v)) == face_A(u,v) at the sampled points. Bilinear faces are pinned by their four corners plus one interior point.

Highest refinement level among the active leaves.

Highest refinement level among the active leaves.

Largest refinement-level difference across any leaf face (0 on a conforming or uniformly refined forest, 1 on a 2:1-balanced adaptive forest). Because FaceNeighbor returns the equal-or-larger neighbour, every level difference is observed from the finer leaf as a coarser leaf neighbour; internal (finer) neighbours contribute nothing from this side.

Largest refinement-level difference across any leaf face (0 on a conforming or uniformly refined forest, 1 on a 2:1-balanced adaptive forest). Because FaceNeighbor returns the equal-or-larger neighbour, every level difference is observed from the finer leaf as a coarser leaf neighbour; internal (finer) neighbours contribute nothing from this side.

Assemble this rank's window of the pre-regrid old field IN DEVICE MEMORY, replacing the base implementation's host-memory migration. The rank's own run is copied device-to-device and the peers' runs are received directly into device memory, so no solution data crosses the host link and the windowed ApplyTransferPlan that follows runs as a device kernel.

Migrate the pre-regrid solution into this rank's old-element window, ready for a windowed ApplyTransferPlan. Call it BEFORE Regrid, which releases the storage the sends read:

Assemble this rank's window of the pre-regrid old field IN DEVICE MEMORY, replacing the base implementation's host-memory migration. The rank's own run is copied device-to-device and the peers' runs are received directly into device memory, so no solution data crosses the host link and the windowed ApplyTransferPlan that follows runs as a device kernel.

Migrate the pre-regrid solution into this rank's old-element window, ready for a windowed ApplyTransferPlan. Call it BEFORE Regrid, which releases the storage the sends read:

GPU implementation of the mortar exchange (see the base class for the algorithm) : traces are staged, reoriented, restricted, and projected in device memory with the SELF_Mortar kernels.

Fills the extBoundary attribute on all sides that participate in a 2:1 nonconforming (mortar) interface.

GPU implementation of the mortar exchange (see the base class for the algorithm) : traces are staged, reoriented, restricted, and projected in device memory with the SELF_Mortar kernels.

Fills the extBoundary attribute on all faces that participate in a 2:1 nonconforming (mortar) interface.

GPU implementation of the vector mortar exchange; the kernels treat the (variable, direction) pairs as 2*nvar independent trace lines.

Fills the extBoundary attribute on all sides participating in a 2:1 nonconforming (mortar) interface; vector analogue of the scalar MortarExchange (see MappedScalar2D_t for the algorithm description).

GPU implementation of the vector mortar exchange; the kernels treat the (variable, direction) pairs as 3*nvar independent trace lines.

Fills the extBoundary attribute on all faces participating in a 2:1 nonconforming (mortar) interface; vector analogue of the scalar MortarExchange (see MappedScalar3D_t for the algorithm description).

Maps receiver-face node indices (i,j) to donor-face node indices (i2,j2) for each of the eight SELF face flips (see the CGNStoSELFflip table above; indices here are 1-based, so the 0-based rule i2 = N-i1 reads i2 = N+2-i). Mortar staging uses this to reorient a small face's trace into the big face's coordinates and to scatter restricted traces back.

GPU implementation of MortarFluxCollect (see the base class for the algorithm and conservation statement). Stages the small sides' boundaryNormal traces in the mortar buffer, then overwrites the big side's integrand on device.

Replaces the big-side boundaryNormal trace on each mortar interface with the L2 projection of the two small sides' boundaryNormal traces.

GPU implementation of MortarFluxCollect (see the base class for the algorithm and conservation statement). Stages the small faces' boundaryNormal traces in the mortar buffer, then overwrites the big face's integrand on device.

Replaces the big-face boundaryNormal trace on each mortar interface with the L2 projection of the four small faces' boundaryNormal traces.

Post the aggregated halo exchange: one MPI_Irecv/MPI_Isend pair per neighboring rank, carrying every (side,variable,component) boundary trace shared with that rank in a single packed device buffer. The boundary array is laid out with the component index outermost, so the pack/unpack kernels treat the vector as 2*nvar scalar variables. Packed buffers are allocated on first use; the shared side tables are built by SideExchange before this is called.

Post the aggregated halo exchange: one MPI_Irecv/MPI_Isend pair per neighboring rank, carrying every (side,variable,component) boundary trace shared with that rank in a single packed device buffer. The boundary array is laid out with the component index outermost, so the pack/unpack kernels treat the vector as 3*nvar scalar variables. Packed buffers are allocated on first use; the shared side tables are built by SideExchange before this is called.

GPU-resident analogue of the base-class mortar message posting; messages are posted on device memory (GPU-aware MPI), following MPIExchangeAsync.

Posts the point-to-point messages required for mortar interfaces whose big and small elements reside on different ranks. The big-side rank sends its edge trace to each remote small-side rank and receives each remote small-side trace; tags are built from the sub-edge global side ids, following the conforming SideExchange tag convention.

GPU-resident analogue of the base-class mortar message posting; messages are posted on device memory (GPU-aware MPI), following MPIExchangeAsync.

Posts the point-to-point messages required for mortar interfaces whose big and small elements reside on different ranks. The big-face rank sends its face trace to each remote small-face rank and receives each remote small-face trace; tags are built from the sub-face global side ids, following the conforming SideExchange tag convention.

GPU-resident mortar message posting for vector data; messages are posted on device memory (GPU-aware MPI), one per variable and physical direction.

Vector analogue of the scalar mortar exchange message posting; each physical direction of each variable is exchanged as its own message.

GPU-resident analogue of the base-class vector mortar message posting; messages are posted on device memory (GPU-aware MPI).

Posts the point-to-point messages required for mortar interfaces whose big and small elements reside on different ranks; vector analogue of the scalar MPIMortarExchangeAsync with one message per direction component.

One-directional messages for MortarFluxCollect on device memory : each remote small side sends its boundaryNormal trace to the big side's rank.

Posts the one-directional messages for MortarFluxCollect : each remote small side sends its boundaryNormal trace to the big side's rank.

Posts the one-directional messages for MortarFluxCollect on device memory : each remote small face sends its boundaryNormal trace to the big face's rank.

Posts the one-directional messages for MortarFluxCollect : each remote small face sends its boundaryNormal trace to the big face's rank.

Newton inverse-map for the 2D SEM element iEl: solve X(xi) = xTarget for the reference coordinate xi in R^2, where X is the high-order interpolant. The Jacobian dX/dxi is taken from geometry%dxds (the covariant basis tensor), interpolated to xi via Lagrange basis.

Newton inverse-map for the 3D SEM element iEl. See NewtonInverse_2D.

Select the geometry buffer to fill this epoch: whichever of the two is not currently active, so the previous epoch's geometry stays intact and readable (AMR Stage 6c).

Select the geometry buffer to fill this epoch: whichever of the two is not currently active, so the previous epoch's geometry stays intact and readable.

L2-normalized Legendre polynomial Ltilde_p(x) = L_p(x)*sqrt((2p+1)/2) on [-1,1], evaluated with the standard three-term recurrence. The normalization gives int_{-1}^{1} Ltilde_p Ltilde_q dx = delta_pq.

L2-normalized Legendre polynomial Ltilde_p(x) = L_p(x)*sqrt((2p+1)/2) on [-1,1], evaluated with the standard three-term recurrence. The normalization gives int_{-1}^{1} Ltilde_p Ltilde_q dx = delta_pq.

.true. if child c's face s is interior to its parent (its neighbour across s is a sibling): the child sits in the half-cube away from parent face s.

The local face directly across an element from face s.

Sibling child index obtained by reflecting c across face s (swap the half-index of the direction normal to s).

Face quadrant (1..4, in face s's trace coordinates) of child c on face s; 0 if c does not touch face s.

The run of elements that rank r-1 owns and that also lies in [first,last]: a..b, empty when b < a. offsetElem is a decomposition's contiguous ownership table, so rank r-1 owns offsetElem(r)+1 .. offsetElem(r+1).

No-stress boundary condition for the BR1 diffusive flux. Reflects the interior gradient so that the normal component of the averaged gradient vanishes at the boundary: sigma_ext = sigma_int - 2 (sigma_int . nhat) nhat This gives zero diffusive flux through the wall and is unconditionally stable for Bassi-Rebay.

Parabolic boundary condition: zero diffusive flux normal to the wall (no-stress for momentum, no-heat-flux for rho*theta).

GPU-accelerated parabolic no-stress / no-heat-flux BC for 2D Entropy-Stable Atmosphere.

No-stress boundary condition for the BR1 diffusive flux. Reflects the interior gradient so that the normal component of the averaged gradient vanishes at the boundary: sigma_ext = sigma_int - 2 (sigma_int . nhat) nhat This gives zero diffusive flux through the wall and is unconditionally stable for Bassi-Rebay.

Parabolic boundary condition: zero diffusive flux normal to the wall (no-stress for momentum, no-heat-flux for rho*theta).

GPU-accelerated parabolic no-stress / no-heat-flux BC for 3D Entropy-Stable Atmosphere. Reflects the normal component of the solution gradient at every wall node so that BR1 averaging gives avgGrad . n = 0 (zero diffusive flux through the wall) for every variable.

For every rank r, the contiguous window [winFirst(r),winLast(r)] of GLOBAL old element indices that rank r's new element range references through the transfer plan. The plan is rank-replicated, so this is a purely local computation: no communication is needed to learn what a peer wants, which is what makes the point-to-point migration cheap here.

For every rank r, the contiguous window [winFirst(r),winLast(r)] of GLOBAL old element indices that rank r's new element range references through the transfer plan. The plan is rank-replicated, so this is a purely local computation: no communication is needed to learn what a peer wants, which is what makes the point-to-point migration cheap here.

Populate the elements and sides arrays for a registered boundary condition. Called after scanning the mesh to determine which faces belong to each bcid.

Post every message of this rank's migration schedule and return without waiting: the receives that fill the remote runs of its window, then the sends that serve the runs of its old range that other ranks' windows need. Receives are posted before sends, as everywhere else in SELF.

PostStepHook is a template routine invoked by the time integrators once immediately after each completed time step (after this%t has advanced by this%dt). The default implementation is a no-op.

PreStepHook is a template routine invoked by the time integrators once immediately before each time step is taken (before any Runge-Kutta stages of that step). The default implementation is a no-op.

PreTendencyHook is a template routine that is used to house any additional calculations that you want to execute at the beginning of the tendency calculation routine. This default PreTendencyHook simply returns back to the caller without executing any instructions

Prolong (interpolate) a parent element's nodal solution onto its four children. uChildren(:,:,:,c) is the field on child c (ordering 1=SW,2=SE,3=NE,4=NW). Exact for the parent's degree-N polynomial representation.

Prolong (interpolate) a parent element's nodal solution onto its eight children. uChildren(:,:,:,:,c) is the field on child c (SELF child ordering = CGNS corner order). Exact for the parent's degree-N polynomial representation.

.true. if child c's side s is interior to its parent (its neighbour across s is a sibling).

The local side directly across an element from side s.

Sibling child index obtained by reflecting c across side s (swap the y-half for the horizontal faces S/N, the x-half for the vertical faces E/W).

Sub-position (1 or 2, in the positive direction of side s) of child c along side s; 0 if c does not touch side s.

Determines whether face sSmall of the small element coincides with quadrant q of face sBig of the big element under the given flip: face_small(F_flip(u,v)) == face_big(quadrant_q(u,v)) at the sampled points.

Reader for HOHQMesh text mesh files in the ISM and ISM-MM formats. The format is auto-detected from the first line: * Line equal to "ISM-MM" (trimmed) => ISM-MM with per-element material name strings and a 4-int count line that includes an unused nEdges field (the ISM-MM writer in HOHQMesh does NOT emit an edge block). * Anything else is treated as plain ISM: the first line is itself the count line "nNodes nElems polyOrder" and there are no per-element material names.

Reader for HOHQMesh 3-D (hexahedral) text mesh files in the ISM and ISM-MM formats, as written by HOHQMesh's WriteISMHexMeshFile (Source/3DSource/Mesh3DOutputMethods.f90). Unlike the 2-D writer, the 3-D writer emits NO format header line: the first line is always the count line "nNodes nElems polyOrder". The two variants differ only in the per-element corner-node line: * ISM : 8 corner-node ids * ISM-MM : 8 corner-node ids followed by a material-name string The variant is auto-detected from the presence of the 9th token.

Recompute the active leaf set by depth-first traversal from the roots. Traversal (rather than a scan of all nodes) is what makes orphaned nodes left behind by coarsening invisible. The root-major DFS order with children visited 1..8 is the Morton (Z-order) curve within each root, which is what makes the contiguous block decomposition of the emitted element list a locality-preserving SFC partition.

Recompute the active leaf set by depth-first traversal from the roots. Traversal (rather than a scan of all nodes) is what makes orphaned nodes left behind by coarsening invisible.

Level-based explicit-stability time step: refinement halves the element scale per level, so a time step dtBase that is stable on the base (level-0) mesh scales to dtBase / 2**MaxLevel on the current forest. Deterministic and exact for the quadtree (child elements are exact half-scale subdivisions); no geometry reduction is needed.

Level-based explicit-stability time step: refinement halves the element scale per level, so a time step dtBase that is stable on the base (level-0) mesh scales to dtBase / 2**MaxLevel on the current forest. Deterministic and exact for the octree (child elements are exact half-scale subdivisions); no geometry reduction is needed.

Build the connectivity of a uniformly (2:1 in each direction) refined 2-D mesh from its base connectivity. Every base element p is split into four children (global id 4*(p-1)+c). The scheme is fully deterministic integer bookkeeping - no coordinate hashing, no flip recomputation, no MPI:

Build the connectivity of a uniformly (2:1 in each direction) refined 3-D mesh from its base connectivity. Every base element p is split into eight children (global id 8*(p-1)+c). The scheme is fully deterministic integer bookkeeping - no coordinate hashing, no flip recomputation, no MPI:

Subdivide one leaf node into eight children. Does nothing (with a warning) if the node is already refined. The caller is responsible for rebuilding the leaf list afterwards (or calling AdaptFromFlags, which does so).

Subdivide one leaf node into four children. Does nothing (with a warning) if the node is already refined. The caller is responsible for rebuilding the leaf list afterwards (or calling AdaptFromFlags, which does so).

Register a boundary condition function with the given bcid and bcname. If the bcid is already registered, the function pointer is updated. The elements and sides arrays are not allocated here; call PopulateBoundaries after scanning the mesh.

GPU regrid: release the device copies of the old boundary-condition element/side arrays, rebuild the model storage and BC maps on the new mesh (base Regrid), then upload the new BC arrays to the device - the same device bookkeeping Init/Free do around the base Init/Free.

Rebind a live model to a new mesh/geometry pair (AMR regrid). The mesh-sized solution storage is reallocated and the boundary-condition registrations and maps are rebuilt for the new mesh, while everything that is not mesh-sized is preserved: the time state (t, dt, entropy, IO counter), the time-integrator selection, configuration flags, and any model-specific parameters (Init is intent(out) and would reset all of these). nvar/nstepped are unchanged - the model solves the same equations on a new mesh.

GPU regrid: release the device copies of the old boundary-condition element/side arrays, rebuild the model storage and BC maps on the new mesh (base Regrid), then upload the new BC arrays to the device - the same device bookkeeping Init/Free do around the base Init/Free.

Rebind a live model to a new mesh/geometry pair (AMR regrid). The mesh-sized solution storage is reallocated and the boundary-condition registrations and maps are rebuilt for the new mesh, while everything that is not mesh-sized is preserved: the time state (t, dt, entropy, IO counter), the time-integrator selection, configuration flags, and any model-specific parameters (Init is intent(out) and would reset all of these). nvar/nstepped are unchanged - the model solves the same equations on a new mesh.

Retire one decomposition from the process-wide lifecycle bookkeeping and, if it was the last one and SELF initialized MPI, finalize.

Base method for reporting the entropy of a model to stdout. Only override this procedure if additional reporting is needed. Alternatively, if you think additional reporting would be valuable for all models, open a pull request with modifications to this base method.

Base method for reporting the entropy of a model to stdout. Only override this procedure if additional reporting is needed. Alternatively, if you think additional reporting would be valuable for all models, open a pull request with modifications to this base method.

Base method for reporting the entropy of a model to stdout. Only override this procedure if additional reporting is needed. Alternatively, if you think additional reporting would be valuable for all models, open a pull request with modifications to this base method.

Method that can be overridden by users to report their own custom metrics after file io

Warn, once, about domain endpoints whose bcid has no registered boundary condition. MapBoundaryConditions establishes the count; ForwardStep calls this.

Warn, once, about mesh boundary edges whose bcid has no registered boundary condition. MapBoundaryConditions establishes the count; ForwardStep calls this.

Warn, once, about mesh boundary faces whose bcid has no registered boundary condition. MapBoundaryConditions establishes the count; ForwardStep calls this.

Reports mesh boundary faces that carry a boundary condition id for which no boundary condition has been registered. The exterior state on such a face is never written by any boundary condition method, so the Riemann solver consumes whatever extBoundary last held.

Method that can be overridden by users to report their own custom metrics after file io

This method can be used to reset all of the boundary elements boundary condition type to the desired value.

This method can be used to reset all of the boundary elements boundary condition type to the desired value.

This method can be used to reset all of the boundary elements boundary condition type to the desired value.

Rebind to a new element count, reusing host pools and device buffers where they fit (AMR Stage 6b), and without uploading zeros the way Init does.

Rebind to a new element count, reusing storage where it fits (AMR Stage 6b). The mortar staging buffer is sized by mesh%nMortars rather than nElem, and nMortars changes with the mesh independently of the element count, so it is invalidated here and lazily re-created at the correct size by the next mortar exchange.

Rebind to a new element count, reusing host pools and device buffers where they fit (AMR Stage 6b), and without uploading zeros the way Init does.

Rebind to a new element count, reusing storage where it fits (AMR Stage 6b). The mortar staging buffer is sized by mesh%nMortars rather than nElem, and nMortars changes with the mesh independently of the element count, so it is invalidated here and lazily re-created at the correct size by the next mortar exchange.

Rebind to a new element count, reusing host pools and device buffers where they fit (AMR Stage 6b). Inherits the Vector2D resize (host pools plus the five device buffers) and adds this class's mortar buffers, which are sized by mesh%nMortars rather than nElem: nMortars changes with the mesh independently of the element count, so a stale buffer would silently under-size the next mortar exchange. They are invalidated here and lazily re-created at the correct size on next use.

Rebind to a new element count, reusing storage where it fits (AMR Stage 6b). See Resize_MappedScalar2D_t for why the mortar staging buffer is invalidated rather than resized.

Rebind to a new element count, reusing host pools and device buffers where they fit (AMR Stage 6b). Inherits the Vector3D resize (host pools plus the five device buffers) and adds this class's mortar buffers, which are sized by mesh%nMortars rather than nElem: nMortars changes with the mesh independently of the element count, so a stale buffer would silently under-size the next mortar exchange. They are invalidated here and lazily re-created at the correct size on next use.

Rebind to a new element count, reusing storage where it fits (AMR Stage 6b). See Resize_MappedScalar3D_t for why the mortar staging buffer is invalidated rather than resized.

Rebind to a new element count, reusing host pools and device buffers where they fit. Deliberately does NOT call UpdateDevice: Init uploads freshly zeroed arrays, which in the adaptive loop is pure waste because the solution transfer overwrites them immediately.

Rebind a live object to a new element count, reusing the existing storage when it fits (AMR Stage 6b). Unlike Init - which is intent(out) and therefore resets the whole object, reallocates, zeroes and, on GPU builds, uploads those zeros - this preserves metadata and equation parsers (neither depends on nElem) and touches storage only when it must grow.

Rebind to a new element count, reusing host pools and device buffers where they fit. Deliberately does NOT call UpdateDevice: Init uploads freshly zeroed arrays, which in the adaptive loop is pure waste because the solution transfer overwrites them immediately.

Rebind a live object to a new element count, reusing the existing storage when it fits (AMR Stage 6b). Unlike Init - which is intent(out) and therefore resets the whole object, reallocates, zeroes and, on GPU builds, uploads those zeros - this preserves metadata and equation parsers (neither depends on nElem) and touches storage only when it must grow.

Rebind a live geometry to a new element count, reusing storage where it fits (AMR Stage 6c). This replaces the Free + Init cycle the adaptive loop performed on a freshly allocated SEMHex every epoch, which threw away exactly the amortization Stage 6b introduced: each member Free released its pools and device buffers, and each Init reallocated, zeroed, rebuilt metadata and equation parsers, and uploaded the zeros.

Rebind a live geometry to a new element count, reusing storage where it fits (AMR Stage 6c). This replaces the Free + Init cycle the adaptive loop performed on a freshly allocated SEMQuad every epoch, which threw away exactly the amortization Stage 6b introduced: each member Free released its pools and device buffers, and each Init reallocated, zeroed, rebuilt metadata and equation parsers, and uploaded the zeros.

Rebind to a new element count, reusing host pools and device buffers where they fit. Deliberately does NOT call UpdateDevice: Init uploads freshly zeroed arrays, which is pure waste in the adaptive loop because the arrays are rewritten immediately after.

Rebind a live object to a new element count, reusing existing storage when it fits (AMR Stage 6b/6c). Unlike Init - which is intent(out), so it resets the object, reallocates, zeroes and reconstructs the 4*nVar equation parsers - this preserves metadata and parsers (neither depends on nElem) and touches storage only when it must grow.

Rebind to a new element count, reusing host pools and device buffers where they fit. Deliberately does NOT call UpdateDevice: Init uploads freshly zeroed arrays, which is pure waste in the adaptive loop because the arrays are rewritten immediately after.

Rebind a live object to a new element count, reusing existing storage when it fits (AMR Stage 6b/6c). Unlike Init - which is intent(out), so it resets the object, reallocates, zeroes and reconstructs the 9*nVar equation parsers - this preserves metadata and parsers (neither depends on nElem) and touches storage only when it must grow.

Rebind to a new element count, reusing host pools and device buffers where they fit. Deliberately does NOT call UpdateDevice: Init uploads freshly zeroed arrays, which is pure waste in the adaptive loop because the fields are rewritten immediately after.

Rebind a live object to a new element count, reusing existing storage when it fits (AMR Stage 6b). Unlike Init - which is intent(out), so it resets the object, reallocates, zeroes and reconstructs the equation parsers - this preserves metadata and parsers (neither depends on nElem) and touches storage only when it must grow. Preserving the parsers also avoids repeating the EquationParser construction that Init performs 2*nVar times as an amdflang workaround.

Rebind to a new element count, reusing host pools and device buffers where they fit. Deliberately does NOT call UpdateDevice: Init uploads freshly zeroed arrays, which is pure waste in the adaptive loop because the fields are rewritten immediately after.

Rebind a live object to a new element count, reusing existing storage when it fits (AMR Stage 6b). Unlike Init - which is intent(out), so it resets the object, reallocates, zeroes and reconstructs the equation parsers - this preserves metadata and parsers (neither depends on nElem) and touches storage only when it must grow. Preserving the parsers also avoids repeating the EquationParser construction that Init performs 3*nVar times as an amdflang workaround.

Determine the energy scale the relative floor is measured against: the pinned value when SetEnergyScale was used, otherwise the largest gate energy over the (rank-local) elements. When comm is present the maximum is taken over the communicator, so every rank applies the same floor and the flags - hence the adapted mesh - do not depend on the decomposition. One small collective per indicator evaluation, i.e. per adaptation epoch, never inside the time-stepping loop.

Determine the energy scale the relative floor is measured against: the pinned value when SetEnergyScale was used, otherwise the largest gate energy over the (rank-local) elements. When comm is present the maximum is taken over the communicator, so every rank applies the same floor and the flags - hence the adapted mesh - do not depend on the decomposition. One small collective per indicator evaluation, i.e. per adaptation epoch, never inside the time-stepping loop.

Ensure energyWeight is associated and sized to nVar. Weights supplied by SetEnergyWeights must match nVar; otherwise the default weights are regenerated: unit weight on the driving variable, or on every variable when ivar is SELF_AMR_ALLVARS. Resolved here rather than at Init because the variable count belongs to the solution field, not to the indicator.

Ensure energyWeight is associated and sized to nVar. Weights supplied by SetEnergyWeights must match nVar; otherwise the default weights are regenerated: unit weight on the driving variable, or on every variable when ivar is SELF_AMR_ALLVARS. Resolved here rather than at Init because the variable count belongs to the solution field, not to the indicator.

Read the incremental-geometry and solution-migration debug switches once. Idempotent, so it is safe (and intended) to call from every path that consults a switch rather than relying on one caller having run first.

Read the Stage 6c geometry and the solution-migration debug switches once. Idempotent, so it is safe (and intended) to call from every path that consults a switch rather than relying on one caller having run first.

Restrict (L2-project) the solution on four children back onto their parent element. Conservative, and the exact left inverse of ProlongToChildren.

Restrict (L2-project) the solution on eight children back onto their parent element. Conservative, and the exact left inverse of ProlongToChildren.

Local Lax-Friedrichs (Rusanov) Riemann flux for linear advection. Entropy-stable: provides symmetric dissipation at element interfaces. Uses the per-face spectral radius lambda = |u.n|, matching upwind (Godunov) for linear advection: tangential faces (u.n=0) contribute zero flux, avoiding spurious dissipation across element interfaces whose face normal is perpendicular to the velocity.

Local Lax-Friedrichs (Rusanov) Riemann flux. Provided as a fallback; the model overrides BoundaryFlux directly with the LMARS solver.

Impedance-matched Riemann flux for acoustic variables 1-4, with zero flux returned for the auxiliary variables 5-7. The acoustic formula is identical to the parent LinearEuler2D model; see riemannflux2d_LinearEuler2D_t for the derivation.

Characteristic-decomposition (impedance-matched) interface flux for linear acoustics with possibly discontinuous sound speed.

Local Lax-Friedrichs (Rusanov) Riemann flux for linear advection. Uses the per-face spectral radius lambda = |u.n|, matching upwind (Godunov) for linear advection: tangential faces (u.n=0) contribute zero flux, avoiding spurious dissipation across element interfaces whose face normal is perpendicular to the velocity.

Local Lax-Friedrichs (Rusanov) Riemann flux.

Impedance-matched (characteristic/Godunov) Riemann flux, identical in form to the 2-D model's. The interface states are resolved with the per-side acoustic impedances Z = rho0*c (each side using its own background density rho0 and sound speed c), so material interfaces in a heterogeneous field are handled with the physically correct transmission/reflection:

Apply boundary conditions for the solution on GPU. Syncs boundary data from device, applies host-side BC dispatch (periodic defaults + registered BCs), then syncs back to device.

Apply boundary conditions for the solution. Periodic boundaries are set as the default; registered boundary conditions overwrite specific endpoints.

Apply registered boundary conditions for the solution. Each boundary condition method loops over its own boundary faces.

Apply registered boundary conditions for the solution.

Set the description of the ivar-th variable

Set the constant-coefficient Laplacian diffusion coefficients (kinematic momentum diffusivity and thermal diffusivity, both in m^2/s) and the dimensionless SIPG jump penalty. Setting nu or kappa > 0 enables the gradient pipeline so that the diffusive flux methods receive solutionGradient.

Set the constant-coefficient Laplacian diffusion coefficients (kinematic momentum diffusivity and thermal diffusivity, both in m^2/s) and the dimensionless SIPG jump penalty. Setting nu or kappa > 0 enables the gradient pipeline so that the diffusive flux methods receive solutionGradient.

Pin the energy scale that normalizes the relative floor, in the units of the gate energy (squared field amplitude times the reference-element area). Two reasons to use this rather than the automatic maximum over elements:

Pin the energy scale that normalizes the relative floor, in the units of the gate energy (squared field amplitude times the reference-element volume). Two reasons to use this rather than the automatic maximum over elements:

Set the per-variable weights w_v >= 0 of the gate energy g_e = sum_v w_v E_tot,e,v, where E_tot,e,v is the exact L2 energy of variable v on the reference element.

Set the per-variable weights w_v >= 0 of the gate energy g_e = sum_v w_v E_tot,e,v, where E_tot,e,v is the exact L2 energy of variable v on the reference element.

Sets the equation parser for the ivar-th variable

Sets the equation parser for the idir direction and ivar-th variable

Sets the equation parser for the idir direction and ivar-th variable

Apply gradient boundary conditions on GPU. Syncs gradient boundary data from device, applies host-side BC dispatch (periodic defaults + registered BCs), then syncs back to device.

Apply boundary conditions for the solution gradient. Periodic boundaries are set as the default; registered boundary conditions overwrite specific endpoints.

Apply registered boundary conditions for the solution gradient. Each boundary condition method loops over its own boundary faces.

Apply registered boundary conditions for the solution gradient.

Initialise a hydrostatically balanced atmosphere with uniform potential temperature theta0, zero velocity, and the geopotential Phi = g*y carried as state variable index 5.

Initialise a hydrostatically balanced atmosphere with uniform potential temperature theta0, zero velocity, and the geopotential Phi = g*z carried as state variable index 6.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Sets the this % interior attribute using the eqn attribute, geometry (for physical positions), and provided simulation time.

Set the name of the ivar-th variable

Five conserved variables: (rho, rhou, rhov, rhotheta, Phi), where Phi = gy is the geopotential. Phi has zero flux (volume and surface) so its tendency is identically zero; it is carried in the state vector solely so that the Souza et al. (2023) non-conservative gravity flux differencing in SourceMethod can read it node-locally.

Six conserved variables: (rho, rhou, rhov, rhow, rhotheta, Phi), where Phi = g*z is the geopotential. Phi has zero flux (volume and surface) so its tendency is identically zero; it is carried in the state vector solely so that the Souza et al. (2023) non-conservative gravity flux differencing in SourceMethod can read it node-locally.

Populate the per-node sigma_x, sigma_y fields. Only nodes inside elements whose material name starts with this%pml_material_prefix receive non-zero damping; all other nodes are forced to zero so the interior solution remains unmodified.

Copy user-supplied physical coordinates into the cloud. xIn must be shape (nPoints, nDim) matching the init dimensions.

Set the relative energy floor of the amplitude gate. An element whose gate energy g satisfies

Set the relative energy floor of the amplitude gate. An element whose gate energy g satisfies

Sets the model % t attribute with the provided simulation time

Update the refine/coarsen thresholds without rebuilding the transform.

Update the refine/coarsen thresholds without rebuilding the transform.

Sets the time integrator method, using a character input

Set the units of the ivar-th variable

Complete the aggregated halo exchange started by SideExchangeStart: wait on the persistent requests and unpack the received traces into extBoundary, applying side flips. The host blocks in MPI_Waitall (which also drives MPI progress) while previously launched device kernels continue to execute. No-op if no exchange is in flight.

Complete the aggregated halo exchange started by SideExchangeStart: wait on the persistent requests and unpack the received traces into extBoundary, applying side flips. The host blocks in MPI_Waitall (which also drives MPI progress) while previously launched device kernels continue to execute. No-op if no exchange is in flight.

Begin the aggregated halo exchange and launch the local (same-rank) side exchange kernel.

Begin the aggregated halo exchange and launch the local (same-rank) side exchange kernel.

Create the smallest 2:1 nonconforming (mortar) mesh: one "big" element of size 2dx x 2dx whose east edge is shared with the west edges of two "small" dx x dx elements. The mesh is conforming everywhere except at the single mortar interface.

Create the smallest 3D 2:1 nonconforming (mortar) mesh: one "big" element of size 2dx x 2dx x 2*dx whose east face is shared with the faces of four "small" dx x dx x dx elements stacked 2x2 in (y,z). The mesh is conforming everywhere except at the single mortar interface.

No source term — upload the zero-initialised host array to device.

No source term — upload the zero-initialised host array to device.

Souza et al. (2023) non-conservative gravity flux differencing on GPU. The geopotential lives at solution(:,:,:,5); the source for rho*v is computed via the SBP-EC two-point form using log-mean density and the contravariant metric.

Souza et al. (2023) non-conservative gravity flux differencing.

Souza et al. (2023) non-conservative gravity flux differencing on GPU. The geopotential lives at solution(:,:,:,:,6); the source for rho*w is computed via the SBP-EC two-point form using log-mean density and the contravariant metric. Fully device-resident.

Souza et al. (2023) non-conservative gravity flux differencing.

Device-resident spatially-dependent linear relaxation (sponge/damping) source term. See sourcemethod_LinearEuler2D_t for the mathematical description; the kernel reads the relaxation rate from solution variable 6 and never leaves the device.

Hu (2001) unsplit PML source term, evaluated per node. In the interior (sigma_x = sigma_y = 0) this leaves the acoustic variables untouched and the auxiliaries integrate q in time but never re-enter the dynamics (since the coupling coefficient sigma_xsigma_y is zero). Inside the PML, the (sigma_x + sigma_y) damping plus the sigma_xsigma_y*phi term produce the correct perfectly-matched behaviour for the linear Euler 2D system.

Spatially-dependent linear relaxation (sponge / damping layer) source.

Device-resident spatially-dependent linear relaxation (sponge/damping) source term. See sourcemethod_LinearEuler3D_t for the mathematical description; the kernel reads the relaxation rate from solution variable 7 and never leaves the device.

Spatially-dependent linear relaxation (sponge / damping layer) source.

This subroutine sets the initial condition for a weak blast wave problem. The initial condition is given by

This subroutine sets the initial condition for a weak blast wave problem. The initial condition is given by

Stage the pre-regrid solution on the DEVICE (Stage 6a), replacing the base implementation's device-to-host copy of the whole field. Regrid subsequently releases solution%interior_gpu, so the field is copied device-to-device into a buffer this model owns; ApplyTransferPlan then reads it from there.

Preserve the current solution ahead of a regrid, so that Regrid may release the storage it lives in. Pair with ApplyTransferPlan, which consumes the staged copy:

Stage the pre-regrid solution on the DEVICE, replacing the base implementation's device-to-host copy of the whole field. Pair with ApplyTransferPlan, which consumes the staged copy:

Preserve the current solution ahead of a regrid, so that Regrid may release the storage it lives in. Pair with ApplyTransferPlan, which consumes the staged copy:

Isoparametric subdivision of one element's geometry node coordinates into its four children. The parent geometry is the degree-nGeo Lagrange interpolant through parentCoords (sampled at geomInterp's control points, i.e. the mesh geometry nodes); each child node coordinate is that interpolant evaluated at the corresponding point of the parent reference square. Exact for any polynomial geometry of degree <= nGeo and for any control-node type, so straight-sided and curved (isoparametric) elements are both handled without approximation.

Isoparametric subdivision of one element's geometry node coordinates into its eight children. The parent geometry is the degree-nGeo Lagrange interpolant through parentCoords (sampled at geomInterp's control points, i.e. the mesh geometry nodes); each child node coordinate is that interpolant evaluated at the corresponding point of the parent reference cube. Exact for any polynomial geometry of degree <= nGeo and for any control-node type, so straight-sided and curved (isoparametric) elements are both handled without approximation.

Arithmetic-mean two-point flux for linear advection. Entropy-conserving with respect to eta(u) = u^2/2.

Entropy-conserving two-point flux function.

Souza et al. (2023, JAMES) entropy-conservative two-point flux for 2-D compressible Euler in (rho, rhov, rhotheta) variables with p = p0(rhoRd*theta/p0)^gamma.

Arithmetic-mean two-point flux for linear advection. Entropy-conserving with respect to eta(u) = u^2/2.

Entropy-conserving two-point flux function (3-D). flux(ivar, d) is the d-th physical component (d=1:x, d=2:y, d=3:z). Override in concrete models. Stub returns zero.

Souza et al. (2023, JAMES) entropy-conservative two-point flux for compressible Euler in (rho, rhov, rhotheta) variables with p = p0(rhoRd*theta/p0)^gamma.

Contravariant EC two-point flux on GPU — fully device-resident.

Contravariant EC two-point flux on GPU — fully device-resident.

Computes pre-projected SCALAR contravariant two-point fluxes for all node pairs and stores them in twoPointFlux%interior(n,i,j,iel,ivar,r).

Computes pre-projected SCALAR contravariant two-point fluxes for all node pairs, following Trixi.jl for curved meshes. Each direction r uses the correct partner AND the correct averaged metric Ja^r.

Souza et al. (2023) entropy-conservative two-point flux on GPU.

Pre-projected scalar contravariant two-point Souza et al. (2023) EC flux. For each node pair (a, b) along reference direction r:

Souza et al. (2023) entropy-conservative two-point flux on GPU. Fully device-resident.

Pre-projected scalar contravariant two-point Souza et al. (2023) EC flux. For each node pair (a, b) along reference direction r:

Create a fully triply-periodic structured hexahedral mesh and store it in SELF's unstructured mesh format. Element geometry and ordering are identical to UniformStructuredMesh; the only difference is connectivity. The faces on the six domain boundaries are wired as interior faces whose neighbor is the element on the opposite side of the domain. This realises the triply periodic box T^3 = [0,Lx] x [0,Ly] x [0,Lz] required by, e.g., the Arnold-Beltrami-Childress (ABC) flow benchmark.

\addtogroup SELF_SupportRoutines @{ \fn UniformPoints Generates a REAL(prec) array of N points evenly spaced between two points.

Produce meshOut, a uniformly refined copy of meshIn: every element is split into eight children. Child geometry is generated by exact isoparametric subdivision of the parent geometry (so the refined mesh represents the identical curved domain), and connectivity / flips are inherited deterministically from the base mesh. meshOut is a fully-formed, conforming Mesh3D_t with 8x the elements, ready for geometry generation and time stepping. Refined globalNodeIDs are left zero: the solver reads flips from sideInfo directly and no 3-D consumer requires them.

Produce meshOut, a uniformly refined copy of meshIn: every element is split into four children. Child geometry is generated by exact isoparametric subdivision of the parent geometry (so the refined mesh represents the identical curved domain), and connectivity / flips are inherited deterministically from the base mesh. meshOut is a fully-formed, conforming Mesh2D_t with 4x the elements, ready for geometry generation and time stepping.

Create a structured mesh and store it in SELF's unstructured mesh format. The mesh is created in tiles of size (tnx,tny). Tiling is used to determine the element ordering.

Create a structured mesh and store it in SELF's unstructured mesh format. The mesh is created in tiles of size (tnx,tny,tnz). Tiling is used to determine the element ordering.

Copy elements, coordinates, and the per-point Lagrange basis cache from host to device. Lazily (re)allocates the cache buffers if their degree changed since the previous call.

Computes a solution update as , where dt is either provided through the interface or taken as the Model's stored time step size (model % dt)

Computes a solution update as , where dt is either provided through the interface or taken as the Model's stored time step size (model % dt)