Finite element bases

Lagrange, BDM and RT bases, tangential trace coordinates and their declared local degrees of freedom.

All API families

Basix reference-element tabulation with explicit executed basis data.

Basix supplies polynomial bases, derivatives and entity transformations. Its native family/degree convention is retained: for example Basix RT degree one is the lowest-order RT element. Reference capabilities do not establish the compatibility, physical conformity or error estimates of a hybrid formulation. The native library is loaded when tabulation is requested; algebraic providers do not require a native element handle.

ReferenceElementSpec dataclass

Declare an element using native Basix family, cell and degree names.

Basix checks supported family/cell/order/variant combinations; this record introduces no degree shifts or new admissibility hypotheses. Family aliases are resolved by Basix's own string_to_family. dof_ordering follows Basix's native old-to-new permutation convention: index j contains the executed slot of canonical DOF j. An empty tuple requests Basix's default order, as does None. Supported permutation/family choices are checked by Basix. Elements use binary64. Descriptors are picklable; create the native element separately inside each worker.

__post_init__() -> None

Freeze the declared permutation and reject malformed specifications.

ReferenceElement dataclass

A process-local native element and read-only executed basis metadata.

basis_matrix contains FiniteElement.coefficient_matrix after native dualization, with rows reordered into the actual executed DOF order, rather than wcoeffs defining only a polynomial span. Its digest identifies the literal binary64 matrix. Replaying coefficients also requires the declared cell, family, variants, order, polyset and backend version. Native handles remain inside their worker; persist metadata and actual basis arrays rather than serializing this object across processes.

NodalReferenceBasis dataclass

An equispaced scalar basis in a caller's literal node order.

permutation[j] identifies the native basis function for caller node j. basis_matrix and its digest capture the actually reordered coefficient rows. The native reference basis and declared permutation must be preserved together when persisting coefficients. Arrays are read-only cache data.

create_reference_element(spec: ReferenceElementSpec) -> ReferenceElement

Create or reuse a read-only Basix element with its native conventions.

Supported shapes, families, variants and orders are inherited from the installed Basix version. Unsupported combinations raise its creation error. The cache is process-local and synchronized during creation; tabulation does not modify the native element or cached basis metadata.

tabulate_reference(element: ReferenceElement, points: ArrayLike, nderiv: int = 0) -> FloatArray

Evaluate native basis values and Cartesian derivatives through nderiv.

The binary64 output has axes (derivative, point, basis, value_component). Scalar elements retain a final axis of size one; vector/tensor values use Basix's flattened component order. Derivatives follow basix.index; any nonnegative derivative order supported by Basix is allowed. Coordinates have shape (points, cell_dimension) and may lie outside the cell. No physical Piola map, cell orientation or global DOF identification is applied.

tabulate_archived_basis(cell: str, degree: int, basis_matrix: ArrayLike, points: ArrayLike, *, components: int = 1, nderiv: int = 0) -> FloatArray

Evaluate an executed scalar/vector polynomial basis from its archived matrix.

basis_matrix contains rows of the actual Basix coefficient_matrix in its native orthonormal Legendre polynomial coordinates. Its rows may include a declared change of basis. The returned array has axes (points, executed_basis) for a scalar and (points, executed_basis, components) for a vector for nderiv=0. Positive derivative orders add a leading derivative axis in Basix's index convention, including values at index zero. Evaluation uses the supplied matrix directly, without dualizing or rebuilding finite-element basis coefficients. Physical maps, continuity, vector blocking and global DOF order belong to the caller.

tabulate_archived_nodal_basis(cell: str, degree: int, basis_matrix: ArrayLike, points: ArrayLike, *, nderiv: int = 0) -> FloatArray

Evaluate the literal executed scalar basis in archived orthonormal coordinates.

Return (point,basis), or (derivative,point,basis) when nderiv is positive. No finite element or independent coefficient matrix is regenerated.

reference_interpolation_points(element: ReferenceElement) -> FloatArray

Return read-only native interpolation points in reference coordinates.

interpolate_reference(element: ReferenceElement, values: ArrayLike) -> FloatArray

Apply native value interpolation to several candidate functions.

Input axes are (interpolation_point, candidate, flattened_value_component). The output has axes (executed_dof, candidate), including a declared DOF permutation. Only value-based functionals are accepted; derivative-based elements require derivative data and are outside this helper's contract. No physical Piola map or entity orientation is inferred.

reference_entity_transformations(element: ReferenceElement) -> dict[str, FloatArray]

Return read-only native entity maps, including H(div)/H(curl) signs.

Matrices act in the intrinsic DOF order of each entity, even when a global dof_ordering is declared. :func:reference_entity_dofs gives those local DOFs' executed indices; embed these maps using those indices instead of applying a global permutation to a small entity matrix. Entity names retain Basix's reference orientation convention. They do not infer PyMHM macroface normals, mesh incidence or a physical trace map.

reference_entity_dofs(element: ReferenceElement) -> tuple[tuple[tuple[int, ...], ...], ...]

Return each entity's intrinsic local order as executed DOF indices.

Axes are (entity_dimension, entity_number, entity_local_dof). Basix already applies the declared old-to-new global permutation to these indices. Pair them with its unchanged entity-local transformations.

reference_base_transformations(element: ReferenceElement) -> FloatArray

Return full DOF maps in the actual executed basis order.

Axes are (transformation, output_dof, input_dof). Basix supplies edge reversals and face rotations/reflections in canonical order. A declared dof_ordering conjugates those maps into the tabulation/basis order. No physical or global trace orientation is selected by this function.

simplex_lagrange_basis(cell: str, degree: int, *, nodes: ArrayLike) -> NodalReferenceBasis

Declare a native equispaced Pk basis using literal barycentric node order.

Positive degree is unrestricted subject to Basix's implementation. Node coordinates must contain its complete equispaced lattice exactly once; this changes only numbering, never solves for a new polynomial basis or merges global nodes by geometric proximity. High-order equispaced nodal interpolation retains its inherent conditioning limitations.

tensor_lagrange_basis(cell: str, degree: int, *, nodes: ArrayLike) -> NodalReferenceBasis

Declare a native equispaced Qk basis in explicit Cartesian node order.

Cells are interval, quadrilateral or hexahedron; degree is positive. Nodes are a complete lattice of shape (basis, reference_dimension). The returned coefficient matrix follows that order, for archive/replay.

tensor_lagrange_tabulation(cell: str, degree: int, points: ArrayLike, *, nodes: ArrayLike, nderiv: int = 1) -> tuple[FloatArray, FloatArray]

Return ordered native Qk values and Cartesian reference gradients.

