mpas

MPAS unstructured spherical centroidal-Voronoi (icosahedral) C-grid on the unit sphere: NCELLS Voronoi cells, NEDGES edges, NVERTICES triangle vertices, with cell valence bounded by MAXEDGES (on a closed sphere the Euler relations give NEDGES = 3*(NCELLS - 2) and NVERTICES = 2*(NCELLS - 2), but the three counts stay independent metaparameters so planar or open meshes can reuse the file). Unlike the rectangular grids there is no coordinate axis: the iteration domains are the three entity interval sets (cells, edges, vertices) plus the ragged connectivity set edges_of_cell (esm-schema IndexSet kind “ragged”; RFC semiring-faq-unified-ir 5.2), which binds by NAME to two keyed factors the consuming document must expose as cell-shaped arrays: nEdgesOnCell (offsets: the valence of each cell) and edgesOnCell (values: the edge ids around each cell, counterclockwise, zero-padded to maxEdges). Geometry is likewise DATA, not formula: mesh files under grids/mpas/mesh/ (section-4.7 subsystem files with literal const arrays) supply the factors nEdgesOnCell, edgesOnCell, edgeSignOnCell, cellsOnEdge, verticesOnEdge, dvEdge, dcEdge, areaCell, latCell, lonCell – and, since the 2026-07 regridding wave WIDENED the mesh contract, the VERTEX-GEOMETRY factors latVertex, lonVertex (vertex = Delaunay-triangle circumcenter, radians, lon in (-pi, pi], same convention as latCell/lonCell) and verticesOnCell (vertex ids around each cell, counterclockwise viewed from outside, zero-padded to maxEdges exactly like edgesOnCell; verticesOnCell[i][k] is the vertex SHARED by consecutive CCW edges edgesOnCell[i][k] and edgesOnCell[i][k+1 cyclic]; on a closed mesh every cell has exactly as many vertices as edges, so nEdgesOnCell doubles as the per-cell vertex count and there is no separate offsets factor) – under exactly those names; a consuming model mounts one as a subsystem and re-exposes the factors a rule needs as bare-name observed aliases (e.g. areaCell := mesh.areaCell), since keyed factors resolve by bare name in the model scope (see problems/divergence_mpas.esm for the pattern). The widened contract makes CELL POLYGON RINGS derivable in-document: the match-less template mpas_cell_rings_deg gathers each cell’s vertex ring into the [cells, maxEdges, ring_coord] lon-lat-DEGREE stack the conservative regridder consumes as src_poly/tgt_poly (regridding/conservative_overlap.esm; ring vertex = [lon, lat] in degrees, CCW), clamping the gather with min(k, nEdgesOnCell(i)) so PADDING SLOTS REPEAT THE CELL’S LAST REAL VERTEX – rectangular (padded) ring storage rather than a ragged ring, chosen because the polygon_intersection_area kernels accept trailing-duplicate padded rings value-identically (verified empirically against the Julia planar and spherical clip kernels: padded and unpadded rings agree at 0 ulp as subject and as clip operand; on a mesh whose MAXEDGES equals every cell’s valence, e.g. the level-0 pentagons, no padding ever occurs). The other closed-form geometry is the resolution scale: cell_area_mean = 4*pi/NCELLS, the mean cell area on the unit sphere (the h^2 measure for convergence sweeps, h ~ sqrt(cell_area_mean)). Orientation conventions are fixed by the mesh data (positive edge normal from cellsOnEdge[e][1] to cellsOnEdge[e][2]; positive tangent v1->v2 with t = k x n; edgeSignOnCell = +1 where the positive normal leaves the cell), per Ringler et al. (2010).

Metaparameters

All sizes are load-time metaparameters (esm-spec §9.7.6): rebind them at the import edge (or through the loader API) and the same files serve every resolution.

Consumer-supplied free names

areaCell cellsOnEdge dcEdge dvEdge edgeSignOnCell edgesOnCell nEdgesOnCell

The mesh geometry arrives as data, not formula: these are the keyed-factor arrays a consuming model must expose under exactly these bare names, normally as observed aliases of a mounted grids/mpas/mesh/ subsystem’s const data (see problems/divergence_mpas.esm). The 2026-07 regridding wave widened the contract with the vertex-geometry factors so cell polygon rings derive in-document via mpas_cell_rings_deg.

