test: enable the value check on the symmetric weak Hessian

docs/upstream-bugs.md

@@ -1,6 +1,7 @@
 # Bugs found in the upstream Python `DiffusionGeometry`
 
-Defects discovered in the reference implementation while porting it to Julia.
+Defects discovered in the reference implementation while porting it to Julia — and,
+in §5, in its *test suite*, where the shipped code is right but nothing checks it.
 
 - **Upstream:** https://github.com/Iolo-Jones/DiffusionGeometry
 - **Revision audited:** `f45b39f` ("update mdg notebook"), clean tree.
@@ -10,14 +11,17 @@ Where a bug has a numerical parity target, the Julia port implements the *correc
 behaviour and `pyparity/gen_fixtures.py` stores a **corrected** reference, so the
 parity gate compares against the right answer rather than pinning the bug. Each
 correction is guarded: if upstream fixes the bug, fixture generation fails loudly
-rather than silently double-correcting. See `PORTING_PLAN.md` §8.
+rather than silently double-correcting. See `PORTING_PLAN.md` §8. A test defect has
+no parity target and so no such guard — §5 is a coverage hole, closed on the Julia
+side only.
 
 | # | Location | Severity | Status in the port |
 |---|---|---|---|
 | 1 | `utils/basis_utils.py::_perm_tables` | **Wrong answers, silent** | Fixed; fixture corrected |
 | 2 | `tensors/vector_fields/vector_field.py::VectorField.from_reconstruction` | **Dead code** (`AttributeError`) | Reimplemented; round-trip gated |
 | 3 | `methods/geodesics.py::geodesic_distances_function` | **Dead code** (`AttributeError`) ×2 | Fixed; fixture corrected |
 | 4 | `utils/basis_utils.py::form_to_ambient_polyvector` (`k == 0` branch) | Cosmetic (unreachable) | N/A |
+| 5 | `tests/test_src/test_hessian.py::test_hessian_02_sym_weak_matrix` | Test defect (coverage hole) | Corrected expectation; check enabled |
 
 ---
 
@@ -207,3 +211,43 @@ a `FormSpace(dg, 0)` directly.
 **In the port.** Not reproduced. `to_ambient(::ScalarFunction)` is defined
 separately as `to_pointwise_basis(f)`, matching the reachable Python behaviour, and
 Julia's `FormSpace` constructor asserts `1 ≤ degree ≤ dim`.
+
+
+---
+
+## 5. The symmetric Hessian's value check is disabled because the *expectation* is wrong
+
+**Where:** `tests/test_src/test_hessian.py`, `test_hessian_02_sym_weak_matrix`.
+
+**What's wrong.** The test assembles the weak symmetric Hessian by hand and then does
+not compare against it:
+
+```python
+# Manual computation: H_sym_weak[i, s, I] = ∑_p μ[p] * φ_i[p] * H[p, j1(s), j2(s), I]
+...
+# NOTE: Disabling strict value check.
+assert hessian_sym_02_weak_computed.shape == (n1 * d_sym, n0)
+# assert np.allclose(hess_comp_sub, hess_sel_flat)
+```
+
+The expectation is missing the off-diagonal multiplicity factor. `hessian_02_sym_weak`
+doubles the `j₁ ≠ j₂` entries, and it is right to: the symmetric basis element for
+`j₁ ≠ j₂` is `dx_{j₁} ⊗ dx_{j₂} + dx_{j₂} ⊗ dx_{j₁}`, so pairing the (symmetric)
+Hessian against it picks up both entries. The same convention is baked into
+`gamma_02_sym`, whose off-diagonal rows and columns carry the matching ×2. Comparing
+the builder against an expectation that drops the factor fails on exactly the
+off-diagonal rows — so the assertion was commented out rather than the expectation
+fixed, and `hessian_02_sym_weak` was left with no value check on either side of the
+port.
+
+**Severity: test defect.** The shipped builder is correct; the coverage was not. But
+it was the *only* weak-operator builder whose values nothing verified — the fixture
+parity gate compares Julia against Python's output, which proves nothing if Python is
+wrong.
+
+**In the port.** The check is enabled, with the factor restored, in
+`test/test_pysrc.jl` → `hessian` → `hessian_02_sym_weak matrix`, and it passes at
+`d = 1..4`. Because that test asserts the very convention in question, a second test
+pins it convention-free: expanding the strong symmetric coefficients into the full
+`(0,2)` basis and pairing them with the full Gram reproduces `hessian_02_weak`, whose
+`⟨H(φ_I), φ_i dx_j ⊗ dx_k⟩` entries involve no symmetric-basis choice at all.