Axes are (point, basis) and (point, basis, reference_dimension). Coordinates may lie outside the cell. Order zero returns an empty gradient coordinate axis; order one returns all Cartesian first derivatives. Mapping into physical coordinates is the caller's explicit geometric operation. Values at literally declared interpolation nodes are the exact Kronecker rows; derivatives retain the executed native tabulation.

nodal_base_transformations(basis: NodalReferenceBasis) -> FloatArray

Conjugate native DOF orientation maps into the declared nodal order.

simplex_lagrange_tabulation(cell: str, degree: int, bary: ArrayLike, *, nodes: ArrayLike, nderiv: int = 2) -> tuple[FloatArray, FloatArray, FloatArray]

Tabulate ordered simplex Pk values, Cartesian gradients and Hessians.

Reference coordinates are bary[:, 1:] with bary[:, 0] = 1-sum(x). Return shapes are (points, basis), (points, basis, dimension) and (points, basis, dimension, dimension). Only derivative orders 0/1/2 are returned here; unrequested gradient/Hessian coordinate axes are empty. These are derivatives on the unit-sum simplex, not derivatives in all independent barycentric variables. Physical derivatives require the caller's affine Jacobian, explicitly. Values at literally declared interpolation nodes are exact Kronecker rows; no proximity threshold or derivative correction is applied.

barycentric_simplex_tabulation(cell: str, degree: int, bary: ArrayLike, *, nodes: ArrayLike, nderiv: int = 2) -> tuple[FloatArray, FloatArray, FloatArray]

Express native simplex derivatives in the canonical barycentric extension.

On the unit-sum hyperplane, the basis is the ordered Basix Pk basis. Its extension is p(lambda_1,...,lambda_d), independent of lambda_0. Returned gradients and Hessians have d+1 coordinate axes; entries involving lambda_0 are zero. Contracting with physical barycentric gradients gives the usual affine derivatives. The extension itself has no physical role.

orthogonal_polynomial_tabulation(cell: str, degree: int, points: ArrayLike, nderiv: int = 0) -> FloatArray

Tabulate Basix's orthonormal polynomial set and Cartesian derivatives.

Values use native reference-cell measure and ordering, with axes (derivative, point, polynomial). Derivative indices follow basix.index. This supplies polynomial spans, rather than dualized finite elements.

legendre_tabulation(points: ArrayLike, degree: int, nderiv: int = 0) -> FloatArray

Evaluate conventional Legendre polynomials on [-1,1] through Basix.

L_j(1)=1; derivatives are with respect to the supplied coordinate, including extrapolation points. Axes are (derivative, *points.shape, degree+1). The affine coordinate map and sqrt(2j+1) normalization convert Basix's orthonormal interval polynomial set into this persisted moment convention.

legendre_values(points: ArrayLike, degree: int) -> FloatArray

Return conventional L_0,...,L_degree at finite real coordinates.

monomial_coefficients(degree: int) -> FloatArray cached

Map orthonormal interval coefficients to the declared power convention.

The coefficients are exact integral moments of x^k against the shifted Legendre polynomials. This changes polynomial coordinates only; Basix owns every repeated value and derivative evaluation.

monomial_tabulation(points: ArrayLike, powers: tuple[tuple[int, ...], ...], *, nderiv: int = 1) -> FloatArray

Tabulate declared component powers using Basix interval polynomial factors.

Preserve the supplied exponent order, including incomplete polynomial subspaces used for MHM moment definitions and archived coefficient vectors. Axes are (derivative, point, monomial); order one contains the value followed by each Cartesian derivative. Any finite real coordinates are allowed. The product construction works in arbitrary coordinate dimension, without imposing a simplex or tensor-cell approximation space on the caller.

bernstein_basis_matrix(dimension: int, degree: int, exponents: tuple[tuple[int, ...], ...]) -> FloatArray

Return the literal executed Bernstein matrix in a declared barycentric order.

Rows are only permuted from the native matrix used by bernstein_tabulation; columns retain Basix's orthonormal reference polynomial coordinates. No interpolation fit or new basis is generated to describe these coefficients.

bernstein_tabulation(points: FloatArray, exponents: tuple[tuple[int, ...], ...]) -> tuple[FloatArray, FloatArray]

Tabulate simplex Bernstein polynomials in an explicit barycentric order.

exponents lists all barycentric multiindices of one degree. The native basis is only permuted; no polynomial fit or numerical basis rotation is performed. Cartesian gradients have axes (point, polynomial, coordinate).

physical_simplex_tabulation(cell: str, degree: int, bary: ArrayLike, *, nodes: ArrayLike, reference_gradients: ArrayLike, nderiv: int = 2) -> tuple[FloatArray, FloatArray, FloatArray]

Map ordered native Pk derivatives using explicitly supplied affine geometry.

reference_gradients has shape (cells, dimension, physical_dimension) and contains the physical gradients of lambda_1 through lambda_dimension. It is the caller's affine Jacobian data; no geometry or orientation is inferred. bary is either shared (points, dimension+1) or per-cell (cells, points, dimension+1). Values retain that shared/per-cell layout; derivatives always have leading cell and point axes. Unrequested derivative coordinate axes are empty. Hessians use the affine chain rule; curved cells require additional geometric derivatives and are outside this adapter.

Discrete inverse constants in physical finite-element coordinates.

laplacian_inverse_bound(gradient: FloatArray, hessian: FloatArray, weights: FloatArray, diameters: FloatArray) -> FloatArray

Return cellwise m=min(1/3,C) with C*h²*||Delta v||²<=||grad v||².

Gradients have axes (cell, point, basis, physical_direction) and Hessians append a second physical direction. Positive quadrature weights integrate the reference cell; its physical volume cancels in the energy quotient. diameters are physical cell diameters. The basis must span the constants with exactly one resolved stiffness kernel per cell.

The spectral quotient applies in two and three dimensions. It describes this supplied polynomial space and quadrature, rather than a mesh-uniform inf-sup certificate. The quadrature must integrate both Gram matrices exactly. Unresolved additional stiffness modes are rejected using the declared relative spectral threshold 1e-12. Piecewise-affine functions have zero Laplacian and return 1/3.

Continuous nodal triangle topology and Basix Pk tabulation.

multiindices(degree: int) -> IntArray

Order barycentric lattice nodes by vertices, oriented edges and interior.

reference_values(degree: int, bary: FloatArray) -> FloatArray

Evaluate native equispaced triangular Pk values in PyMHM node order.

bary is a finite real array of shape (points, 3); evaluation is allowed outside the reference triangle on the unit-sum barycentric hyperplane. Coordinates and returned values use binary64, with the same node order as :func:multiindices. Basix tabulates values without allocating derivatives.