Index sets

Geometry templates

Rules

Complete auto-applied rewrite rules on spatial D: imported stencil + boundary conditions in one makearray (esm-spec §9.6.8). Import a rule and every matching derivative in your model lowers through it.

fv_divergence_cell

op:div order:1 bc:none axes:cells

Source: grids/mpas/rules/fv_divergence_cell.esm

Rewrites div(F) — match pattern:

{"op": "div", "args": ["F"]}

Definition:

Finite-volume cell divergence of an edge-normal flux on the MPAS Voronoi grid (the TRiSK discrete divergence, Ringler et al. 2010 eq. 21): div_i = (1/A_i) * sum_{k=1..nEdgesOnCell(i)} edgeSignOnCell(i,k) * F(e) * dvEdge(e) with e = edgesOnCell(i,k) — the exact Gauss/divergence theorem on cell i with the flux integral along each edge approximated by the midpoint rule, F being the flux component along the edge’s positive normal (which points from cellsOnEdge(e,1) to cellsOnEdge(e,2)) and edgeSignOnCell flipping it outward from cell i. Matches the open-namespace rewrite-target {op: div, args: [F]} (esm-spec 9.6.8: div is author sugar with no evaluable core, and a rule MAY match it directly and emit a one-shot scheme), match-scoped with where: {F: {shape: [edges]}} (esm-spec 9.6.1) so it fires ONLY on the divergence of a bare field declared over this grid’s edges index set — the axis-less div pattern would otherwise fire on every div in a consuming component. The shape constraint is the mechanism that lets two MPAS meshes coexist in one document: import this rule twice under distinct index-set renames (esm-spec 9.7.7), and each instance’s where shape follows its renamed edge set (e.g. [meshA.edges] vs [meshB.edges]), so each div lowers by its own mesh’s rule with no first-declared-wins tie-break. Because the constraint requires a bare variable of exactly that shape, a consumer that divergences a COMPOUND inline flux (e.g. div(u*h)) must first bind the flux to a declared edge-shaped observed (the divergence_mpas problem’s flux observed is exactly that pattern). The body is a single aggregate over the cells axis whose contracted index k runs over the RAGGED index set edges_of_cell (esm-spec 4.3.1 { from, of }; RFC semiring-faq-unified-ir 5.2): per-cell dynamic bound k <= nEdgesOnCell(i) plus a gather through the edgesOnCell values factor — this is the connectivity-based stencil form that replaces the makearray region form of the rectangular grids. The division by areaCell(i) sits inside the summand, which equals (1/A_i)sum since A_i is constant over k. There are no boundary regions and no makearray because a closed sphere has no boundary: every edge separates exactly two cells (bc:none; a rule for a limited-area MPAS mesh would instead overwrite the boundary-cell entries, and would be a different rule). The free names nEdgesOnCell, edgesOnCell, edgeSignOnCell, dvEdge, areaCell are the grid’s keyed-factor contract (see grids/mpas/grid.esm): the consuming model must expose them as bare-name arrays, normally as observed aliases of a grids/mpas/mesh/ subsystem’s const data. Accuracy: exactly conservative for any F (sum_i A_idiv_i = 0, each edge cancels between its two cells) and exactly annihilates discrete curls (div of F_e = (eta(v2)-eta(v1))/dvEdge(e) is zero to rounding — the TRiSK div-curl identity; problems/divergence_mpas.esm asserts it); order tag is the guaranteed first-order bound on arbitrarily distorted Voronoi meshes, though the midpoint flux is second-order on centroidal/symmetric meshes (Ringler et al. 2010, section 3.2).

References

• Ringler, T., Thuburn, J., Klemp, J. B., & Skamarock, W. C. (2010). A unified approach to energy conservation and potential vorticity dynamics for arbitrarily-structured C-grids. Journal of Computational Physics, 229(9), 3065-3090. Eq. (21) and the div-curl mimetic identities of section 3.
• Thuburn, J., Ringler, T. D., Skamarock, W. C., & Klemp, J. B. (2009). Numerical representation of geostrophic modes on arbitrarily structured C-grids. Journal of Computational Physics, 228(22), 8321-8335.

Conformance fixtures

• ast: fv_divergence_cell — status active
• convergence: divergence_mpas — status active
• simulation: divergence_mpas — status active

