Merge pull request 'test: Java golden-value oracles for the five DCE math cores + P1-2/P1-3 fixes' (#27) from tests/java-golden-oracles into main

Reviewed-on: #27
2026-05-30 08:01:01 +00:00
parent 18b9c61492 ba83974525
commit 03e8000ac6
9 changed files with 580 additions and 61 deletions
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -282,7 +282,7 @@ Expected results: `test-fast` + `quality-gates` green on every push, 0 skipped,
 ## Release state

 Current release: **v0.10.0** (tag on `main`, released 2026-05-26 — the "reviewer-ready" release).
-Phases 1–9a complete, Phase 8b-Lite CGAL API surface complete (all 5 DCE models reachable via `<CGAL/Discrete_*.h>`), Phase 9b block-FD HyperIdeal Hessian shipped (~96× speed-up).  v0.10.0 added: 100 % Doxygen public-API coverage (396/396 symbols, 0 warnings), a 14-gate structural quality suite (4 required in CI), the reviewer materials package (`doc/reviewer/`), `output_uv_map` for 4 of 5 DCE entries, six new roadmap phases + three RESEARCH phases, and a six-mode compile-time workflow matrix.  Numbers at release: **259/259 tests pass, 0 skipped** (23 non-CGAL + 236 CGAL).  Next planned milestones: Phase 9c (4g-polygon, genus g > 1) and Phase 9b-analytic (Schläfli identity).  See `doc/release-policy.md` for the version-tag policy, `CHANGELOG.md` for full release notes, and `doc/roadmap/phases.md` for the phase plan.
+Phases 1–9a complete, Phase 8b-Lite CGAL API surface complete (all 5 DCE models reachable via `<CGAL/Discrete_*.h>`), Phase 9b block-FD HyperIdeal Hessian shipped (~96× speed-up).  v0.10.0 added: 100 % Doxygen public-API coverage (396/396 symbols, 0 warnings), a 14-gate structural quality suite (4 required in CI), the reviewer materials package (`doc/reviewer/`), `output_uv_map` for 4 of 5 DCE entries, six new roadmap phases + three RESEARCH phases, and a six-mode compile-time workflow matrix.  Numbers (single source of truth = `doc/api/tests.md`): **272/272 tests pass, 0 skipped** (26 non-CGAL + 246 CGAL).  Next planned milestones: Phase 9c (4g-polygon, genus g > 1) and Phase 9b-analytic (Schläfli identity).  See `doc/release-policy.md` for the version-tag policy, `CHANGELOG.md` for full release notes, and `doc/roadmap/phases.md` for the phase plan.

 ## Phase 8 CGAL-package decisions (frozen 2026-05-19)

@@ -302,14 +302,18 @@ The 2026-05-21 audit corrected four mis-labels (now fixed): `InversiveDistanceFu

 ### Java↔C++ math-correctness audit (2026-05-29)

-Full line-by-line audit of all math-critical headers against `de.varylab.discreteconformal.*` lives in **`doc/reviewer/java-port-audit.md`** (read it before re-investigating any of these).  All 11 findings are resolved or noted; the four that needed code changes are now ✅ FIXED and **all CGAL tests pass** (240 after the Euclidean holonomy/τ end-to-end tests + the spherical edge-DOF closed-form oracle landed on top — see `doc/api/tests.md`):
+Full line-by-line audit of all math-critical headers against `de.varylab.discreteconformal.*` lives in **`doc/reviewer/java-port-audit.md`** (read it before re-investigating any of these).  All 11 findings are resolved or noted; the four that needed code changes are now ✅ FIXED and **all CGAL tests pass** (246 after the Euclidean holonomy/τ end-to-end tests, the spherical edge-DOF closed-form oracle, and the Java golden-value oracle suites — incl. three *full-mesh* oracles driving the real `EuclideanCyclicFunctional`/`SphericalFunctional` on a shared tetrahedron, one of them an edge-DOF gradient oracle — landed on top; see `doc/api/tests.md`):

 - **Finding 3** (`spherical_functional.hpp`) — spherical edge-DOF now uses Java's *replacement* parameterization (`Λ = λ_e` when the edge is a DOF, via helper `spher_eff_lambda`); edge gradient is `α_opp⁺ + α_opp⁻ − θ_e` (dropped the extra `−(S_f⁺+S_f⁻)/2` term).  Vertex-only path is bit-for-bit unchanged.
 - **Finding 4** (`spherical_hessian.hpp`) — added an always-compiled `throw std::logic_error` edge-DOF guard (mirrors Finding 2 for Euclidean).
 - **Finding 6** (`period_matrix.hpp`) — `compute_period_matrix` now calls the faithful `normalizeModulus` (matches Java oracle: `0 ≤ Re ≤ ½`, `Im ≥ 0`, `|τ| ≥ 1`); `reduce_to_fundamental_domain` retained for the canonical SL(2,ℤ) domain.
 - **Finding 9** (`inversive_distance_functional.hpp`) — degenerate-face gradient now uses limiting angles instead of skipping (mirrors Finding 1); the genuinely-non-real `l*sq <= 0` skip is kept.

-**Open follow-ups (tests only, not bugs):** the edge-DOF paths for Findings 3 & 4 are verified for compilation + the unchanged vertex-only case but **not yet against a Java oracle at the solution level** — missing-test items 4, 7, 10 in the audit doc (golden-value oracle tests, period-matrix sign-convention test, inversive-distance degenerate test).  Build/test reminder: the CGAL build dir used for this audit is **`build-cgal`** at the *repo root* (not under `code/`), target `conformallab_cgal_tests`, filter `ctest -R '^cgal\.'`.
+**Java golden-value oracles (2026-05-29, P0):** five suites now pin the C++ pure-math core bit-for-bit (1e-12) against the compiled Java library (openjdk 17, real `Clausen.Л` / `HyperIdealUtility` / `DiscreteEllipticUtility.normalizeModulus`): `HyperIdealGoldenJava` (Clausen/Л/ImLi₂, ζ₁₃/₁₄/₁₅/ζ, both tetrahedron-volume formulas), `EuclideanGoldenJava` (angle formula + 2·Л energy), `SphericalGoldenJava` (law-of-cosines angles + β relations + Л energy), and `PeriodMatrix.NormalizeModulus_GoldenJava` (audit missing-test item 7, ✅ done).  Oracle harnesses live in `/tmp/oracle/*.java` (recipe in the audit doc).
+
+**Full-mesh oracles (2026-05-29):** `EuclideanGoldenJava.FullMeshGradientAndEnergy_Tetrahedron` and `SphericalGoldenJava.FullMeshGradientAndEnergy_Tetrahedron` drive the *real* `EuclideanCyclicFunctional` / `SphericalFunctional` on a shared tetrahedron (`/tmp/oracle/{tet.obj,EucMeshOracle.java,SphereMeshOracle.java}`) and pin both the per-vertex gradient (`Θ−Σα`) and `ΔE = E(x)−E(0)` to 1e-12 — the energy cross-checks C++'s Gauss-Legendre *path integral* against Java's *closed-form* functional (audit missing-test item 5, ✅ done).  **Finding surfaced:** the spherical oracle must call Java's raw `conformalEnergyAndGradient`, NOT `evaluate()` — the latter pre-runs a 1-D Brent maximization over the global-scale gauge (`maximizeInNegativeDirection`), which C++ deliberately factors into the Newton solver's `spherical_gauge_shift` instead (both correct; different factoring).
+
+**Open follow-ups (tests only, not bugs):** (a) the spherical edge-DOF *gradient* is now Java-oracle'd at the solution level (`SphericalGoldenJava.FullMeshEdgeDofGradient_Tetrahedron` — vertex + edge components vs raw `conformalEnergyAndGradient`, locking Finding 3); the only remaining edge-DOF gap is *Newton-to-convergence* with edge DOFs, blocked by the Finding-4 spherical-Hessian guard (needs an FD-Hessian or guard relaxation first) — a solver feature, not a correctness gap; (b) inversive-distance degenerate-face robustness (item 10) has no Java oracle because the functional is research with no Java parent — only a NaN-free limiting-angle regression test applies.  Build/test reminder: the CGAL build dir used for this audit is **`build-cgal`** at the *repo root* (not under `code/`), target `conformallab_cgal_tests`, filter `ctest -R '^cgal\.'`.

 ## Documentation map

@@ -422,7 +426,7 @@ Every turn re-sends the whole conversation, so context accumulation is the large

 ### Token hygiene — command discipline (Tier 2)

- **Scope + filter together.** Path-scope searches (`grep -rn "newton_" code/include/`, not repo-wide). Filter test output (`ctest -R "cgal.NewtonSolver" --output-on-failure`, or `--gtest_filter` for one suite) instead of dumping all 259 results.
+- **Scope + filter together.** Path-scope searches (`grep -rn "newton_" code/include/`, not repo-wide). Filter test output (`ctest -R "cgal.NewtonSolver" --output-on-failure`, or `--gtest_filter` for one suite) instead of dumping all 272 results.
 - **Never read raw logs into context.** Redirect to a file, then `grep`/`tail` it. With `run_in_background`, read the output file selectively rather than pulling it whole.
 - **Don't re-read a file you just edited.** `Edit` errors on stale state — the harness tracks it; a `Read`-back after a successful edit is wasted context.
 - **Reference by line number** (`newton_solver.hpp:147`) instead of re-pasting code blocks.
--- a/code/include/layout.hpp
+++ b/code/include/layout.hpp
@@ -176,6 +176,14 @@ struct HolonomyData {
    std::vector<Eigen::Vector2d> translations;     ///< Euclidean / spherical translation per cut edge.
    std::vector<MobiusMap>       mobius_maps;      ///< Hyperbolic Möbius isometry per cut edge (Phase 7).
    std::vector<std::size_t>     cut_edge_indices; ///< Index (in the cut-graph edge list) of each holonomy entry.
+
+    /// Residual rotation |arg(a)| (radians) of the Euclidean deck isometry
+    /// z ↦ a·z + b per cut edge.  Zero for a perfectly flat cone metric;
+    /// a non-zero value signals under-convergence or a genuine cone-angle
+    /// defect, in which case `translations[i]` (the isometry's translation
+    /// part b) is only an approximate lattice generator.  Empty for the
+    /// hyperbolic (Möbius) path.
+    std::vector<double>          residual_rotation;
 };

 // ── Internal helpers ──────────────────────────────────────────────────────────
@@ -490,11 +498,23 @@ inline void set_root_huv_2d(
 // different locations.  The rigid motion identifying those two copies is the
 // deck transformation of the generator loop (cut edge + tree path) — i.e. the
 // holonomy.  For a flat cone metric (Θ ≡ 2π) the linear part is trivial, so the
-// holonomy is the pure translation given by the displacement of the shared
-// edge's midpoint between the two developments.  That translation is exactly
-// the lattice generator ω consumed by compute_period_matrix().
+// holonomy is the pure translation that identifies the two developments of the
+// shared edge.  We recover it as a full rigid motion z ↦ a·z + b fitted to the
+// correspondence (Bs,Bt) ↦ (As,At): the shared edge endpoints S,T as developed
+// in face B map to the same endpoints as developed in face A.  Because both
+// faces use the same edge length, |a| = 1 and the fit is an exact
+// orientation-preserving isometry.  Its translation part b is the lattice
+// generator ω consumed by compute_period_matrix(); the rotation magnitude
+// |arg(a)| is returned as a convergence diagnostic.  For a perfectly flat cone
+// metric (Θ ≡ 2π) a = 1 and b reduces to the midpoint displacement midA − midB
+// (the previous formula), so converged flat tori are bit-for-bit unchanged.
+struct EuclideanHolonomyResult {
+    std::vector<Eigen::Vector2d> translations;       ///< ω_i = translation part b
+    std::vector<double>          residual_rotation;  ///< |arg(a_i)| per cut edge
+};
+
 template <typename EdgeLenFn>
-inline std::vector<Eigen::Vector2d> euclidean_holonomy(
+inline EuclideanHolonomyResult euclidean_holonomy(
    const ConformalMesh& mesh,
    const CutGraph&      cut,
    EdgeLenFn&&          edge_len)
@@ -565,9 +585,10 @@ inline std::vector<Eigen::Vector2d> euclidean_holonomy(
    for (auto f : mesh.faces())
        if (!face_done[static_cast<std::size_t>(f.idx())]) develop(f);

-    // ── Holonomy translation per cut edge ─────────────────────────────────────
-    std::vector<Eigen::Vector2d> omega;
-    omega.reserve(cut.cut_edge_indices.size());
+    // ── Holonomy isometry per cut edge ─────────────────────────────────────────
+    EuclideanHolonomyResult out;
+    out.translations.reserve(cut.cut_edge_indices.size());
+    out.residual_rotation.reserve(cut.cut_edge_indices.size());
    for (std::size_t ce : cut.cut_edge_indices) {
        Edge_index e = *std::next(mesh.edges().begin(), static_cast<std::ptrdiff_t>(ce));
        Halfedge_index h  = mesh.halfedge(e);
@@ -575,7 +596,8 @@ inline std::vector<Eigen::Vector2d> euclidean_holonomy(
        if (mesh.is_border(h) || mesh.is_border(ho)
            || !face_done[static_cast<std::size_t>(mesh.face(h).idx())]
            || !face_done[static_cast<std::size_t>(mesh.face(ho).idx())]) {
-            omega.push_back(Eigen::Vector2d::Zero());
+            out.translations.push_back(Eigen::Vector2d::Zero());
+            out.residual_rotation.push_back(0.0);
            continue;
        }
        // Face A = face(h) places edge e as halfedge h:  S=source(h), T=target(h)
@@ -585,12 +607,15 @@ inline std::vector<Eigen::Vector2d> euclidean_holonomy(
        //   target(ho)=S  →  S=pos of source(next(ho)), T=pos of source(ho)
        C Bt = hpos[static_cast<std::size_t>(ho.idx())];
        C Bs = hpos[static_cast<std::size_t>(mesh.next(ho).idx())];
-        C midA = 0.5 * (As + At);
-        C midB = 0.5 * (Bs + Bt);
-        C w    = midA - midB;     // pure translation for a flat cone metric
-        omega.push_back(Eigen::Vector2d(w.real(), w.imag()));
+        // Fit the deck isometry g(z) = a·z + b with (Bs,Bt) ↦ (As,At).
+        // |a| = 1 exactly (shared edge length); b is the translation part.
+        C edgeB = Bt - Bs;
+        C a = (std::abs(edgeB) > 1e-300) ? (At - As) / edgeB : C(1.0, 0.0);
+        C b = As - a * Bs;
+        out.translations.push_back(Eigen::Vector2d(b.real(), b.imag()));
+        out.residual_rotation.push_back(std::abs(std::arg(a)));
    }
-    return omega;
+    return out;
 }

 } // namespace detail
@@ -719,7 +744,9 @@ inline Layout2D euclidean_layout(
    if (cut && holonomy) {
        holonomy->cut_edge_indices = cut->cut_edge_indices;
        holonomy->mobius_maps.clear();
-        holonomy->translations = detail::euclidean_holonomy(mesh, *cut, edge_len);
+        auto eh = detail::euclidean_holonomy(mesh, *cut, edge_len);
+        holonomy->translations      = std::move(eh.translations);
+        holonomy->residual_rotation = std::move(eh.residual_rotation);

        // Preserve the per-cut-edge seam UV in halfedge_uv (texture atlas), as
        // before, so HalfedgeUV-based tests still see the seam-crossing layout.
--- a/code/include/period_matrix.hpp
+++ b/code/include/period_matrix.hpp
@@ -114,11 +114,18 @@ inline std::complex<double> reduce_to_fundamental_domain(std::complex<double> ta
 // is_in_fundamental_domain — check membership in F with tolerance tol.
 // ─────────────────────────────────────────────────────────────────────────────
 /// `true` iff `τ` lies inside the standard SL(2,ℤ) fundamental domain
-/// with tolerance `tol`.
+/// `F = { Im τ > 0, −½ ≤ Re τ < ½, |τ| ≥ 1 }` with tolerance `tol`.
+///
+/// This is the half-open domain produced by `reduce_to_fundamental_domain`
+/// (the right boundary `Re τ = +½ ≡ −½` is excluded via `T`).  It is NOT the
+/// mirror-folded `0 ≤ Re τ ≤ ½` domain produced by `normalizeModulus` (used
+/// inside `compute_period_matrix`); a τ with `Re τ < 0` is a legitimate member
+/// of `F` here but would be folded to `Re ≥ 0` by `normalizeModulus`.
 inline bool is_in_fundamental_domain(std::complex<double> tau, double tol = 1e-9)
 {
    if (tau.imag() <= 0.0) return false;
-    if (std::abs(tau.real()) > 0.5 + tol) return false;
+    if (tau.real() < -0.5 - tol) return false;  // left boundary closed
+    if (tau.real() >  0.5 - tol) return false;  // right boundary Re = +½ excluded (≡ −½)
    if (std::abs(tau) < 1.0 - tol) return false;
    return true;
 }
--- a/code/tests/cgal/test_euclidean_functional.cpp
+++ b/code/tests/cgal/test_euclidean_functional.cpp
@@ -26,6 +26,7 @@
 #include "euclidean_geometry.hpp"
 #include "euclidean_functional.hpp"
 #include "euclidean_hessian.hpp"
+#include "clausen.hpp"
 #include <gtest/gtest.h>
 #include <cmath>
 #include <vector>
@@ -290,3 +291,117 @@ TEST(EuclideanFunctional, GradientCheck_MixedPinnedVertices)
    EXPECT_TRUE(gradient_check_euclidean(mesh, x, maps))
        << "Gradient check failed for mixed pinned/variable vertices";
 }
+
+// ─────────────────────────────────────────────────────────────────────────────
+// Golden-value oracle — pin the Euclidean angle formula and the 2·Л(α) energy
+// term bit-for-bit against the upstream Java reference (EuclideanCyclicFunctional
+// .triangleEnergyAndAlphas, lines 341-361), captured by running the compiled Java
+// library (openjdk 17) with the real de.varylab…Clausen.Л on these exact edge
+// lengths.  Companion to HyperIdealGoldenJava and SphericalGoldenJava: locks the
+// absolute angle/energy values against an independent implementation, catching
+// silent index/sign drift the curl-free path-integral gradient check cannot see.
+//
+// To regenerate: /tmp/oracle/EucOracle.java.  Values are Java printf %.17g.
+// ─────────────────────────────────────────────────────────────────────────────
+TEST(EuclideanGoldenJava, AngleAndLobachevskyEnergyFromLengths)
+{
+    auto check = [](double l12, double l23, double l31,
+                    double a1_g, double a2_g, double a3_g, double L_g) {
+        auto fa = euclidean_angles_from_lengths(l12, l23, l31);
+        EXPECT_TRUE(fa.valid);
+        EXPECT_NEAR(fa.alpha1, a1_g, 1e-12);
+        EXPECT_NEAR(fa.alpha2, a2_g, 1e-12);
+        EXPECT_NEAR(fa.alpha3, a3_g, 1e-12);
+        const double Lterm = 2.0 * Lobachevsky(fa.alpha1)
+                           + 2.0 * Lobachevsky(fa.alpha2)
+                           + 2.0 * Lobachevsky(fa.alpha3);
+        EXPECT_NEAR(Lterm, L_g, 1e-12);
+    };
+    check(1.0, 1.2, 0.9,
+          1.3637649752769678, 0.82416964552030680, 0.95365803279251860,
+          1.9456273836230942);
+    check(1.0, 1.0, 1.0,
+          1.0471975511965979, 1.0471975511965979, 1.0471975511965979,
+          2.0298832128193070);
+}
+
+// ─────────────────────────────────────────────────────────────────────────────
+// FULL-MESH golden oracle — the strongest cross-check: drives the REAL upstream
+// EuclideanCyclicFunctional (openjdk 17) on a tetrahedron loaded from a shared
+// OBJ (identical topology + geometry to make_tetrahedron()), and pins BOTH the
+// per-vertex gradient G_v = Θ_v − Σα AND the energy difference ΔE = E(x) − E(0)
+// bit-for-bit against it.
+//
+// This closes audit missing-test item 5 (full-mesh energy + gradient at a known
+// x).  It is genuinely independent of the C++ implementation in two ways:
+//   • the gradient is the upstream library's own analytic gradient (not an FD
+//     check, which only proves curl-freeness of the C++ self-consistent energy);
+//   • the energy is Java's CLOSED-FORM functional value, whereas C++ computes it
+//     as a Gauss-Legendre PATH INTEGRAL of its gradient — two different methods
+//     that must agree on ΔE (the initialEnergy constant and the φ·λ⁰ term cancel
+//     in the difference, and there are no edge DOFs here).
+//
+// Setup parity (verified against UnwrapUtility.prepareInvariantDataEuclidean):
+//   closed mesh, ALL 4 vertices variable (no pin), Θ_v = 2π, no edge DOFs,
+//   λ°_e = 2·log(|p_i − p_j|), per-vertex u(P) = 0.10·X − 0.07·Y + 0.13·Z.
+//
+// To regenerate: /tmp/oracle/{tet.obj,EucMeshOracle.java}.  Values are Java %.17g.
+// ─────────────────────────────────────────────────────────────────────────────
+TEST(EuclideanGoldenJava, FullMeshGradientAndEnergy_Tetrahedron)
+{
+    constexpr double TWO_PI = 2.0 * 3.14159265358979323846264338328;
+
+    auto mesh = make_tetrahedron();   // same 4 vertices as /tmp/oracle/tet.obj
+    auto maps = setup_euclidean_maps(mesh);
+    compute_euclidean_lambda0_from_mesh(mesh, maps);
+
+    // All four vertices are DOFs (Java prepareInvariantData makes every interior
+    // vertex variable on a closed mesh); Θ_v = 2π; no edge DOFs.
+    int idx = 0;
+    for (auto v : mesh.vertices()) {
+        maps.v_idx[v]   = idx++;
+        maps.theta_v[v] = TWO_PI;
+    }
+
+    auto u_of = [](const Point3& p) {
+        return 0.10 * p.x() - 0.07 * p.y() + 0.13 * p.z();
+    };
+
+    std::vector<double> x(static_cast<std::size_t>(idx), 0.0);
+    for (auto v : mesh.vertices())
+        x[static_cast<std::size_t>(maps.v_idx[v])] = u_of(mesh.point(v));
+
+    auto G = euclidean_gradient(mesh, x, maps);
+
+    // Java golden gradients keyed by vertex position (order-independent lookup).
+    struct GoldRow { double X, Y, Z, G; };
+    const GoldRow gold[4] = {
+        { 1,  1,  1, 3.5277511803984396},
+        { 1, -1, -1, 3.2837611905358440},
+        {-1,  1, -1, 2.3437485291237100},
+        {-1, -1,  1, 3.4111097143011780},
+    };
+    for (auto v : mesh.vertices()) {
+        const auto& p = mesh.point(v);
+        const double g = G[static_cast<std::size_t>(maps.v_idx[v])];
+        bool matched = false;
+        for (const auto& row : gold) {
+            if (std::abs(p.x() - row.X) < 1e-9 &&
+                std::abs(p.y() - row.Y) < 1e-9 &&
+                std::abs(p.z() - row.Z) < 1e-9) {
+                EXPECT_NEAR(g, row.G, 1e-12)
+                    << "gradient mismatch at (" << p.x() << "," << p.y()
+                    << "," << p.z() << ")";
+                matched = true;
+                break;
+            }
+        }
+        EXPECT_TRUE(matched) << "unexpected vertex position";
+    }
+
+    // ΔE = E(x) − E(0).  C++ path integral vs Java closed-form functional value.
+    std::vector<double> x0(static_cast<std::size_t>(idx), 0.0);
+    const double dE = euclidean_energy(mesh, x, maps)
+                    - euclidean_energy(mesh, x0, maps);
+    EXPECT_NEAR(dE, 0.15962619236187336, 1e-12);
+}
--- a/code/tests/cgal/test_phase7.cpp
+++ b/code/tests/cgal/test_phase7.cpp
@@ -338,12 +338,52 @@ TEST(PeriodMatrix, ComputePeriodMatrix_UnitSquare)
 TEST(PeriodMatrix, ComputePeriodMatrix_ReducedTau_InFD)
 {
    // ω_1 = (1, 0), ω_2 = (0.5, 0.25)  →  τ = 0.5 + 0.25i
-    // |τ| = sqrt(0.25 + 0.0625) ≈ 0.559 < 1 → needs S step
+    // |τ| = sqrt(0.25 + 0.0625) ≈ 0.559 < 1 → needs S step.
+    // compute_period_matrix reduces with normalizeModulus (Finding 6), whose
+    // mirror-folded target domain is { 0 ≤ Re ≤ ½, Im > 0, |τ| ≥ 1 } — note the
+    // RIGHT boundary Re = +½ is CLOSED here.  This is NOT the half-open SL(2,ℤ)
+    // domain of is_in_fundamental_domain (−½ ≤ Re < ½), which excludes Re = +½
+    // because +½ ≡ −½ under T.  For this input normalizeModulus lands exactly on
+    // τ = ½ + i, so we must check the normalizeModulus domain, not the SL(2,ℤ)
+    // one (asserting is_in_fundamental_domain here would wrongly fail on +½).
    HolonomyData hol;
    hol.translations = { Eigen::Vector2d(1.0, 0.0), Eigen::Vector2d(0.5, 0.25) };
    PeriodData pd = compute_period_matrix(hol, /*reduce=*/true);
    EXPECT_TRUE(pd.in_fundamental_domain);