reference_basis(degree: int, bary: FloatArray) -> tuple[FloatArray, FloatArray, FloatArray]

Return Basix Pk values and derivatives in the canonical barycentric extension.

Shapes are (points, basis), (points, basis, 3) and (points, basis, 3, 3). The extension is p(lambda_1,lambda_2): lambda_0 derivatives are zero, while the other entries are Cartesian reference derivatives. Contract with all three physical barycentric gradients to obtain physical derivatives. Input coordinates must sum to one; off-triangle points on that hyperplane are supported. The node order remains :func:multiindices; independent off-hyperplane extensions are not part of the finite-element field.

nodal_space(mesh: TriangleMesh, degree: int) -> tuple[IntArray, FloatArray]

Enumerate continuous nodal unknowns without merging geometrically close nodes.

tabulate(mesh: TriangleMesh, degree: int, bary: FloatArray, *, backend: Literal['portable', 'basix'] = 'basix') -> tuple[IntArray, FloatArray, FloatArray, FloatArray, FloatArray]

Return DOFs, coordinates, values, physical gradients and physical Hessians.

Basix supplies equispaced Pk values and Cartesian derivatives in :func:multiindices order. Points must sum to one. The affine geometry, continuous topology and nodal coefficient order are explicit PyMHM data. backend='portable' is a compatibility spelling for this same Basix execution; there is no separate polynomial implementation.

element_tabulate(mesh: TriangleMesh, degree: int, bary: FloatArray, *, backend: Literal['portable', 'basix'] = 'basix') -> tuple[IntArray, FloatArray, FloatArray, FloatArray, FloatArray]

Tabulate Pk at separate barycentric quadrature points in each triangle.

bary has shape (cells, points, 3). All returned basis arrays retain the leading cell axis, so material intersections do not require evaluating every element's basis at every other element's quadrature points. backend follows :func:tabulate; tabulation preserves the declared continuous topology and nodal coefficient order.

trace_coupling(coarse: TriangleMesh, cell: int, fine: TriangleMesh, skeleton: SkeletonSpace, degree: int) -> FloatArray

Integrate signed nodal traces over independently subdivided macrofaces.

scalar_operators(mesh: TriangleMesh, degree: int, *, diffusion: Any = 1.0, source: Any = 0.0, reaction: Any = 0.0, advection: Any = (0.0, 0.0), skew_advection: bool = False, order: int = 6, element_backend: Literal['portable', 'basix'] = 'basix') -> tuple[Any, Any, FloatArray]

Assemble continuous Pk scalar diffusion, reaction and advective operators.

With skew_advection=False the convection term is (beta.grad(u),v). Its antisymmetric part is selected when skew_advection=True. A variable conservative velocity also requires half its divergence in reaction. Basix tabulates local Pk values and derivatives in the same nodal coefficient order. element_backend='portable' is a compatibility spelling for this same execution.

Quadratic Brezzi–Douglas–Marini vectors on affine triangles.

The nine face degrees of freedom are integrals of the oriented normal component against Legendre polynomials of degrees zero through two. Three cell moments complete the element. Contravariant Piola transforms preserve normal moments and divergence; odd edge moments include parameter reversal.

bdm2_dofs(mesh: TriangleMesh) -> IntArray

Return twelve conforming vector DOFs per triangle, with cell moments last.

bdm2_basis(mesh: TriangleMesh, barycentric: FloatArray) -> tuple[FloatArray, FloatArray]

Evaluate global-oriented BDM2 bases and divergences in every triangle.

The outputs have shapes (cells, points, 12, 2) and (cells, points, 12). Barycentric coordinates have shape (points, 3) or (cells, points, 3) for separate material-intersection rules.

bdm2_evaluate(mesh: TriangleMesh, coefficients: FloatArray, barycentric: FloatArray) -> tuple[FloatArray, FloatArray]

Evaluate a BDM2 vector or row-wise tensor and its physical divergence.

Coefficients have shape (3*faces+3*cells,) for a vector or (3*faces+3*cells, rows) for a row-wise tensor. Tensor output axes are (cell, point, row, column); divergence retains the row axis.

validate_bdm2_trace(skeleton: SkeletonSpace, refinement: int) -> None

Require a representable piecewise quadratic trace on aligned fine edges.

bdm2_trace_map(mesh: TriangleMesh, cell: int, fine: TriangleMesh, skeleton: SkeletonSpace) -> FloatArray

Map scalar signed macro normal traces to fine boundary BDM2 moments.

The result has one row per boundary face and Legendre moment (degrees 0–2), and one column per scalar macroface basis. Vector traces reuse this matrix independently for each component. Breaks must align with fine edges, as checked by validate_bdm2_trace before local assembly.

Polynomial H(div) families with independently enriched interior triangle modes.

BDM(k,n) retains normal degree k and every interior mode of BDM(k+n). The choices n=0,1,2 realize the BDM, BDM-plus and BDM-double-plus families. All degrees of freedom are physical normal moments or reference cell moments; contravariant Piola mapping preserves their duality.

BDMFamily dataclass

Triangle vector polynomial family with normal degree k and interior degree k+n.

degree is at least one; enrichment is a nonnegative integer. The complete polynomial degree is their sum, divergence has degree one less, and the normal trace retains degree k. Each fine face has k+1 moments and each cell has (k+n)^2-1 interior moments. Removing higher face modes leaves the full zero-normal bubble subspace, rather than truncating component monomials.

polynomial_degree: int property

Return the complete interior vector polynomial degree.

interior_size: int property

Return the number of cell moments per vector.

local_size: int property

Return the total number of vector basis functions in one triangle.

reference_coefficients: FloatArray property

Return the executed native-BDM to normal/interior moment transform.

Rows follow the native Basix candidate basis; columns follow three face blocks of retained normal degree and every complete interior bubble. Persist this matrix together with the native coefficient matrix when replaying physical coefficients across numerical-library versions.

__post_init__() -> None

Validate the finite polynomial family without imposing a tested-degree cap.

size(mesh: TriangleMesh) -> int

Return the global conforming vector dimension on a local mesh.

dofs(mesh: TriangleMesh) -> IntArray

Return global-oriented face moments followed by private cell moments.

basis(mesh: TriangleMesh, barycentric: FloatArray) -> tuple[FloatArray, FloatArray]

Evaluate Piola vector bases and divergence at common or cellwise points.

Barycentric coordinates may have shape (q,3) or (cells,q,3), allowing separate material-intersection quadrature without changing the basis. Returned arrays have shape (cells,q,dofs,2) and (cells,q,dofs).

