formoniq: FEEC at any polynomial degree

[luiswirth]

Stacked on #130, which it needs for the tensor layer. Review that one
first; the base will retarget to main once it merges.

One commit, passing fmt, clippy, test and doc.

What it is

P⁻_r Λ^k is Sym^r ⊗ Λ^k in barycentric coordinates, so the trimmed
family needs no construction of its own: it is the tensor the algebra
already builds, read at a degree.

Barycentric coordinates are what make the coefficient factor homogeneous.
Σ λ_i = 1 supplies the missing variable, so
dim P_r(ℝ^n) = binom(n+r, n) = dim Sym^r(ℝ^(n+1)), no graded layer sits
above a homogeneous factor, and a polynomial degree is one number exactly
as a grade is.

The two structural operators are one operation in opposite directions.
d moves a degree from the coefficients to the blade, the Koszul
operator moves it back, and dκ + κd = (r+k) id on the homogeneous
pieces is what cuts out the trimmed spaces. Whitney forms fall out as
k! κ(dλ_σ), so the familiar deletion formula
W_σ = k! Σ_i (-1)^i λ_{σ_i} dλ_{σ∖σ_i} is a consequence of what κ does
to a blade rather than a definition.

The geometric decomposition

P⁻_r Λ^k(K) = ⊕_{f ⊆ K} int_f, a direct sum over every subsimplex, each
summand the piece whose support is all of f and no smaller face. That
is what makes the global space a subspace of HΛ^k rather than a direct
sum of unrelated cell spaces.

Two cells sharing a face give that face's degrees of freedom the same
global index with no communication and no orientation convention: a dof
is named by data intrinsic to f (a subsimplex of f and a monomial in
f's own barycentric coordinates, both in positions within f), and a
skeleton stores every simplex colex-sorted, so both cells see the same
order. Nothing here consults a CellOrdering.

Assembly

Through an explicit local-to-global map, so one routine serves every
degree and the first-order case, where a dof is a k-simplex, is that map
rather than a different mechanism.

The mass matrix is integrated, not quadratured. The integrand is
polynomial in the barycentric coordinates and their moments are the
closed form n! α! / (|α| + n)!, so there is no rule to be exact to and
no order to check. HodgeMassElmat is this at r = 1.

TrimmedComplex implements the same HilbertComplex, so every problem
written against that interface — source, eigenvalue, heat, wave, Dirac —
runs at any degree with no solver touched. tests/higher_order.rs
exercises all five at r = 2, each against a law rather than a stored
number.

Convergence

Hodge-Laplace eigenvalues on [0, π]^n under the absolute boundary
condition, against the exact bottom of the spectrum:

dim 3, grade 1	ndofs	error	rate
r=1	4913	1.2e-2	1.91
r=2	3665	4.4e-4	3.50
r=3	10345	1.8e-6	5.75

Rates come out at 2, 4 and 6 across dimensions and grades, matching the
predicted O(h^(2r)).

Compared against the exact eigenvalue rather than a Richardson
extrapolation on purpose: at r = 3 the error reaches roundoff within
two refinements, and a self-convergence estimate cannot tell that from a
bug. Getting the exact value right needed deriving it: the Hodge
Laplacian on the box splits into scalar Laplacians with a Dirichlet
factor on each direction in I and a Neumann factor on each outside it,
so the smallest nonzero eigenvalue is max(k, 1) — not the specific
eigenform's λ = n, which the first run wrongly compared against.

Assembly costs 2.1 ms, 13 ms and 77 ms for r = 1, 2, 3 on a fixed 3D
mesh; space setup 193 µs, 559 µs and 1.22 ms.

Ties to trusted code

The r = 1 case is identical, not merely close, and three tests say so:
the trimmed basis is pointwise the Whitney basis, the reference d is
exactly the simplicial coboundary with entries in {−1, 0, 1}, and the
assembled dof counts, mass matrices and differential equal the Whitney
complex's on real meshes.

Ported to the new tensor layer

PolyForm is Sym^r ⊗ Λᵏ over the barycentric ℝⁿ⁺¹, which is now
stated in the slots rather than in a tensor-wide dimension: both slots
carry cell_dim + 1 explicitly. That the coefficient and blade factors
live over the same enlarged space, and that at_bary leaves it by
pulling back along unit_difbarys, is visible in the types.

What the tests caught

• The discrete d was scattered wrongly. It filled each column from
  a single cell, but a basis function's derivative is supported on all
  cells containing its dof, so the scatter now averages over the
  multiplicity. Caught by the r=1-equals-Whitney test.
• The local-to-global map capped meshes at 64 vertices, by routing a
  Simplex through a 64-bit Combination of global vertex indices.
  Caught by the p-convergence run, not by a test — no small-mesh test can
  reach it.
• Q_r ⊗ H was Kroneckered by hand under a comment asserting the
  order matched the tensor's strides, with nothing enforcing it.
  Reversing it yields a symmetric positive definite matrix of the right
  shape and the wrong meaning, so no shape check catches it. Now one
  function and one test.

Not done

A higher-order L² error against an analytic solution needs an interpolant
into P⁻_r Λ^k. The convergence study uses eigenvalues, which need none,
but a source-problem convergence table would want it.

[luiswirth]

Parking this. #130 is in, and the next work is architectural: gramian moves above multialgebra so a metric is the Sym²V* tensor it is, and simplicial splits into topology and a new regge crate carrying the metric geometry — invariants 5 and 1 each becoming a crate boundary.

Both touch the layers this builds on, so rebasing this first would mean doing it twice. Revisiting once the ladder settles. The work is complete and green as it stands; nothing here is in doubt.