-    EXPECT_TRUE(is_in_fundamental_domain(pd.tau, 1e-9));
+    const double tol = 1e-9;
+    EXPECT_GT(pd.tau.imag(), 0.0);
+    EXPECT_GE(pd.tau.real(), 0.0 - tol);          // 0 ≤ Re   (mirror fold)
+    EXPECT_LE(pd.tau.real(), 0.5 + tol);          // Re ≤ ½   (closed right edge)
+    EXPECT_GE(std::abs(pd.tau), 1.0 - tol);       // |τ| ≥ 1
+    // Concretely: τ = ½ + i.
+    EXPECT_NEAR(pd.tau.real(), 0.5, 1e-12);
+    EXPECT_NEAR(pd.tau.imag(), 1.0, 1e-12);
+}
+
+// ─────────────────────────────────────────────────────────────────────────────
+// Golden-value oracle — pin normalizeModulus (the τ reduction used by
+// compute_period_matrix, Finding 6) bit-for-bit against the upstream Java
+// reference (de.varylab.discreteconformal.util.DiscreteEllipticUtility.
+// normalizeModulus), captured by calling the compiled Java method (openjdk 17)
+// on these exact τ.  This locks the sign/fold conventions of the SL(2,ℤ)+mirror
+// reduction (0 ≤ Re ≤ ½, Im ≥ 0, |τ| ≥ 1) against an independent implementation,
+// catching drift the existing in-FD membership checks cannot (they only assert
+// the result lies in F, not that it is the SAME representative Java picks).
+//
+// To regenerate: /tmp/oracle/TauOracle.java.  Values are Java printf %.17g.
+// ─────────────────────────────────────────────────────────────────────────────
+TEST(PeriodMatrix, NormalizeModulus_GoldenJava)
+{
+    auto chk = [](double re, double im, double re_g, double im_g) {
+        C n = conformallab::normalizeModulus(C(re, im));
+        EXPECT_NEAR(n.real(), re_g, 1e-12);
+        EXPECT_NEAR(n.imag(), im_g, 1e-12);
+    };
+    chk(0.3, 0.5,  0.11764705882352933, 1.4705882352941178);  // |τ|<1 → S + folds
+    chk(-0.4, 1.3, 0.40000000000000000, 1.3000000000000000);  // Re<0 mirror fold
+    chk(2.7, 0.8,  0.41095890410958880, 1.0958904109589043);  // large Re → T
+    chk(0.1, 2.0,  0.10000000000000000, 2.0000000000000000);  // already in F
+    chk(-1.6, 0.9, 0.41237113402061850, 0.92783505154639180); // T + S + mirror
 }

 // ════════════════════════════════════════════════════════════════════════════