evaluate(mesh: TriangleMesh, coefficients: FloatArray, barycentric: FloatArray) -> tuple[FloatArray, FloatArray]

Evaluate a conforming vector or row-wise tensor and its divergence.

validate_trace(skeleton: SkeletonSpace, refinement: int) -> None

Require normal degrees no greater than k and fine-edge-aligned segments.

trace_map(mesh: TriangleMesh, cell: int, fine: TriangleMesh, skeleton: SkeletonSpace) -> FloatArray

Map signed scalar macro tractions to the oriented fine normal moments.

Raviart–Thomas vectors and their canonical moments on affine triangles.

RTField dataclass

Immutable canonical Raviart--Thomas coefficients for repeated point evaluation.

Construction validates and copies the complete real coefficient vector. Subsequent calls validate only their points and incident cells, and use the same Piola arithmetic as :func:rt_evaluate_points. This avoids scanning an entire global field for every small evaluation batch. The original mesh is retained; its geometry arrays follow :class:TriangleMesh's read-only contract.

__init__(mesh: TriangleMesh, coefficients: Any, degree: int) -> None

Copy coefficients into bytes-backed read-only storage, preserving precision.

__reduce__() -> tuple[Any, tuple[Any, ...]]

Reestablish immutable storage after serialization, including spawn workers.

evaluate(points: FloatArray, cell_indices: IntArray) -> tuple[FloatArray, FloatArray]

Evaluate one incident cell per point, preserving one-sided physical values.

rt_degree(degree: int) -> int

Validate a nonnegative mathematical Raviart--Thomas order.

rt_interior_tests(degree: int, points: FloatArray) -> FloatArray

Evaluate canonical P_(m-1) interior tests: monomials through RT2, orthonormal above.

rt_moment_coefficients(degree: int) -> FloatArray cached

Construct the canonical RT dual basis from exact reference moments.

rt_dofs(mesh: TriangleMesh, degree: int) -> IntArray

Return oriented face moments followed by cell-interior moments per triangle.

rt_basis(mesh: TriangleMesh, degree: int, barycentric: FloatArray) -> tuple[FloatArray, FloatArray]

Evaluate globally oriented RT_m basis vectors and physical divergence.

Mathematical order m corresponds to Basix/DOLFINx RT degree m+1. The outputs have shapes (cells, points, (m+1)(m+3), 2) and (cells, points, DOFs). Barycentric points may have shape (points,3) or (cells,points,3).

rt_evaluate(mesh: TriangleMesh, coefficients: FloatArray, degree: int, barycentric: FloatArray) -> tuple[FloatArray, FloatArray]

Evaluate an RT vector field and its divergence in every fine triangle.

rt_evaluate_points(mesh: TriangleMesh, coefficients: FloatArray, degree: int, points: FloatArray, cell_indices: IntArray) -> tuple[FloatArray, FloatArray]

Evaluate one specified incident RT cell at each physical point.

Points have shape (n,2), with an integer cell index for each point. No geometric nearest-neighbor guess or interface averaging is performed. Points must belong to their declared triangle, up to barycentric roundoff. Values and divergence use the same canonical moments and Piola orientation as :func:rt_evaluate; storage scales with the number of requested points.

rt_interpolate(mesh: TriangleMesh, field: Any, degree: int, order: int = 6) -> FloatArray

Interpolate a single-valued physical vector field using canonical RT moments.

This routine requires a well-defined normal trace. The MHM moment reconstruction supplies its own averaged one-sided interior-face data.

Basix scalar Pk topology, tabulation and affine tetrahedral operators.

tetrahedron_quadrature(order: int = 5) -> tuple[FloatArray, FloatArray]

Return positive Duffy-product quadrature with weights normalized to unit volume.

tetra_nodal_space(mesh: TetraMesh, degree: int) -> tuple[IntArray, FloatArray]

Enumerate continuous Pk nodes using exact topological identities.

tetra_basis(degree: int, bary: FloatArray) -> tuple[FloatArray, FloatArray]

Return Basix Pk values and canonical four-axis barycentric derivatives.

The extension and unit-sum point convention follow :func:pymhm.fem.scalar.tetrahedron_topology.tetra_values_gradients for every degree.

tetra_face_basis(degree: int, bary: FloatArray, *, opposite_vertex: int) -> FloatArray

Restrict the executed nodal Pk basis to one exact reference face.

bary has four unit-sum coordinates with coordinate opposite_vertex exactly zero. Entries whose integer nodal weight on that vertex is nonzero have identically zero trace and are excluded through the declared topology. Values on the face support retain the native volume table unchanged. The volume basis, its coefficients and physical derivative maps are unaffected.

tetra_tabulate(mesh: TetraMesh, degree: int, bary: FloatArray, *, backend: Literal['portable', 'basix'] = 'basix') -> tuple[IntArray, FloatArray, FloatArray, FloatArray]

Return nodal DOFs/coordinates, basis values and physical gradients.

Basix supplies equispaced Pk values/Cartesian derivatives in the literal :func:tetra_indices node order. Barycentric points must sum to one. PyMHM supplies the affine cell geometry and continuous DOF topology. backend='portable' is a compatibility spelling for this same Basix execution; there is no separate cardinal-polynomial implementation.

tetra_element_tabulate(mesh: TetraMesh, degree: int, bary: FloatArray, *, backend: Literal['portable', 'basix'] = 'basix') -> tuple[IntArray, FloatArray, FloatArray, FloatArray, FloatArray]

Tabulate continuous Pk values, physical gradients and physical Hessians.

backend follows :func:tetra_tabulate. Basix tabulates both derivative orders to the library before applying the affine Jacobian.

tetra_operators(mesh: TetraMesh, degree: int = 2, *, diffusion: Any = 1.0, source: Any = 0.0, order: int = 5, element_backend: Literal['portable', 'basix'] = 'basix') -> tuple[Any, Any, FloatArray]

Assemble scalar diffusion, consistent mass and source with positive quadrature.

Basix tabulates equispaced Pk in the declared nodal order. Integration and affine geometry retain their conventions. element_backend='portable' is a compatibility spelling for this same Basix execution. Compiled binary64 diffusion integration compensates both Cartesian contractions and quadrature sums; mass and source retain their formulas.

tetra_boundary_nodes(mesh: TetraMesh, degree: int) -> IntArray

Return sorted continuous-Pk nodal coordinates on the exterior boundary.

Exterior incidence is tested in affine barycentric coordinates against the opposite local vertex, with the existing geometric tolerance 1e-12. This identifies nodes; it neither projects boundary data nor chooses a PDE sign.

