Correct and consolidate Form.h's entity-domain mapping

[jorgensd]

Correct and consolidate Form.h's entity-domain mapping

Form's constructor maps integration entities in this->mesh() to the
corresponding cells of an argument's or coefficient's mesh. That mapping has
five defects, all currently masked because only facet integrals reach the code.
Fixing them removes an arbitrary restriction, and mixed-dimensional ridge
integrals fall out as a consequence.

cpp/dolfinx/fem/Form.h                           +89 -60
cpp/dolfinx/fem/utils.cpp                         +1 -22
cpp/dolfinx/fem/utils.h                           +1 -21
cpp/test/fem/form.cpp                            +52  -1
python/doc/source/release_notes.md               +12
python/test/unit/fem/test_assemble_submesh.py   +113
.github/workflows/fenicsx-refs.env                +1  -1   ← temporary, see below

1. Entity dimension conflated with the argument mesh's codimension

In both the argument and the coefficient block, before:

int codim = tdim - mesh0->topology()->dim();          // from the ARGUMENT mesh
auto c_to_f = topology.connectivity(tdim, tdim - 1);  // hard-coded facets

codim describes the argument mesh; the connectivity has to be indexed by the
integral's entity dimension. They coincide only for facet integrals, so the
code reads as if it handles the general case and does not. Any other integral
type would index the facet list with a non-facet local index.

The entity dimension is now computed once and used for the lookup:

const int edim = integral_entity_dim(type, tdim);
...
std::shared_ptr<const graph::AdjacencyList<std::int32_t>> c_to_e
    = topology.connectivity(tdim, edim);

The same conflation remains in pack.h on the Expression path — out of
scope here, noted below.

2. Three copies of the IntegralType → entity dimension switch

It was spelled out inside impl::entity_permutations in Form.h, and again in
fem/utils.h and fem/utils.cpp. They have to agree and nothing enforced it.
Replaced by one

constexpr int integral_entity_dim(IntegralType type, int tdim);

beside the IntegralType enum. It lives in Form.h rather than utils.h
because of an include cycle — fem/utils.h → Function.h → assembler.h →
the assemble_*_impl.h headers — which makes a helper defined in utils.h
invisible to the assemblers that include it. Form.h is upstream of that cycle.
This accounts for the deletions in utils.{h,cpp}. It throws on an unrecognised
type rather than returning a plausible-looking dimension.

3. An arbitrary codimension restriction

compute_facet_domains threw "Codimension > 1 not supported." from a branch
that is already general: it works entirely through the passed-in adjacency, and
ridge uses the same 2-wide (cell, local_entity) row layout as
exterior_facet. The codim == 1 / else throw pair collapses to one
codim >= 1 branch, and the function is renamed compute_entity_domains to
match what it does. Lifting the restriction is deletion, not addition.

4. A comment asserting the opposite of what happens

The old comment claimed the local-entity column is reused "since
create_submesh preserves the local facet index". It does not — on the submesh
that entity is the cell, so the column has no meaning. It is harmless only
because no consumer reads column 1: the vector and matrix assemblers read
entities0(f, 0), and coefficient packing takes
submdspan(entities, full_extent, 0). The comment now says that.

5. A precondition the mapping never checked

compute_entity_domains maps an integration entity to a cell of the other
mesh, so it can only express the case where the integral's entity dimension
equals that mesh's dimension. Nothing checked it. Form now throws naming both
dimensions.

This is not reachable from Python: FFCx rejects both mismatched combinations
while compiling the kernel, because the quadrature points have the wrong
dimension for the element.

ridge integral  + 2D submesh data → "Point dim (1) does not match element dim (2)."
vertex integral + 1D submesh data → "Point dim (0) does not match element dim (1)."

So the check is a backstop for callers assembling a Form directly in C++ with
their own kernel, and it documents the precondition — it is not fixing a bug
users can hit. It is covered by a new Catch2 test that builds such a Form and
matches the message, so it is not carried unexercised.

The argument and coefficient blocks were otherwise identical and are now one
shared map_entities helper, reducing each call site to a single line.

What this enables

With the mapping correct, a form may integrate over the ridges of a parent mesh
with an argument or coefficient on a codimension-2 submesh. The other two pieces
already exist:

• FFCx generates the kernels — FEniCS/ffcx#879
  adds the closure-dofs gather that pulls a submesh's coordinate dofs out of the
  parent cell's coordinate_dofs buffer, indexed by quadrature_permutation.