--- a/code/tests/cgal/test_spherical_functional.cpp
+++ b/code/tests/cgal/test_spherical_functional.cpp
@@ -27,6 +27,8 @@
 #include "mesh_builder.hpp"
 #include "spherical_functional.hpp"
 #include "spherical_hessian.hpp"
+#include "spherical_geometry.hpp"
+#include "clausen.hpp"
 #include <gtest/gtest.h>
 #include <cmath>
 #include <vector>
@@ -416,3 +418,237 @@ TEST(SphericalFunctional, GaugeFix_AlreadyAtGaugeReturnsTNearZero)
    EXPECT_NEAR(t, 0.0, 1e-5)
        << "Gauge shift from the symmetric point should be ~0; got " << t;
 }
+
+// ─────────────────────────────────────────────────────────────────────────────
+// Golden-value oracle — pin the spherical law-of-cosines angle formula, the β
+// half-angle relations, and the Lobachevsky energy term bit-for-bit against the
+// upstream Java reference (SphericalFunctional.triangleEnergyAndAlphas, lines
+// 401-453), captured by running the compiled Java library (openjdk 17) with the
+// real de.varylab…Clausen.Л on these exact arc lengths.  Unlike the Schläfli
+// path-integral gradient check (which only verifies curl-freeness), this locks
+// the absolute angle/β/energy values against an independent implementation,
+// catching any silent index/sign/convention drift in the spherical port.
+//
+// To regenerate: /tmp/oracle/SphereOracle.java (recipe in doc/reviewer/
+// java-port-audit.md).  Values are Java printf %.17g.  Tolerance 1e-12.
+//
+// Index map (C++ spherical_angles ↔ Java): with l12=lij, l23=ljk, l31=lki the
+// returned alpha1/alpha2/alpha3 are Java αi/αj/αk (angle opposite ljk/lki/lij).
+// ─────────────────────────────────────────────────────────────────────────────
+TEST(SphericalGoldenJava, AngleBetaEnergyFromLengths)
+{
+    auto check = [](double lij, double ljk, double lki,
+                    double ai_g, double aj_g, double ak_g,
+                    double bi_g, double bj_g, double bk_g, double L_g) {
+        auto fa = spherical_angles(lij, ljk, lki);
+        EXPECT_TRUE(fa.valid);
+        EXPECT_NEAR(fa.alpha1, ai_g, 1e-12);
+        EXPECT_NEAR(fa.alpha2, aj_g, 1e-12);
+        EXPECT_NEAR(fa.alpha3, ak_g, 1e-12);
+
+        const double ai = fa.alpha1, aj = fa.alpha2, ak = fa.alpha3;
+        const double bi = 0.5 * (PI + ai - aj - ak);
+        const double bj = 0.5 * (PI - ai + aj - ak);
+        const double bk = 0.5 * (PI - ai - aj + ak);
+        EXPECT_NEAR(bi, bi_g, 1e-12);
+        EXPECT_NEAR(bj, bj_g, 1e-12);
+        EXPECT_NEAR(bk, bk_g, 1e-12);
+
+        const double Lterm =
+            Lobachevsky(ai) + Lobachevsky(aj) + Lobachevsky(ak) +
+            Lobachevsky(bi) + Lobachevsky(bj) + Lobachevsky(bk) +
+            Lobachevsky(0.5 * (PI - ai - aj - ak));
+        EXPECT_NEAR(Lterm, L_g, 1e-12);
+    };
+
+    // Scalene spherical triangle (all Δ > 0, Δ_ijk < 2π).
+    check(1.0, 1.2, 0.9,
+          1.5305813141072122, 0.99684981272794600, 1.1246069519101605,
+          1.2753586015294494, 0.74162710015018330, 0.86938423933239760,
+          1.3576100185550408);
+    // Equilateral spherical triangle.
+    check(0.7, 0.7, 0.7,
+          1.1225596283199812, 1.1225596283199812, 1.1225596283199812,
+          1.0095165126349062, 1.0095165126349060, 1.0095165126349060,
+          1.6806828584976297);
+}
+
+// ─────────────────────────────────────────────────────────────────────────────
+// FULL-MESH golden oracle (spherical) — drives the REAL upstream SphericalFunctional
+// on the shared tetrahedron (scaled so every arc length stays < π) and pins BOTH
+// the per-vertex gradient G_v = Θ_v − Σα AND ΔE = E(x) − E(0) bit-for-bit.
+//
+// IMPORTANT — the oracle calls Java's `conformalEnergyAndGradient` (the RAW
+// θ−Σα gradient), NOT `evaluate()`: `evaluate()` first runs a 1-D Brent
+// maximization over the global-scale gauge direction (`maximizeInNegativeDirection`),
+// which the C++ `spherical_gradient` deliberately does NOT — C++ factors that
+// gauge into the Newton solver's `spherical_gauge_shift` instead.  Comparing
+// against `evaluate()` would (wrongly) drive every component to ~0.  The raw
+// gradient is the piece that corresponds 1:1 to the C++ functional.
+//
+// Setup parity (UnwrapUtility.prepareInvariantDataHyperbolicAndSpherical, scale):
+//   closed mesh, ALL 4 vertices variable, Θ_v = 2π, no edge DOFs,
+//   λ°_e = 2·log(SCALE·|p_i − p_j|)  (Java uses chord length × scale, whereas the
+//   C++ compute_lambda0_from_mesh helper assumes unit-sphere vertices and uses
+//   ARC length — so we set λ° directly here to match Java exactly),
+//   per-vertex u(P) = 0.10·X − 0.07·Y + 0.13·Z, SCALE = 0.2.
+//
+// To regenerate: /tmp/oracle/{tet.obj,SphereMeshOracle.java}.  Values are Java %.17g.
+// ─────────────────────────────────────────────────────────────────────────────
+TEST(SphericalGoldenJava, FullMeshGradientAndEnergy_Tetrahedron)
+{
+    constexpr double TWO_PI = 2.0 * 3.14159265358979323846264338328;
+    constexpr double SCALE  = 0.2;
+
+    auto mesh = make_tetrahedron();   // same 4 vertices as /tmp/oracle/tet.obj
+    auto maps = setup_spherical_maps(mesh);
+
+    // λ°_e = 2·log(SCALE · chord), matching Java's prepareInvariantData(scale).
+    for (auto e : mesh.edges()) {
+        auto h  = mesh.halfedge(e);
+        const auto& p1 = mesh.point(mesh.source(h));
+        const auto& p2 = mesh.point(mesh.target(h));
+        const double dx = p1.x() - p2.x(), dy = p1.y() - p2.y(), dz = p1.z() - p2.z();
+        const double chord = std::sqrt(dx*dx + dy*dy + dz*dz);
+        maps.lambda0[e] = 2.0 * std::log(SCALE * chord);
+    }
+
+    int idx = 0;
+    for (auto v : mesh.vertices()) {
+        maps.v_idx[v]   = idx++;
+        maps.theta_v[v] = TWO_PI;
+    }
+
+    auto u_of = [](const Point3& p) {
+        return 0.10 * p.x() - 0.07 * p.y() + 0.13 * p.z();
+    };
+
+    std::vector<double> x(static_cast<std::size_t>(idx), 0.0);
+    for (auto v : mesh.vertices())
+        x[static_cast<std::size_t>(maps.v_idx[v])] = u_of(mesh.point(v));
+
+    auto G = spherical_gradient(mesh, x, maps);
+
+    struct GoldRow { double X, Y, Z, G; };
+    const GoldRow gold[4] = {
+        { 1,  1,  1, 2.7671034786104927},
+        { 1, -1, -1, 2.5216834054857546},
+        {-1,  1, -1, 1.5401761310866633},
+        {-1, -1,  1, 2.6518569531861100},
+    };
+    for (auto v : mesh.vertices()) {
+        const auto& p = mesh.point(v);
+        const double g = G[static_cast<std::size_t>(maps.v_idx[v])];
+        bool matched = false;
+        for (const auto& row : gold) {
+            if (std::abs(p.x() - row.X) < 1e-9 &&
+                std::abs(p.y() - row.Y) < 1e-9 &&
+                std::abs(p.z() - row.Z) < 1e-9) {
+                EXPECT_NEAR(g, row.G, 1e-12)
+                    << "gradient mismatch at (" << p.x() << "," << p.y()
+                    << "," << p.z() << ")";
+                matched = true;
+                break;
+            }
+        }
+        EXPECT_TRUE(matched) << "unexpected vertex position";
+    }
+
+    std::vector<double> x0(static_cast<std::size_t>(idx), 0.0);
+    const double dE = spherical_energy(mesh, x, maps)
+                    - spherical_energy(mesh, x0, maps);
+    EXPECT_NEAR(dE, 0.16409141487397116, 1e-12);
+}
+
+// ─────────────────────────────────────────────────────────────────────────────
+// FULL-MESH EDGE-DOF golden oracle (spherical, Finding 3) — the missing solution-
+// level Java oracle for the edge-DOF gradient path.  Makes two opposite edges of
+// the shared tetrahedron variable (replacement parameterization: λ_e = x[e_idx],
+// no u-shift) and pins BOTH the vertex gradient (Θ−Σα) AND the edge gradient
+// (α_opp⁺ + α_opp⁻ − θ_e, θ_e = π) bit-for-bit against the REAL upstream
+// SphericalFunctional.conformalEnergyAndGradient (raw gradient, not evaluate()).
+//
+// This is a pure GRADIENT oracle (no ΔE): with an edge DOF, x = 0 means λ_e = 0
+// ⇒ spherical arc length = π (degenerate), so the path-integral-from-origin
+// energy reference is ill-defined — the gradient at a fixed non-degenerate x is
+// the unambiguous Finding-3 quantity.  It does NOT touch the spherical Hessian,
+// so it is unaffected by the Finding-4 edge-DOF Hessian guard.
+//
+// To regenerate: /tmp/oracle/SphereEdgeOracle.java.  Values are Java %.17g.
+// ─────────────────────────────────────────────────────────────────────────────
+TEST(SphericalGoldenJava, FullMeshEdgeDofGradient_Tetrahedron)
+{
+    constexpr double TWO_PI = 2.0 * 3.14159265358979323846264338328;
+    constexpr double PI_C   = 3.14159265358979323846264338328;
+    constexpr double SCALE  = 0.2;
+
+    auto mesh = make_tetrahedron();
+    auto maps = setup_spherical_maps(mesh);
+
+    for (auto e : mesh.edges()) {
+        auto h  = mesh.halfedge(e);
+        const auto& p1 = mesh.point(mesh.source(h));
+        const auto& p2 = mesh.point(mesh.target(h));
+        const double dx = p1.x()-p2.x(), dy = p1.y()-p2.y(), dz = p1.z()-p2.z();
+        maps.lambda0[e] = 2.0 * std::log(SCALE * std::sqrt(dx*dx + dy*dy + dz*dz));
+        maps.theta_e[e] = PI_C;
+        maps.e_idx[e]   = -1;
+    }
+
+    int idx = 0;
+    for (auto v : mesh.vertices()) { maps.v_idx[v] = idx++; maps.theta_v[v] = TWO_PI; }
+
+    auto pos = [&](Vertex_index v){ return mesh.point(v); };
+    auto same = [](const Point3& p, double X, double Y, double Z) {
+        return std::abs(p.x()-X)<1e-9 && std::abs(p.y()-Y)<1e-9 && std::abs(p.z()-Z)<1e-9;
+    };
+    // Edge connects positions (aX,aY,aZ)–(bX,bY,bZ) in either order?
+    auto connects = [&](Edge_index e, const double a[3], const double b[3]) {
+        auto h = mesh.halfedge(e);
+        const Point3& s = pos(mesh.source(h));
+        const Point3& t = pos(mesh.target(h));
+        return (same(s,a[0],a[1],a[2]) && same(t,b[0],b[1],b[2])) ||
+               (same(s,b[0],b[1],b[2]) && same(t,a[0],a[1],a[2]));
+    };
+    const double A[3]={1,1,1}, B[3]={1,-1,-1}, C[3]={-1,1,-1}, D[3]={-1,-1,1};
+
+    // Make edges A–B and C–D variable (replacement parameterization).
+    int edof = idx;  // edge DOFs start after the 4 vertex DOFs
+    Edge_index eAB, eCD;
+    for (auto e : mesh.edges()) {
+        if (connects(e, A, B)) { maps.e_idx[e] = edof++; eAB = e; }
+        else if (connects(e, C, D)) { maps.e_idx[e] = edof++; eCD = e; }
+    }
+    ASSERT_EQ(edof, 6);  // 4 vertex + 2 edge DOFs
+
+    auto u_of = [](const Point3& p){ return 0.10*p.x() - 0.07*p.y() + 0.13*p.z(); };
+
+    std::vector<double> x(6, 0.0);
+    for (auto v : mesh.vertices())
+        x[static_cast<std::size_t>(maps.v_idx[v])] = u_of(mesh.point(v));
+    // Edge DOF value = λ⁰_e + 0.1 (same perturbation as the Java oracle).
+    x[static_cast<std::size_t>(maps.e_idx[eAB])] = maps.lambda0[eAB] + 0.1;
+    x[static_cast<std::size_t>(maps.e_idx[eCD])] = maps.lambda0[eCD] + 0.1;
+
+    auto G = spherical_gradient(mesh, x, maps);
+
+    // Vertex gradients keyed by position.
+    struct VRow { double X, Y, Z, G; };
+    const VRow vg[4] = {
+        { 1,  1,  1, 2.5036751008175546},
+        { 1, -1, -1, 2.2303362601050205},
+        {-1,  1, -1, 1.8463710314817632},
+        {-1, -1,  1, 2.7499719487341570},
+    };
+    for (auto v : mesh.vertices()) {
+        const auto& p = mesh.point(v);
+        const double g = G[static_cast<std::size_t>(maps.v_idx[v])];
+        bool m = false;
+        for (auto& r : vg)
+            if (same(p, r.X, r.Y, r.Z)) { EXPECT_NEAR(g, r.G, 1e-12); m = true; break; }
+        EXPECT_TRUE(m);
+    }
+    // Edge gradients: A–B and C–D (Java golden values).
+    EXPECT_NEAR(G[static_cast<std::size_t>(maps.e_idx[eAB])], -0.35189517043413690, 1e-12);
+    EXPECT_NEAR(G[static_cast<std::size_t>(maps.e_idx[eCD])], -0.44101986058895950, 1e-12);
+}
--- a/code/tests/test_hyper_ideal_utility.cpp
+++ b/code/tests/test_hyper_ideal_utility.cpp
@@ -4,10 +4,12 @@
 // Port of de.varylab.discreteconformal.functional.HyperIdealUtilityTest (Java/JUnit).

 #include "hyper_ideal_utility.hpp"