Canonical arbitrary-order Raviart--Thomas moments on affine tetrahedra.

RT_m = [P_m]^3 + x times homogeneous P_m. Face moments use the declared Pk normal tests; interior moments test [P_(m-1)]^3. The interior test basis is Bernstein, and its degrees of freedom are moments, not expansion coefficients.

RTTetraFamily dataclass

Bases: HDiv3DFamily

RT_m tetrahedral flux with complete P_m divergence and normal trace.

This family shares canonical physical face orientation and affine Piola mapping with :class:HDiv3DFamily. Its interior coordinates are the canonical volume moments, without bubble orthonormalization.

coefficients: FloatArray property

Return the executed RT candidate-to-moment-dual basis for archival replay.

__init__(degree: int = 0) -> None

Select a nonnegative RT polynomial order.

tabulate(points: FloatArray, *, coefficients: FloatArray | None = None) -> tuple[FloatArray, FloatArray, FloatArray]

Evaluate flux, divergence and complete pressure monomials.

rt3d_interior_tests(points: FloatArray, degree: int) -> FloatArray

Tabulate the componentwise P_(m-1) canonical volume moment tests.

H(div) families with independent face order and complete interior enrichment.

Tetrahedra retain all zero-normal BDM_(p+1) bubbles with normal Pk, k<=p+1. Prisms use horizontal [P_(p+1)(triangle) tensor P_p(z)]^2 and vertical P_p(triangle) tensor P_(p+1)(z), restricted to triangular Pk and rectangular Qk normals, k<=p. Divergence is the complete Pp or W_pp pressure space.

HDiv3DFamily dataclass

Mixed polynomial family with separate normal and pressure degrees.

Tetrahedral pressure degree p uses all zero-normal bubbles of [P_(p+1)]^3; normal degree k may range from zero to p+1. Prisms use the tensor family declared in this module, with zero<=k<=p. Default k=p=1 reproduces the 18-mode tetrahedron and 27-mode prism; p=2,k=1 gives 32 tetrahedral modes. Fixed normal order and increasing pressure degree enrich only interiors.

face_sizes: tuple[int, ...] property

Return complete triangular Pk or rectangular Qk normal moment counts.

local_size: int property

Return the number of local flux modes, including zero-normal interior bubbles.

pressure_size: int property

Return the complete divergence-space dimension.

interior_size: int property

Return the number of cell-private flux modes.

coefficients: FloatArray property

Return the read-only component-monomial coefficient matrix of this basis.

Rows are grouped by vector component. The original normal-P1 tetrahedral p=1,2 and prismatic p=1 families retain their lexicographic monomial candidates. Other orders use the stable Bernstein candidates declared in pymhm.fem.hdiv.moments_3d.candidates. Columns are face moments followed by uniquely oriented interior bubbles. The family degrees and candidate convention are part of the coefficient-vector contract. Save this matrix alongside persisted flux DOFs and pass it to :meth:tabulate when replaying those coordinates.

interior_moment_seeds: tuple[tuple[int, tuple[int, ...]], ...] property

Return the ordered (component, monomial exponent) interior seed tests.

Bubble functions are dual to these tests before positive-diagonal Cholesky orthonormalization. The final interior DOFs are their L2 moments against the resulting orthonormal bubble functions; face lifts are L2 orthogonal to every bubble.

__post_init__() -> None

Validate the supported mathematical spaces.

tabulate(points: FloatArray, *, coefficients: FloatArray | None = None) -> tuple[FloatArray, FloatArray, FloatArray]

Evaluate reference flux, analytic divergence and pressure monomials.

coefficients is an optional archived reference basis matrix with exactly the same row/column convention as :attr:coefficients. Supplying it preserves the executed basis when replaying stored DOFs.

reference_vertices(kind: CellKind) -> FloatArray

Return the reference vertices in the declared tetrahedral or prism order.

reference_faces(kind: CellKind) -> tuple[tuple[int, ...], ...]

Return triangle corners or tensor-ordered quadrilateral corners of each face.

cell_quadrature(kind: CellKind, order: int) -> tuple[FloatArray, FloatArray]

Return positive reference volume quadrature, with physical reference weights.

face_quadrature(corners: int, order: int) -> tuple[FloatArray, FloatArray]

Return positive unit-triangle or unit-square face quadrature.

face_shape(uv: FloatArray, corners: int) -> FloatArray

Evaluate affine triangle or bilinear tensor-square geometric shape functions.

face_polynomials(uv: FloatArray, corners: int, degree: int = 1) -> FloatArray

Evaluate complete Pk triangle or Qk square normal moment tests.

Degrees zero and one preserve the monomial convention 1, u, v, uv. Higher orders use orthogonal polynomials normalized in the mean-square face norm (Dubiner on triangles, tensor Legendre on quadrilaterals). Their first function is one. Ordering follows total degree and descending first index, with no singular evaluation at triangle vertices.

face_size(corners: int, degree: int) -> int

Return the dimension of a complete Pk or tensor Qk face space.

Moment-defined BDM and tensor-prismatic spaces of arbitrary polynomial order.

Complete BDM spaces are restricted through their normal polynomials, rather than by discarding high face moments. All zero-normal bubbles are retained. The construction uses explicit Nedelec interior tests and positive triangular factors; no singular-vector orientation defines persisted coordinates.

candidates(kind: CellKind, degree: int, points: FloatArray) -> tuple[FloatArray, FloatArray]

Tabulate stable Bernstein vector candidates and their divergence.

Rows in a persisted matrix are grouped by component. Within a component, simplex barycentric exponents are lexicographic; prism factors enumerate the triangle first and then the interval. Horizontal prism degrees are (degree+1, degree), and vertical degrees are (degree, degree+1).

tetrahedral_candidate_coefficients(degree: int) -> FloatArray cached

Preserve archived component/Bernstein rows through native BDM interpolation.

interior_tests(kind: CellKind, degree: int, points: FloatArray) -> FloatArray

Evaluate ordered interior moment tests for pressure degree degree.

Tetrahedral BDM_(degree+1) uses the first-kind Nedelec space of order degree-1. Prisms use triangular Nedelec tests tensor P_degree(z) in horizontal components and P_degree(x,y) tensor P_(degree-1)(z) vertically. Negative polynomial orders denote the empty space.

coefficients(kind: CellKind, pressure_degree: int, normal_degree: int) -> FloatArray cached

Construct uniquely oriented face lifts and all normalized interior bubbles.