• The preceding commit on this branch makes the ridge quadrature_permutation
  real. Ridge kernels were previously handed an unconditionally-zero
  permutation, so that gather would have read the wrong closure ordering on any
  reflected edge.

Nothing else needed changing. EntityMap is a pure index map with no notion of
codimension; create_submesh, compute_integration_domains (whose
exterior_facet/vertex/ridge paths already share one dim-parameterised
branch), coefficient packing and all three assemblers are already
codimension-generic.

IntegralType::vertex still throws, with a message naming the supported types.

Tests

Python — python/test/unit/fem/test_assemble_submesh.py, 8 new cases.
Nothing previously put a codim-2 submesh in a form; test_submesh_full creates
one but checks only topology and geometry.

test_mixed_dom_codim_2 integrates over every edge of a unit cube with data on
the edge submesh, and compares against integrating over the submesh directly:

coefficient          submesh dx = 83.685432181533   parent dr = 83.685432181533
grad(g)·grad(g)      submesh dx = 222.087665379149  parent dr = 222.087665379149
SpatialCoordinate    submesh dx = 343.500436540654  parent dr = 343.500436540654
g * x[0]             submesh dx = 42.105345728328   parent dr = 42.105345728328

The grad and SpatialCoordinate cases go through the FFCx gather and so
through the ridge permutations. The test is not vacuous: 58% of the cube's
edges are reflected relative to their parent cell's local orientation, and the
forms report needs_facet_permutations (asserted). Degree 2 is the
orientation-sensitive case — an edge's dofs are [v0, v1, midpoint], so a wrong
permutation swaps the endpoints.

test_mixed_dom_codim_2_arguments covers arguments rather than coefficients,
for both a matrix and a vector. Both run over GhostMode.none and
GhostMode.shared_facet, and in parallel — the "first incident cell" choices in
fem/utils.h and mesh/utils.h are made independently and are rank-local,
reconciled only by the permutation.

Two negative tests pin that the mismatched combinations are rejected, without
asserting which layer does the rejecting (see §5).

Six of these eight fail without this change, with
"Integral type not supported.".

C++ — cpp/test/fem/form.cpp, 1 new case. Builds a Form directly with a
dummy kernel, a ridge integral and a facet submesh, and uses CHECK_THROWS_WITH
so it pins §5's message rather than any exception.

A release-notes entry is added under v0.12.0 (draft).

Verification

• Python: 2861 passed, 69 skipped, 5 xfailed across fem/ and mesh/.
• C++: 90 test cases, 87 passed, 3 skipped, 29807 assertions.
• codim-2 tests also pass under mpirun -n 3.
• clang-format, ruff check and ruff format --check clean.

Before merging

.github/workflows/fenicsx-refs.env currently pins ffcx_ref to
dokken/manifold-fixes so CI can run the grad and SpatialCoordinate cases.
This is a temporary commit and must be reverted to main once
FEniCS/ffcx#879 has merged.

Required first, in order:

1. FEniCS/ffcx#882 — Fix vertex mesh
   compilation with geometry tables (dokken/point-coordinate-tables). Removes
   the assert ttype != "ones" that blocked point-mesh coordinate elements.
   Not functionally required by this PR's tests — a tet parent with an
   interval submesh has a P1 coordinate element, so the all-ones case never
   arises. It is in the chain only because
   FEniCS/ffcx#879 carries it.
2. FEniCS/ffcx#879 — Fix coordinate
   element access of codim >=1 assembly (dokken/manifold-fixes). The real
   functional dependency; without it the grad(g)·grad(g) and
   SpatialCoordinate cases cannot compile.
3. FEniCS/dolfinx#3904 — the
   ridge-permutations PR. The base this is stacked on, not a sibling.

Out of scope

• Mixed-dimensional vertex integrals — same pattern, and FFCx would need the
  0-dimensional quadrature case first.
• The third copy of the codim logic on the Expression path in pack.h,
  carrying the same conflation plus a parent-to-submesh-only restriction. It is
  already codim-generic and correct for codim 2.
• Two-sided ridge integrals — get_cell_entity_pairs is
  static_assert(num_cells == 1); patch assembly is a separate feature.

AI disclosure

This PR was generated with the help of CLAUDE (OPUS 5). I have reviewed,
modified and tested the code.