+#include "hyper_ideal_geometry.hpp"
 #include "clausen.hpp"

 #include <gtest/gtest.h>
 #include <cmath>
+#include <complex>

 using conformallab::calculateTetrahedronVolume;
 using conformallab::calculateTetrahedronVolumeWithIdealVertexAtGamma;
@@ -93,3 +95,43 @@ TEST(HyperIdealUtilityTest, CompareGeneralAndIdealFormulaCase2) {
    double V  = calculateTetrahedronVolume(bi, bj, bk, ai, aj, ak);
    EXPECT_NEAR(Ve, V, 1e-12);
 }
+
+// ─────────────────────────────────────────────────────────────────────────────
+// Golden-value oracle tests — pin the C++ pure-math core bit-for-bit against the
+// upstream Java reference (de.varylab.discreteconformal.functional.{Clausen,
+// HyperIdealUtility}), captured by running the compiled Java library (openjdk 17)
+// on these exact inputs.  Unlike the FD gradient checks (which only verify
+// curl-freeness, since the C++ energy is the path-integral of its own gradient),
+// these lock the absolute values of the math-critical helpers against an
+// independent implementation, catching any silent convention/formula drift.
+//
+// To regenerate: see doc/reviewer/java-port-audit.md (oracle harness recipe).
+// Values are Java's Double.toString output (shortest round-trip).  Tolerance is
+// 1e-12 (well above the ~1e-15 inter-platform libm divergence for these ranges).
+// ─────────────────────────────────────────────────────────────────────────────
+TEST(HyperIdealGoldenJava, ClausenLobachevskyImLi2) {
+    EXPECT_NEAR(conformallab::clausen2(0.7),  0.954448086482735,  1e-12);
+    EXPECT_NEAR(conformallab::clausen2(2.5),  0.4335982032355327, 1e-12);
+    EXPECT_NEAR(conformallab::clausen2(-1.3), -0.9897032532295984, 1e-12);
+    EXPECT_NEAR(Lobachevsky(0.9),             0.4121734067662043, 1e-12);
+    EXPECT_NEAR(conformallab::ImLi2(std::complex<double>(0.3, 0.4)),
+                                              0.46136289181910894, 1e-12);
+}
+
+TEST(HyperIdealGoldenJava, ZetaFamily) {
+    EXPECT_NEAR(conformallab::zeta13(0.5, 0.7, 0.9), 2.663195966482385,  1e-12);
+    EXPECT_NEAR(conformallab::zeta14(0.4, 0.8),      1.8262295633065202, 1e-12);
+    EXPECT_NEAR(conformallab::zeta15(0.6),           2.216976794676588,  1e-12);
+    EXPECT_NEAR(conformallab::zeta(0.5, 0.7, 0.9),   1.6156519307269948, 1e-12);
+}
+
+TEST(HyperIdealGoldenJava, TetrahedronVolumes) {
+    // Generic non-degenerate generalized hyperbolic tetrahedron.
+    EXPECT_NEAR(calculateTetrahedronVolume(0.6, 0.7, 0.8, 0.5, 0.9, 0.4),
+                2.9074633382516435, 1e-12);
+    // One ideal vertex at gamma, non-symmetric angles (Java test2 inputs).
+    EXPECT_NEAR(calculateTetrahedronVolumeWithIdealVertexAtGamma(
+                    0.6623267054958116, 1.437248992086214, 1.0420169560077686,
+                    0.6896178197389236, 0.5195634857410114, 0.6304500578493993),
+                2.0459750326926214, 1e-12);
+}
--- a/doc/api/tests.md
+++ b/doc/api/tests.md
@@ -7,14 +7,14 @@ Pure-math tests, only Eigen required. Covers Java utilities ported in Phase 1–
 | File | What it tests |
 |---|---|
 | `test_clausen.cpp` | Clausen Cl₂, Lobachevsky Л, ImLi₂ — values at known points |