Moment equations use orthogonal face tests and QR-conditioned interior tests. The triangular change of tests is also applied to the right-hand sides, preserving the declared Nedelec moment coordinates and their positive-diagonal bubble orientation. This avoids an ill-conditioned monomial moment inversion even when wider accumulation is unavailable.

Tensor-product Raviart--Thomas elements with independent interior enrichment.

On each rectangle, normal traces have degree k and interior vector functions have RT order s=k+n. The pressure is Q_s. Face DOFs are integral Legendre moments; cell DOFs are coefficients of zero-normal-trace polynomial bubbles. The construction implements the face/interior separation used by the mixed MHM family of Duran et al. (2019), without copying a reference implementation.

tensor_rt_moment_coefficients(degree: int, enriched_degree: int) -> FloatArray cached

Express declared face lifts and bubble coordinates in native rectangular RT.

tensor_rt_basis(mesh: CartesianMacroMesh, degree: int, enrichment: int, points: FloatArray) -> tuple[FloatArray, FloatArray, FloatArray]

Evaluate globally oriented RT vectors, divergence and modal Q_s pressure.

points are reference coordinates in [0,1]^2, either shared (q,2) or cell dependent (cells,q,2). Vector/divergence arrays include a leading cell axis. The pressure basis is cell independent and ordered with x degree varying fastest. The Piola transform preserves all normal moments even on anisotropic rectangles.

tensor_rt_dofs(mesh: CartesianMacroMesh, degree: int, enrichment: int) -> IntArray

Map conforming face moments and independent cell bubbles to global DOFs.

tensor_rt_trace_map(mesh: CartesianMacroMesh, cell: int, fine: CartesianMacroMesh, skeleton: SkeletonSpace, degree: int) -> FloatArray

Integrate each oriented macro flux against fine-edge Legendre moments.

tensor_rt_operators(mesh: CartesianMacroMesh, degree: int, enrichment: int, permeability: Any, source: Any, order: int) -> tuple[Any, Any, FloatArray]

Assemble mixed inverse-permeability mass and divergence moments.

One-sided polynomial traces and local scalar boundary integration.

PreparedScalarTrace dataclass

Read-only face geometry for repeated one-sided nodal trace evaluations.

Intervals retain the original fine-edge order and endpoint tolerance. A sorted interval index visits only overlapping pieces; the original last assignment wins at shared endpoints, including tolerance overlaps. Arrays are snapshots of the geometry, while each evaluation reads current field coefficients. No process-global geometry cache or field averaging is used.

selections(parameter: FloatArray) -> Iterator[tuple[FloatArray, IntArray, IntArray]]

Yield overlapping sample indices in the original fine-boundary order.

evaluate(coefficients: FloatArray, parameter: FloatArray) -> FloatArray

Evaluate without rebuilding edge geometry or changing endpoint ownership.

edge_basis(degree: int, parameter: FloatArray) -> FloatArray

Evaluate edge nodes ordered as endpoints followed by interior nodes.

edge_pieces(macro: TriangleMesh, fine: TriangleMesh, face: int, degree: int) -> Iterator[tuple[FloatArray, IntArray]]

Yield oriented macro parameters and nodal DOFs of matching fine edges.

prepare_scalar_trace(macro: TriangleMesh, fine: TriangleMesh, face: int, degree: int) -> PreparedScalarTrace

Prepare immutable oriented edge pieces and their exact tolerance intervals.

scalar_trace(macro: TriangleMesh, fine: TriangleMesh, face: int, degree: int, coefficients: FloatArray, parameter: FloatArray) -> FloatArray

Evaluate a continuous local trace without averaging across macrofaces.

diffusive_boundary_matrix(macro: TriangleMesh, fine: TriangleMesh, faces: tuple[int, ...], degree: int, velocity: Any, order: int) -> Any

Assemble +half beta.n boundary mass for prescribed diffusive flux.

On these exterior faces the retained multiplier is -K grad(u).n, rather than the interior Robin multiplier (-K grad(u)+beta*u/2).n.

strong_boundary_dofs(macro: TriangleMesh, fine: TriangleMesh, faces: tuple[int, ...], degree: int, value: Any) -> tuple[IntArray, FloatArray]

Collect unique local nodal Dirichlet DOFs and evaluate the boundary data.

Scalar skeletal moments on the union of incident physical fine-edge partitions.

incident_face_breaks(skeleton: SkeletonSpace, local_meshes: Sequence[TriangleMesh] | Mapping[int, TriangleMesh], face: int, degree: int) -> FloatArray

Return integration cuts without enriching the skeletal approximation space.

The partition contains the skeletal breaks and all fine-edge endpoints on either incident side. Coincident physical endpoints are merged at the same floating-point threshold used by the edge trace assembly.

integrate_dirichlet_trace(skeleton: SkeletonSpace, local_meshes: Sequence[TriangleMesh], values: Any, degree: int, order: int, *, faces: Sequence[int] | None = None) -> FloatArray

Integrate scalar pressure data against Lambda on every exterior macroface.

Fine-edge cuts resolve prescribed data that are piecewise polynomial on the local boundary. They do not change Lambda's degrees or breakpoints. Additional discontinuities of an arbitrary callback still require an integration partition that resolves those discontinuities. faces selects exterior macrofaces; omitted faces have zero moments and the callback is never evaluated there.

Basix tensor Qk tables and scalar operators on Cartesian rectangles.

Rows retain the declared x-fastest nodal ordering. Product quadratures have unit-sum weights; physical cell measures are applied by the operators.

quadrilateral_quadrature(order: int = 4) -> tuple[FloatArray, FloatArray]

Return product Gauss points on [0,1]^2 and unit-sum integration weights.

cardinal_polynomials(degree: int) -> tuple[FloatArray, ...] cached

Export the reordered Basix interval basis in ascending power coordinates.

This representation conversion supports coefficient archives; evaluation and derivatives use Basix directly. The orthonormal native polynomial j equals sqrt(2j+1) times the shifted conventional Legendre polynomial.

qk_basis(degree: int, points: Any) -> tuple[FloatArray, FloatArray]

Evaluate equidistant Qk basis and reference gradients at unit-square points.

Nodes and basis functions are ordered first in x, then in y. Basix interval factors supply values and derivatives, including at nodal points; their tensor product retains the historical coefficient order.

qk_space(mesh: CartesianMacroMesh, degree: int) -> tuple[IntArray, FloatArray]

Return the shared nodal DOFs and coordinates of a continuous Qk space.

grid_resolves_material(mesh: CartesianMacroMesh, field: CartesianCellField) -> bool

Determine whether every material line coincides with a fine-grid line.