test/test_pysrc.jl

@@ -8,6 +8,10 @@
 #
 # `d` is the package's exterior-derivative accessor, so the ambient dimension is
 # spelled `dim` throughout to avoid shadowing it.
+#
+# One test here has no upstream counterpart: `hessian` → `hessian_02_sym_weak agrees
+# with the full (0,2) Hessian`, which corroborates the value check upstream disabled.
+# It is marked as port-added where it appears.
 
 using Combinatorics: combinations
 # `coeffs` and `np_reshape` are internal; the tests need them to poke at the
@@ -164,13 +168,52 @@ end
 # ── test_hessian.py ───────────────────────────────────────────────────────────
 @testset "hessian" begin
     for_each_geom() do dg, (dim, n, n0, n1)
-        d_sym = size(get_symmetric_basis_indices(dim), 1)
+        sym_idx = get_symmetric_basis_indices(dim)
+        d_sym = size(sym_idx, 1)
+        u, mu = function_basis(dg), measure(dg)
+        H = hessian_functions(dg.cache)
 
         # Upstream asserts only the shape here — its value check is commented out
-        # (`NOTE: Disabling strict value check`). Ported as-is; see the note in the
-        # summary rather than silently strengthening it.
-        @testset "hessian_02_sym_weak shape" begin
-            @test size(weak(hessian(dg))) == (n1 * d_sym, n0)
+        # (`NOTE: Disabling strict value check`). The check is enabled below; the
+        # expectation upstream compares against is the thing that is wrong, not the
+        # builder. See docs/upstream-bugs.md §5.
+        @testset "hessian_02_sym_weak matrix" begin
+            H_weak = weak(hessian(dg))
+            @test size(H_weak) == (n1 * d_sym, n0)
+
+            # H_sym_weak[i, s, I] = w_s · Σ_p μ_p φ_i(p) H[p, j₁(s), j₂(s), I], with
+            # w_s = 2 off the diagonal: the symmetric basis element for j₁ ≠ j₂ is
+            # dx_{j₁} ⊗ dx_{j₂} + dx_{j₂} ⊗ dx_{j₁}, so it picks up both entries of
+            # the (symmetric) Hessian. Upstream's expectation drops that factor.
+            i_sel, s_sel, I_sel = sample_indices(n1), sample_indices(d_sym), sample_indices(n0)
+            manual = zeros(length(i_sel) * length(s_sel), length(I_sel))
+            for (a, i) in enumerate(i_sel), (b, s) in enumerate(s_sel)
+                j1, j2 = sym_idx[s, 1], sym_idx[s, 2]
+                w = j1 == j2 ? 1.0 : 2.0
+                for (c, I) in enumerate(I_sel)
+                    manual[(a - 1) * length(s_sel) + b, c] =
+                        w * sum(u[p, i] * H[p, j1, j2, I] * mu[p] for p in 1:n)
+                end
+            end
+            rows = vec([(i - 1) * d_sym + s for s in s_sel, i in i_sel])
+            @test isapprox(H_weak[rows, I_sel], manual; atol=1e-6)
+        end
+
+        # No upstream counterpart — added by the port. The ×2 above is a claim about
+        # the symmetric basis, so pin it without reference to that basis: expand the
+        # strong symmetric coefficients into the full (0, 2) basis, pair them with the
+        # full Gram, and the result must be the full weak Hessian
+        # ⟨H(φ_I), φ_i dx_j ⊗ dx_k⟩ — which carries no such convention.
+        @testset "hessian_02_sym_weak agrees with the full (0,2) Hessian" begin
+            c_sym = gram_inv(tensor02sym_space(dg)) * weak(hessian(dg))
+            c_full = reduce(hcat, [expand_symmetric_tensor_coeffs(c_sym[:, I], n1, dim)
+                                   for I in 1:n0])
+
+            lhs = gram(tensor02_space(dg)) * c_full
+            rhs = hessian_02_weak(u, H, mu, n1)
+            rows = flat_idx(sample_indices(n1), dim * dim)
+            cols = sample_indices(n0)
+            @test isapprox(lhs[rows, cols], rhs[rows, cols]; atol=1e-6)
         end
 
         # Strong/weak relation: G · (G⁻¹ · H_weak) == H_weak.