-| `test_hyper_ideal_utility.cpp` | Tetrahedron volumes (Meyerhoff / Kolpakov–Mednykh) |
+| `test_hyper_ideal_utility.cpp` | Tetrahedron volumes (Meyerhoff / Kolpakov–Mednykh) + Java golden-value oracle (Clausen/Л/ImLi₂, ζ₁₃/₁₄/₁₅/ζ, both volume formulas) |
 | `test_matrix_utility.cpp` | Matrix helpers |
 | `test_surface_curve_utility.cpp` | Surface curve utilities |
 | `test_discrete_elliptic_utility.cpp` | Discrete elliptic functions |
 | `test_p2_utility.cpp` | P2 projective utilities |
 | `test_hyper_ideal_visualization_utility.cpp` | Poincaré disk projection, circumcircle |

-**Total: 23 tests, 0 skipped.**
+**Total: 26 tests, 0 skipped.**

 ---

@@ -30,7 +30,9 @@ All tests have CTest prefix `cgal.` (set via `TEST_PREFIX "cgal."` in CMakeLists
 | `ConformalMeshValidity` | `test_conformal_mesh.cpp` | 1 | CGAL validity for all factory meshes |
 | `HyperIdealFunctional` | `test_hyper_ideal_functional.cpp` | 7 | FD gradient checks + Hessian symmetry |
 | `SphericalFunctional` | `test_spherical_functional.cpp` | 13 | Angle formula + gradient + gauge-fix + cross-module Hessian check + closed-form edge-DOF oracle |
+| `SphericalGoldenJava` | `test_spherical_functional.cpp` | 3 | Java golden oracles: (1) law-of-cosines angles + β relations + Л energy term; (2) full-mesh tetrahedron — per-vertex gradient + ΔE; (3) full-mesh **edge-DOF** gradient (vertex + edge components, Finding 3) vs raw `conformalEnergyAndGradient` |
 | `EuclideanFunctional` | `test_euclidean_functional.cpp` | 12 | Angle formula + gradient + cross-module Hessian check |
+| `EuclideanGoldenJava` | `test_euclidean_functional.cpp` | 2 | Java golden oracles: (1) angle formula + 2·Л(α) energy term; (2) full-mesh tetrahedron — per-vertex gradient + ΔE vs real `EuclideanCyclicFunctional` |
 | `EuclideanHessian` | `test_euclidean_hessian.cpp` | 9 | Cotangent Laplacian structure, FD agreement, PSD, null space |
 | `SphericalHessian` | `test_spherical_hessian.cpp` | 8 | Derivative correctness, NSD at equilibrium |
 | `NewtonSolver` | `test_newton_solver.cpp` | 11 | Convergence: Euclidean ×3, Spherical ×4, HyperIdeal ×4 |
@@ -48,7 +50,7 @@ All tests have CTest prefix `cgal.` (set via `TEST_PREFIX "cgal."` in CMakeLists
 | `HalfedgeUV` | `test_phase7.cpp` | 4 | Size = number of half-edges, seam consistency, boundary half-edges = 0 |
 | `PriorityBFS` | `test_phase7.cpp` | 3 | Success, no seam on open meshes, all vertices placed |
 | `NormaliseEuclidean` | `test_phase7.cpp` | 2 | UV centroid = 0, halfedge_uv centroid = 0 |
-| `PeriodMatrix` | `test_phase7.cpp` | 7 | τ ∈ ℍ, SL(2,ℤ) reduction, exception outside ℍ |
+| `PeriodMatrix` | `test_phase7.cpp` | 8 | τ ∈ ℍ, SL(2,ℤ) reduction, exception outside ℍ, Java `normalizeModulus` golden oracle |
 | `FundamentalDomain` | `test_phase7.cpp` | 7 | Genus-1 parallelogram CCW, generators, g > 1 empty |
 | `TilingCopy/Neighbourhood` | `test_phase7.cpp` | 4 | Translation correct, tile count |
 | `HolonomyEndToEnd` | `test_phase7.cpp` | 3 | Full pipeline τ on tori of revolution (4×4, hex 6×6, 8×8) vs analytic modulus i·√(R²−r²)/r |
@@ -67,7 +69,7 @@ All tests have CTest prefix `cgal.` (set via `TEST_PREFIX "cgal."` in CMakeLists
 | `NewtonPhase9a` | `test_newton_phase9a.cpp` | 7 | Phase 9a-Newton — convergence for the two new circle-packing solvers |
 | `CGALPhase8bLite` | `test_cgal_phase8b_lite.cpp` | 17 | Phase 8b-Lite — CGAL entries for all 5 DCE models + `output_uv_map` (Euclidean, Spherical, HyperIdeal, Inversive-Distance) + CP-Euclidean throws-clearly + pipe-operator chaining |

-**Total: 240 tests, 0 skipped.**
+**Total: 246 tests, 0 skipped.**

 ---

--- a/doc/reviewer/java-port-audit.md
+++ b/doc/reviewer/java-port-audit.md
@@ -124,18 +124,26 @@ cycles, not unweighted BFS).
  requires integration of holomorphic differentials — not implemented") so it is
  consistent with the documented contract, not a regression.

-### 8. ℹ️ NOTE — Holonomy translation assumes zero residual rotation
+### 8. ✅ FIXED — Holonomy now fits the full deck isometry (was: midpoint only)
 `detail::euclidean_holonomy` (`layout.hpp`) develops each face along a dual tree
-that never crosses a cut edge, then reads each generator's translation as the
-**midpoint displacement** `ω = midA − midB` of the shared cut edge between its two
-independent developments. This equals the true deck-translation **only when the
-two developments differ by a pure translation** (linear part = identity), i.e. for
-a perfectly flat cone metric (Θ ≡ 2π). When Newton has not fully converged (small
-residual cone-angle defect), there is a residual rotation and the midpoint
-displacement is a first-order approximation rather than the exact ω. A fully
-robust version would fit the rigid motion mapping edge `(Bs,Bt)↦(At,As)` and read
-its translation part. Acceptable for converged flat tori (the `torus_8x8` test
-passes); worth tightening if higher-genus or under-converged inputs are used.
+that never crosses a cut edge.  It **previously** read each generator's
+translation as the **midpoint displacement** `ω = midA − midB` of the shared cut
+edge between its two independent developments — exact **only** when the two
+developments differ by a pure translation (linear part = identity), i.e. a
+perfectly flat cone metric (Θ ≡ 2π).  Under residual cone-angle defect that was
+a first-order approximation.
+
+**Fix (2026-05-29):** it now fits the full orientation-preserving deck isometry
+`g(z) = a·z + b` to the correspondence `(Bs,Bt) ↦ (As,At)` (`a = (At−As)/(Bt−Bs)`,
+`b = As − a·Bs`).  Because both faces use the same edge length, `|a| = 1` exactly,
+so the fit is an exact isometry; `ω = b` is its translation part and `|arg(a)|`
+is exposed per cut edge as `HolonomyData::residual_rotation` (a convergence
+diagnostic).  For a flat cone metric `a = 1` and `b` reduces to `midA − midB`, so
+converged flat tori are **bit-for-bit unchanged** — all 240 CGAL tests still pass,
+including `HolonomyEndToEnd.{Torus4x4,TorusHex6x6,Torus8x8}`.  Note: when
+`|arg(a)|` is non-negligible the holonomy is a genuine rotation and *no* lattice
+translation is well-defined; the diagnostic now surfaces that instead of silently
+averaging.

 ---

@@ -200,37 +208,75 @@ Audited against the literature and found faithful:
   (This also fixed the Euclidean holonomy extraction — develop across the dual
   tree only — see `layout.hpp` / `cut_graph.hpp`.)

-4. ⚠️ **PARTLY DONE (2026-05-29)** — added `EdgeGradient_RegularTetClosedForm`:
-   an *independent closed-form* oracle that pins each edge-DOF gradient to
-   `π/3` (= 2·(2π/3) − π) on the regular spherical tetrahedron, where the corner
-   angle 2π/3 follows from the spherical law of cosines.  This locks the
-   Finding-3 formula `α_opp⁺ + α_opp⁻ − θ_e` against a value the path-integral FD
-   check cannot detect (the energy is the integral of G, so FD only proves
-   curl-freeness).  **Still open:** (i) Newton-to-convergence *with* edge DOFs is
-   blocked by the Finding-4 spherical-Hessian guard (`throw` on edge DOFs), so a
-   converged-metric test needs an FD-Hessian or guard relaxation first; (ii) a
-   live-Java golden-value oracle still requires running the upstream library.
+4. ✅ **DONE (2026-05-29)** — the Finding-3 spherical edge-DOF gradient is now
+   locked at two levels: (a) `EdgeGradient_RegularTetClosedForm` — an independent
+   *closed-form* oracle pinning each edge-DOF gradient to `π/3` (= 2·(2π/3) − π)
+   on the regular spherical tetrahedron; and (b)
+   `SphericalGoldenJava.FullMeshEdgeDofGradient_Tetrahedron` — a *live-Java*
+   full-mesh oracle (item 5) on a generic asymmetric configuration that pins
+   BOTH the vertex gradient AND the edge gradient `α_opp⁺ + α_opp⁻ − θ_e` to
+   1e-12 vs the real raw `conformalEnergyAndGradient`.  The vertex-DOF spherical
+   path is likewise locked by `SphericalGoldenJava.FullMeshGradientAndEnergy_…`.
+   **Still open (solver feature, not a correctness gap):** Newton-to-convergence
+   *with* edge DOFs is blocked by the Finding-4 spherical-Hessian guard (`throw`
+   on edge DOFs), so a converged-metric test needs an FD-Hessian or guard
+   relaxation first.

-5. **Cross-check vs Java numeric oracles** — there is no test that pins C++
-   functional/gradient values against recorded Java outputs on a fixed mesh. A
-   handful of golden-value tests (energy + gradient at a known `x`) would catch
-   any future silent divergence from the reference far more directly than the
-   self-consistent FD checks (which only verify curl-freeness, since the C++
-   energy is itself the path-integral of its own gradient).
+5. ✅ **DONE (2026-05-29)** — golden-value oracles now span both levels:
+   • *pure-math / per-triangle* (1e-12 vs compiled upstream, openjdk 17):
+   `HyperIdealGoldenJava` (Clausen/Л/ImLi₂, ζ₁₃/₁₄/₁₅/ζ, both tetrahedron-volume
+   formulas), `EuclideanGoldenJava` (angle formula + 2·Л energy),
+   `SphericalGoldenJava` (law-of-cosines angles + β relations + Л energy),
+   `PeriodMatrix.NormalizeModulus_GoldenJava` (item 8).
+   • *full-mesh energy + gradient at a known `x`* — the originally-fragile part,
+   now solved via a shared OBJ tetrahedron loaded identically in both languages:
+   `EuclideanGoldenJava.FullMeshGradientAndEnergy_Tetrahedron` and
+   `SphericalGoldenJava.FullMeshGradientAndEnergy_Tetrahedron` drive the REAL
+   `EuclideanCyclicFunctional` / `SphericalFunctional` and pin every per-vertex
+   gradient `Θ−Σα` (matched by vertex position) AND `ΔE = E(x)−E(0)` to 1e-12.
+   The ΔE check is doubly independent: C++ computes energy as a Gauss-Legendre
+   PATH INTEGRAL of its gradient, Java as a CLOSED-FORM sum — the initialEnergy
+   constant and the φ·λ⁰ term cancel in the difference (no edge DOFs here).
+   • *full-mesh EDGE-DOF gradient* —
+   `SphericalGoldenJava.FullMeshEdgeDofGradient_Tetrahedron` makes two opposite
+   tetrahedron edges variable (replacement parameterization λ_e = x[e_idx]) and
+   pins BOTH the vertex gradient AND the edge gradient `α_opp⁺+α_opp⁻−θ_e` to
+   1e-12 vs the real raw `conformalEnergyAndGradient`, closing the Finding-3
+   solution-level gap (item 4).  It is a pure gradient oracle (no ΔE): with an
+   edge DOF, x=0 ⇒ λ_e=0 ⇒ arc length π is degenerate, so the path-integral-from-
+   origin energy reference is ill-defined; the fixed-x gradient is the
+   unambiguous quantity.  It never touches the Hessian, so the Finding-4 edge-DOF
+   Hessian guard does not block it.  Harness: `/tmp/oracle/SphereEdgeOracle.java`.
+   **Subtlety surfaced:** the spherical oracle must call Java's raw
+   `conformalEnergyAndGradient`, not `evaluate()`; `evaluate()` first runs a 1-D
+   Brent maximization over the global-scale gauge (`maximizeInNegativeDirection`),
+   which the C++ `spherical_gradient` deliberately omits (C++ factors that gauge
+   into the Newton solver's `spherical_gauge_shift`).  Harnesses:
+   `/tmp/oracle/{tet.obj,EucMeshOracle.java,SphereMeshOracle.java,SphereEdgeOracle.java}`.
+   The *edge-DOF* full-mesh gradient oracle is now also done (item 4,
+   `SphereEdgeOracle.java`).

 6. **CP-Euclidean `clausen` vs `clausen2`** — a single test asserting
   `|clausen(x) − clausen2(x)| < 1e-12` across a sweep would document Finding 5
   and guard against an approximation regression.

-7. **Period-matrix domain convention (Finding 6)** — production now uses
-   `normalizeModulus` (folds `Re ≥ 0`). Add a torus whose true τ has `Re(τ) < 0`
-   before reduction and assert the resulting `Re(τ) ≥ 0` convention, so any future
-   switch back to `reduce_to_fundamental_domain` is caught. Pair with a golden τ
-   recorded from the Java `DiscreteEllipticUtility` on the same mesh.
+7. ⚠️ **PARTLY DONE (2026-05-29)** — the `Re ≥ 0` fold convention is now locked at
+   the function level by item 8 (`NormalizeModulus_GoldenJava`, including `Re < 0`
+   inputs).  **Still open:** an *end-to-end* torus whose true τ has `Re(τ) < 0`
+   before reduction, asserting the production pipeline yields `Re(τ) ≥ 0` (so any
+   future switch back to `reduce_to_fundamental_domain` is caught at the mesh
+   level, not just the function level).

-8. **`normalizeModulus` vs Java oracle (Finding 6)** — the faithful
-   `normalizeModulus` is only checked against ad-hoc values. Add a few recorded
-   Java `normalizeModulus` outputs as golden values to lock the `Re ≥ 0` fold.
+8. ✅ **DONE (2026-05-29)** — `PeriodMatrix.NormalizeModulus_GoldenJava` pins the
+   C++ `normalizeModulus` against five recorded outputs of the real Java
+   `DiscreteEllipticUtility.normalizeModulus` (1e-12), including `Re < 0` mirror
+   folds, large-`Re` `T`-steps, and `|τ| < 1` `S`-inversions — locking the
+   `0 ≤ Re ≤ ½` fold convention.  (Generated via `/tmp/oracle/TauOracle.java`.)
+   This also surfaced a real interaction: `compute_period_matrix` can land
+   exactly on `Re = +½` (the closed right edge of `normalizeModulus`'s mirror-
+   folded domain), which the half-open SL(2,ℤ) `is_in_fundamental_domain`
+   (`−½ ≤ Re < ½`, P1-3) correctly *excludes*; `ComputePeriodMatrix_ReducedTau_InFD`
+   was updated to assert the `normalizeModulus` domain rather than the SL(2,ℤ) one.

 9. **Higher-genus cut graph (Finding 7)** — only genus-1 is exercised end-to-end.
   Add a genus-2 mesh and assert `compute_cut_graph` returns exactly `2g = 4` cut