rectangle_intersection_quadrature(origins: FloatArray, spacing: FloatArray, field: CartesianCellField, order: int) -> tuple[FloatArray, FloatArray]

Integrate each rectangle separately on its exact intersections with pixels.

Reference points have shape (cells, points, 2); weights have shape (cells, points) and sum to one on each rectangle. Zero-weight padding combines variable intersection counts without extrapolating the material.

quadrilateral_operators(mesh: CartesianMacroMesh, degree: int = 1, *, permeability: Any = 1.0, source: Any = 0.0, order: int = 4) -> tuple[Any, Any, FloatArray]

Assemble sparse Qk diffusion, mass and load using bounded cell batches.

CartesianCellField coefficients are integrated over exact intersections with material pixels, retaining a uniform-rule fast path on aligned grids. Arbitrary callback discontinuities are not inferred: those require an aligned grid or a declared quadrature-sensitivity study.

quadrilateral_trace_coupling(coarse: CartesianMacroMesh, cell: int, fine: CartesianMacroMesh, skeleton: SkeletonSpace, degree: int) -> FloatArray

Integrate signed Qk traces, splitting at fine and skeleton breakpoints.

Scalar polynomial spaces and oriented skeleton coordinates on planar macroedges.

FaceSpace dataclass

Discontinuous Legendre or continuous nodal polynomials on one macroface.

Breaks parameterize the oriented face on [0, 1]. Each segment has its own degree. The default is one constant multiplier on the entire face. continuous=True shares endpoints between segments on this macroface, with no continuity imposed across different macrofaces. Degrees must then be positive; interior interpolation nodes are equally spaced.

size: int property

Return the number of scalar face unknowns.

__post_init__() -> None

Reject gaps, reversed intervals and inconsistent polynomial degrees.

uniform(degree: int = 0, subdivisions: int = 1, *, continuous: bool = False) -> FaceSpace classmethod

Construct an equally partitioned face with the same degree per segment.

evaluate(parameter: Any) -> FloatArray

Evaluate Basix scalar bases at oriented coordinates in [0, 1].

Continuous segment nodes are ordered left endpoint, right endpoint, then increasing interior nodes. Discontinuous modes use conventional unnormalized Legendre polynomials, whose degree-j squared mass on a segment of length L is L/(2j+1). Internal breaks belong to the segment on their right; a face endpoint belongs to its adjacent segment.

constant_coefficients() -> FloatArray

Represent the constant one exactly in this face's declared basis.

quadrature(order: int = 4) -> tuple[FloatArray, FloatArray]

Return Gauss points and weights, separately on every segment.

Normalize the reference rule to its exact mass two before mapping each subinterval. Extended accumulation prevents a rounded weight-sum defect from becoming an artificial volume change in nearly incompressible displacement boundary data. No field value is clipped. Returned coordinates and weights retain binary64 storage.

SkeletonSpace dataclass

Independent face partitions and degrees with component-interleaved DOFs.

size: int property

Return the total number of trace coefficients.

__init__(mesh: TriangleMesh, faces: tuple[FaceSpace, ...] | None = None, components: int = 1) -> None

Assign unique face numbers and validate the geometry-to-space association.

dofs(face: int) -> IntArray

Return the global DOFs belonging to a face.

cell_dofs(cell: int) -> IntArray

Concatenate trace DOFs in the cell's counterclockwise side order.

interface_pairing(test_space: Any, trial_space: Any, cell: int, *, order: int = 6) -> FloatArray

Compatibility entry point for the dimension-independent physical trace pairing.

Oriented triangular-face spaces, coupling and physical scalar boundary moments.

Trace coefficients represent densities in the macroface's canonical normal. Bernstein coefficients are independent or shared within each physical macroface.

TriangularSkeleton

Piecewise Pk flux densities, discontinuous or C0 within each macroface.

subdivisions counts uniform edge divisions, giving its square many subtriangles by default. Explicit face_partitions instead prescribe a conforming nonuniform triangulation of each original macroface. Each subtriangle has (k+1)(k+2)/2 local Bernstein modes. continuous=True shares vertex and edge coefficients inside a macroface and requires k >= 1. Sequences select independent face degrees and continuity flags. No degrees of freedom are shared between macrofaces, including their common edges. Flux uses the canonical normal of TetraMesh.

__init__(mesh: TetraMesh, subdivisions: Any = 1, *, degree: Any = 0, face_partitions: tuple[FloatArray, ...] | None = None, continuous: Any = False) -> None

Create uniform or explicit conforming triangular partitions with Pk modes.

Explicit vertices use barycentric coordinates in mesh.faces[face] ordering. Every original face must be covered without holes, overlaps or hanging internal edges. subdivisions and partitions retain the uniform-grid configuration; face_partition returns the active geometry.

face_partition(face: int) -> FloatArray

Return active subtriangles in original macroface barycentric coordinates.

face_weights(face: int) -> FloatArray

Return positive subtriangle area fractions summing to one on a macroface.

basis(face: int, barycentric: FloatArray) -> FloatArray

Evaluate the Bernstein Pk density basis on a canonical subtriangle.

Degrees zero and one retain their constant and barycentric coordinate conventions. All modes integrate to the same fraction of subface area, and coefficients equal to one represent a constant density exactly.

subtriangle_dofs(face: int, segment: int) -> IntArray

Return global coefficients in the subtriangle's local Bernstein ordering.

C0 vertex and edge coefficients are shared only within this macroface. The map transports local orientations to the declared canonical numbering.

evaluate(face: int, barycentric: FloatArray) -> FloatArray

Evaluate the complete macroface basis in original barycentric coordinates.

Columns follow dofs(face). At partition interfaces the first incident subtriangle supplies the DG one-sided value; C0 values agree on both sides.

dofs(face: int) -> IntArray

Return the global coefficients for one validated macroface index.

cell_dofs(cell: int) -> IntArray

Concatenate the four oriented macroface coefficient blocks of a cell.

tetra_trace_coupling(coarse: TetraMesh, cell: int, fine: TetraMesh, skeleton: TriangularSkeleton, degree: int) -> FloatArray

Integrate signed nodal Pk traces against aligned triangular Bernstein modes.

tetra_boundary_data(skeleton: TriangularSkeleton, dirichlet: Any, neumann: dict[int, Any], order: int) -> tuple[FloatArray, dict[int, float]]

Integrate weak pressure moments and project prescribed outward normal fluxes.

The quadrature floor integrates the Bernstein Gram matrix exactly even when the independently selected trace degree exceeds the local degree. Nonpolynomial boundary data can require a higher caller-supplied order.

Basix RT spaces with declared tensor moments and contravariant Piola maps.