Convergence — case divergence_mpas

Error norms read from the committed golden tests/conformance/convergence/divergence_mpas/golden/errors.json (binding: julia; nothing recomputed for display).

Expected order: 1 (± 0.25).

fv_gradient_edge

op:grad order:1 bc:none axes:edges

Source: grids/mpas/rules/fv_gradient_edge.esm

Rewrites grad(phi) — match pattern:

{"op": "grad", "args": ["phi"]}

Definition:

Finite-volume edge-normal gradient of a cell scalar on the MPAS Voronoi grid (the TRiSK discrete gradient, Ringler et al. 2010 eq. 22): grad_e = (phi(cellsOnEdge(e,2)) - phi(cellsOnEdge(e,1))) / dcEdge(e) — the component, along each edge’s positive normal, of the gradient of the cell field phi, formed as the two-point difference of the flanking cell values divided by the center-to-center arc dcEdge(e). The positive normal points from cellsOnEdge(e,1) to cellsOnEdge(e,2) (the MPAS orientation shared with edgeSignOnCell and with the divergence rule’s flux convention), so the difference is taken head-minus-tail in that direction and dcEdge is the great-circle distance the difference spans. This is the exact MIMETIC DUAL of grids/mpas/rules/fv_divergence_cell: that rule maps an edge-normal flux to a cell divergence over cells; this rule maps a cell scalar to an edge-normal gradient over edges, and the two are negative adjoints of one another (see Accuracy). Matches the open-namespace rewrite-target {op: grad, args: [phi]} (esm-spec 9.6.8, 4.2: grad is author sugar with no evaluable core, and a rule MAY match it directly on the bare op and emit a one-shot scheme rather than routing through the D-sugar identity grad(f,dim:x)=D(f,wrt:x)), match-scoped with where: {phi: {shape: [cells]}} (esm-spec 9.6.1) so it fires ONLY on the gradient of a bare field declared over this grid’s cells index set — the axis-less grad pattern would otherwise fire on every grad in a consuming component. The shape constraint is the mechanism that lets two MPAS meshes coexist in one document: import this rule twice under distinct index-set renames (esm-spec 9.7.7), and each instance’s where shape follows its renamed cell set (e.g. [meshA.cells] vs [meshB.cells]), so each grad lowers by its own mesh’s rule with no first-declared-wins tie-break. Because the constraint requires a bare variable of exactly that shape, a consumer that takes the gradient of a COMPOUND inline cell expression (e.g. grad(a*b)) must first bind it to a declared cell-shaped observed. The body is a single aggregate over the edges axis with output index e and NO contracted index — a pure connectivity map, one two-point stencil per edge — because on the MPAS C-grid the gradient is intrinsically edge-local: every edge separates exactly two cells and the edge-normal derivative reads only those two, unlike the divergence whose body contracts the ragged edges_of_cell set. There is no makearray and no boundary region because a closed sphere has no boundary: every edge has both a cellsOnEdge(e,1) and a cellsOnEdge(e,2) (bc:none; a rule for a limited-area MPAS mesh would instead special-case the boundary edges that have only one interior cell, and would be a different rule). The free names cellsOnEdge, dcEdge are the grid’s keyed-factor contract (see grids/mpas/grid.esm): the consuming model must expose them as bare-name arrays, normally as observed aliases of a grids/mpas/mesh/ subsystem’s const data (problems/gradient_mpas.esm is the pattern). Accuracy: exactly annihilates spatially-constant fields (phi(c2)-phi(c1) = 0 for phi = const, hence grad_e = 0 at every edge to rounding — problems/gradient_mpas.esm asserts this), which is the gradient’s null-space mimetic property and the discrete statement that grad(const) = 0. The rule is the negative discrete adjoint of fv_divergence_cell under the TRiSK inner products <a,b>_cells = sum_i areaCell(i)a_ib_i and <u,v>_edges = sum_e dvEdge(e)*dcEdge(e)u_ev_e: summation by parts on the closed sphere carries no boundary term, so sum_i areaCell(i)phi_idiv(F)_i = -sum_e dvEdge(e)*dcEdge(e)grad(phi)_eF_e for any cell field phi and edge field F, the discrete Gauss/Green identity the two operators share (Ringler et al. 2010, section 3). The order tag is the guaranteed first-order bound on arbitrarily distorted Voronoi meshes, though the two-point difference is second-order on centroidal/symmetric meshes where the edge midpoint bisects the cell-center arc (Ringler et al. 2010, section 3.2). A convergence-order sweep on refined MPAS meshes is a possible follow-up, out of scope for the exactness gate that ships with this rule.