Normal DOFs are integral tensor-Legendre moments. Scalar Qk pressure uses an ordinary pullback; divergence carries the inverse geometric determinant.

tensor_legendre_values(degree: int, points: FloatArray) -> FloatArray

Evaluate tensor Legendre polynomials, with the last coordinate varying fastest.

mapped_rt_dofs(mesh: HexMesh, degree: int) -> IntArray

Map shared integral normal moments and independent cell-interior bubbles.

mapped_reference_coefficients(degree: int) -> FloatArray cached

Express declared hexahedral face lifts and bubble coefficients in native RT.

mapped_rt_basis(mesh: HexMesh, degree: int, points: FloatArray) -> tuple[FloatArray, FloatArray, FloatArray]

Tabulate oriented physical RT_k vectors, their exact Piola divergence and scalar Q_k.

Sparse scatter of user-integrated element forms in declared coefficient maps.

assemble_element_blocks(blocks: Any, rows: Any, columns: Any, shape: tuple[int, int]) -> Any

Scatter element matrices to a CSC operator without selecting a PDE.

blocks[e,i,j] pairs test coefficient rows[e,i] with trial coefficient columns[e,j]. Repeated global coordinates sum, including repeated coordinates within an element. Square matrices declare identical trial/test maps and equal dimensions explicitly. Real/complex entries keep their numeric dtype. The caller owns integration, orientation signs and local basis transformations; none are inferred during this scatter.

DG mass, curl and oriented scalar/tangential trace operators.

TangentialTraceSpace

Scalar TM or two-component tangential density on actual macrofaces.

In 3D canonical faces use deterministic orthonormal tangent frames. Coefficients represent a density in the declared canonical tangent frame; the opposite incident cell receives the opposite sign. In scalar 2D pairings there is one coordinate per face basis. No normal coordinate is stored by the vector 3D pairing.

__init__(base: SkeletonSpace | TriangularSkeleton) -> None

Wrap scalar face spaces retaining their independent degree and partition.

dofs(face: int) -> IntArray

Return component-interleaved global tangent coefficients of one face.

cell_dofs(cell: int) -> IntArray

Concatenate face blocks in actual macrocell incidence order.

CurlOperators dataclass

DG mass blocks and magnetic evolution C=-curl, paired with electric -C transpose.

electric_mass: Any cached property

Return the block-diagonal SPD electric mass matrix.

magnetic_mass: Any cached property

Return the block-diagonal SPD magnetic mass matrix.

frequency_bound() -> float

Bound the mass-scaled curl using positive weighted Gram row sums.

For S=L_mu^-1 C L_epsilon^-T, any strictly positive x gives rho(S.T S) <= max_i (abs(S.T S) x)_i / x_i. A few positive power iterations tighten this Gershgorin bound; no Ritz estimate is accepted as an upper bound. The elementary 1/infinity-norm bound is retained if it is smaller. Roundoff padding covers sparse dot products.

scalar_basis(degree: int, bary: FloatArray) -> tuple[FloatArray, FloatArray]

Use cardinal simplex functions without merging DG element nodes.

quadrature(mesh: Any, material: Any, order: int) -> tuple[FloatArray, FloatArray, Any]

Return physical-weight simplex/rectangle rules split at material interfaces.

physical_points(mesh: Any, reference: FloatArray) -> FloatArray

Map common or elementwise reference coordinates to physical DG points.

physical_basis(mesh: Any, degree: int, bary: FloatArray) -> tuple[FloatArray, FloatArray]

Tabulate common or per-cell simplex/rectangle rules with physical gradients.

cell_basis(mesh: Any, cell: int, degree: int, points: FloatArray) -> FloatArray

Evaluate one incident DG polynomial at physical facet quadrature points.

facet_rule(mesh: Any, face: int, order: int) -> tuple[FloatArray, FloatArray]

Integrate a fine edge or triangular face with positive physical weights.

mass_blocks(mesh: Any, degree: int, material: Any, components: int, order: int) -> FloatArray

Integrate SPD mass blocks with component-interleaved element-local nodes.

derivatives(mesh: Any, degree: int, order: int) -> tuple[Any, ...]

Assemble central DG derivatives with incident traces at macro boundaries.

trace_coupling(skeleton: TangentialTraceSpace, cell: int, fine: Any, degree: int, order: int) -> Any

Integrate signed electric/trace products with a degree-dependent quadrature floor.

assemble_local(mesh: Any, skeleton: TangentialTraceSpace, cell: int, degree: int, permittivity: Any, permeability: Any, order: int) -> CurlOperators

Assemble paired scalar/vector 2D or vector/vector 3D DG curl operators.

Electric and magnetic labels identify the two mass-weighted coordinate groups; their differential pairing is C=-curl and -C.T. Coefficients are interleaved by component within each independent element-local node. This is central-flux DG, not a conforming Nedelec space. The tangential density represents H cross n and loads the electric equation. The paired operators implement Lanteri et al. (2018), DOI:10.1137/16M110037X, Sec. 5; its energy identity is distinct from discrete Gauss-law or error estimates, which are not established by this assembly.

field_values(value: Any, points: FloatArray, components: int) -> FloatArray

Evaluate a real scalar/vector Maxwell field with explicit component count.

volume_load(local: CurlOperators, value: Any, components: int, material: Any, order: int) -> FloatArray

Integrate a physical field against DG tests with its optional mass material.

tangential_rules(skeleton: TangentialTraceSpace, faces: Iterable[int], order: int) -> Iterator[tuple[int, FloatArray, FloatArray, FloatArray, FloatArray, IntArray]]

Yield physical face rules, canonical frames and independent trace coordinates.

The basis and weights retain every explicit segment or triangular partition. Scalar edge pairings use the unit scalar frame; vector face pairings use the declared canonical orthonormal tangent frame. The quadrature floor integrates the polynomial trace mass exactly; nonpolynomial material and boundary data require independent order checks.

tangential_mass(skeleton: TangentialTraceSpace, face_weights: Mapping[int, Any], order: int) -> Any

Assemble a positive scalar-weighted mass in canonical tangent coordinates.

tangential_load(skeleton: TangentialTraceSpace, value: Any, faces: Iterable[int], order: int) -> FloatArray

Integrate a scalar/vector field's canonical tangential moments on chosen faces.

A callable receives (physical_points, canonical_normals). Only the tangential projection contributes in 3D. Each face is added once in the caller's explicit order; this routine does not infer boundary conditions.

tangential_moments(operators: Iterable[CurlOperators], fields: Iterable[FloatArray], trace_size: int) -> FloatArray

Accumulate signed incident tangential field moments in explicit cell order.