References

• Ringler, T., Thuburn, J., Klemp, J. B., & Skamarock, W. C. (2010). A unified approach to energy conservation and potential vorticity dynamics for arbitrarily-structured C-grids. Journal of Computational Physics, 229(9), 3065-3090. Eq. (22) and the discrete Green/adjoint (grad = -div^T) identities of section 3.
• Thuburn, J., Ringler, T. D., Skamarock, W. C., & Klemp, J. B. (2009). Numerical representation of geostrophic modes on arbitrarily structured C-grids. Journal of Computational Physics, 228(22), 8321-8335.

Conformance fixtures

• ast: fv_gradient_edge — status active
• simulation: gradient_mpas — status active

fv_laplacian_cell

op:laplacian order:1 bc:none axes:cells

Source: grids/mpas/rules/fv_laplacian_cell.esm

Rewrites laplacian(phi) — match pattern:

{"op": "laplacian", "args": ["phi"]}

Definition:

Finite-volume cell Laplacian of a cell scalar on the MPAS Voronoi grid (the TRiSK discrete Laplacian, Ringler et al. 2010, the div-of-grad composition of eqs. 21-22): lap_i = (1/A_i) * sum_{k=1..nEdgesOnCell(i)} edgeSignOnCell(i,k) * dvEdge(e) * (phi(cellsOnEdge(e,2)) - phi(cellsOnEdge(e,1))) / dcEdge(e) with e = edgesOnCell(i,k) – the exact Gauss/divergence theorem on cell i applied to the edge-normal gradient of phi, i.e. div(grad phi) with the flux integral along each edge approximated by the midpoint rule. This rule is the mimetic COMPOSITION of grids/mpas/rules/fv_divergence_cell with grids/mpas/rules/fv_gradient_edge: the divergence’s per-edge flux F(e) is replaced INLINE by the two-point edge-normal gradient grad_e = (phi(c2) - phi(c1)) / dcEdge(e) that fv_gradient_edge would emit (c1 = cellsOnEdge(e,1), c2 = cellsOnEdge(e,2); the positive edge normal points c1 -> c2, the MPAS orientation shared by edgeSignOnCell and by both parent rules), so the summand is the outward-signed edge-normal gradient flux edgeSignOnCell(i,k)dvEdge(e)grad_e and the sum is the net flux through cell i’s boundary, per unit area A_i. Emitting the composed one-shot stencil (rather than lowering div and grad on separate fixpoint passes) is the open-namespace option esm-spec 9.6.8 explicitly grants: ‘A rule author MAY instead match laplacian directly and emit a one-shot stencil.’ Matches the open-namespace rewrite-target {op: laplacian, args: [phi]} (esm-spec 9.6.8, 4.2: laplacian is author sugar with no evaluable core – laplacian(u) = sum_i D(D(u,x_i),x_i) – and a rule MAY match it directly on the bare op and emit a scheme rather than routing through the nested-D identity), match-scoped with where: {phi: {shape: [cells]}} (esm-spec 9.6.1) so it fires ONLY on the Laplacian of a bare field declared over this grid’s cells index set – the axis-less laplacian pattern would otherwise fire on every laplacian in a consuming component. The shape constraint is the mechanism that lets two MPAS meshes coexist in one document: import this rule twice under distinct index-set renames (esm-spec 9.7.7), and each instance’s where shape follows its renamed cell set (e.g. [meshA.cells] vs [meshB.cells]), so each laplacian lowers by its own mesh’s rule with no first-declared-wins tie-break. Because the constraint requires a bare variable of exactly that shape, a consumer that takes the Laplacian of a COMPOUND inline cell expression (e.g. laplacian(ab)) must first bind it to a declared cell-shaped observed. The body is a single aggregate over the cells axis whose contracted index k runs over the RAGGED index set edges_of_cell (esm-spec 4.3.1 { from, of }; RFC semiring-faq-unified-ir 5.2): per-cell dynamic bound k <= nEdgesOnCell(i) plus a gather through the edgesOnCell values factor – the same connectivity-based stencil form as fv_divergence_cell, since the Laplacian IS a divergence (its output is a cell aggregate contracting the cell’s ragged edge ring), only with the edge flux built in-place from phi and dcEdge instead of read from a declared edge field. The division by areaCell(i) sits inside the summand, which equals (1/A_i)sum since A_i is constant over k, exactly as in fv_divergence_cell. There are no boundary regions and no makearray because a closed sphere has no boundary: every edge separates exactly two cells, so every cellsOnEdge(e,1)/cellsOnEdge(e,2) gather and every edgeSignOnCell entry the ragged bound reaches is well defined (bc:none; a rule for a limited-area MPAS mesh would instead overwrite the boundary-cell entries – one-sided/Neumann/Dirichlet closure of the edges with only one interior cell – and would be a different rule). The free names edgeSignOnCell, edgesOnCell, dvEdge, cellsOnEdge, dcEdge, areaCell are the grid’s keyed-factor contract (see grids/mpas/grid.esm): the consuming model must expose them as bare-name arrays, normally as observed aliases of a grids/mpas/mesh/ subsystem’s const data (problems/laplacian_mpas.esm is the pattern) – the union of fv_divergence_cell’s factors (edgeSignOnCell, edgesOnCell, dvEdge, areaCell) and fv_gradient_edge’s factors (cellsOnEdge, dcEdge), since the composed stencil needs both the divergence’s cell-boundary connectivity and the gradient’s flanking-cell gather. Accuracy: exactly ANNIHILATES spatially-constant fields (phi(c2)-phi(c1) = 0 for phi = const at every edge, hence every summand edgeSignOnCelldvEdge0/dcEdge = 0 and lap_i = 0 at every cell to rounding – the Laplacian’s null-space mimetic property, the discrete statement that laplacian(const) = 0; problems/laplacian_mpas.esm asserts it). Exactly CONSERVATIVE for any cell field phi (sum_i A_ilap_i = 0: each interior edge e appears in exactly its two cells cellsOnEdge(e,1) and cellsOnEdge(e,2) with opposite edgeSignOnCell and the same dvEdgegrad_e, so the two cells’ contributions cancel and the global integral of the discrete Laplacian is zero – the discrete divergence-theorem closure inherited from fv_divergence_cell; problems/laplacian_mpas.esm asserts it for a non-constant field). And SYMMETRIC NEGATIVE-SEMIDEFINITE under the TRiSK inner products <a,b>_cells = sum_i areaCell(i)a_ib_i and <u,v>_edges = sum_e dvEdge(e)*dcEdge(e)u_ev_e: because lap = div(grad) and fv_gradient_edge is the negative discrete adjoint of fv_divergence_cell (grad = -div^T under those inner products, the closed-sphere summation-by-parts identity of section 3), the Laplacian satisfies <phi, lap psi>_cells = -sum_e dvEdge(e)*dcEdge(e)grad(phi)_egrad(psi)_e, which is symmetric in phi<->psi and equals -||grad phi||^2 <= 0 at psi = phi – the discrete Dirichlet energy, the property that makes the TRiSK Laplacian a stable diffusion operator. The order tag is the guaranteed first-order bound on arbitrarily distorted Voronoi meshes, though the composed midpoint flux is second-order on centroidal/symmetric meshes where the edge midpoint bisects the cell-center arc (Ringler et al. 2010, section 3.2). A convergence-order sweep on refined MPAS meshes is a possible follow-up, out of scope for the exactness gate that ships with this rule.

References

• Ringler, T., Thuburn, J., Klemp, J. B., & Skamarock, W. C. (2010). A unified approach to energy conservation and potential vorticity dynamics for arbitrarily-structured C-grids. Journal of Computational Physics, 229(9), 3065-3090. Eqs. (21)-(22) and the discrete div-grad adjoint (grad = -div^T) / Dirichlet-energy identities of section 3.
• Thuburn, J., Ringler, T. D., Skamarock, W. C., & Klemp, J. B. (2009). Numerical representation of geostrophic modes on arbitrarily structured C-grids. Journal of Computational Physics, 228(22), 8321-8335.

Conformance fixtures

• ast: fv_laplacian_cell — status active
• simulation: laplacian_mpas — status active