diff --git a/code/tests/cgal/test_euclidean_functional.cpp b/code/tests/cgal/test_euclidean_functional.cpp index b66662f..f1f4063 100644 --- a/code/tests/cgal/test_euclidean_functional.cpp +++ b/code/tests/cgal/test_euclidean_functional.cpp @@ -125,7 +125,16 @@ TEST(EuclideanFunctional, AngleSumEqualsPi) } // ════════════════════════════════════════════════════════════════════════════ -// Degenerate triangle → valid = false +// Degenerate triangle — limiting angles (Finding-F, java-port-audit item 1) +// +// The BPS-energy convex C¹ extension assigns the *limiting* angles when the +// triangle inequality is violated: the corner OPPOSITE the over-long edge +// gets π, the other two get 0. These tests lock that behaviour in for +// both euclidean_angles_from_lengths() and the gradient accumulation. +// +// Before the java-port-audit Finding 1 fix, the degenerate-face code +// returned {0,0,0} and the gradient skipped the face entirely, producing +// the wrong gradient on near-flip configurations. // ════════════════════════════════════════════════════════════════════════════ TEST(EuclideanFunctional, DegenerateTriangleReturnsFalse) @@ -135,6 +144,88 @@ TEST(EuclideanFunctional, DegenerateTriangleReturnsFalse) EXPECT_FALSE(fa.valid); } +TEST(EuclideanFunctional, DegenerateTriangle_LimitingAngles_L23TooLong) +{ + // l23 = 10 >> l12 + l31 = 2 → α₁ = π (at v1, opposite l23), α₂=α₃=0. + auto fa = euclidean_angles_from_lengths(1.0, 10.0, 1.0); + EXPECT_FALSE(fa.valid); + EXPECT_NEAR(fa.alpha1, PI, 1e-12) << "corner opposite over-long l23 must be π"; + EXPECT_NEAR(fa.alpha2, 0.0, 1e-12); + EXPECT_NEAR(fa.alpha3, 0.0, 1e-12); +} + +TEST(EuclideanFunctional, DegenerateTriangle_LimitingAngles_L31TooLong) +{ + // l31 too long → α₂ = π (at v2, opposite l31). + auto fa = euclidean_angles_from_lengths(1.0, 1.0, 10.0); + EXPECT_FALSE(fa.valid); + EXPECT_NEAR(fa.alpha2, PI, 1e-12) << "corner opposite over-long l31 must be π"; + EXPECT_NEAR(fa.alpha1, 0.0, 1e-12); + EXPECT_NEAR(fa.alpha3, 0.0, 1e-12); +} + +TEST(EuclideanFunctional, DegenerateTriangle_LimitingAngles_L12TooLong) +{ + // l12 too long → α₃ = π (at v3, opposite l12). + auto fa = euclidean_angles_from_lengths(10.0, 1.0, 1.0); + EXPECT_FALSE(fa.valid); + EXPECT_NEAR(fa.alpha3, PI, 1e-12) << "corner opposite over-long l12 must be π"; + EXPECT_NEAR(fa.alpha1, 0.0, 1e-12); + EXPECT_NEAR(fa.alpha2, 0.0, 1e-12); +} + +TEST(EuclideanFunctional, DegenerateTriangle_GradientPicksUpPiCorner) +{ + // Build a single triangle. Force a degenerate effective-length by + // setting a large negative lambda0 on two edges so l12 >> l23 + l31. + // + // DOFs: all vertices free. x = 0 (no conformal scaling). + // lambda0: e_opp_v3 (i.e. l12) is huge; the other two are near-zero. + // Expected: the gradient at v3 (opposite l12) picks up −π from the + // degenerate face; G_v3 = Θ_v3 − π = 2π − π = π. + auto mesh = make_triangle(); + auto maps = setup_euclidean_maps(mesh); + + // Identify edges: h0=halfedge(face), source(h0)=v1, source(next(h0))=v2, etc. + auto f = *mesh.faces().begin(); + auto h0 = mesh.halfedge(f); + auto h1 = mesh.next(h0); + auto h2 = mesh.next(h1); + + // Assign large lambda0 to the edge opposite v3 (= edge of h0, i.e. e12). + Edge_index e12 = mesh.edge(h0); + Edge_index e23 = mesh.edge(h1); + Edge_index e31 = mesh.edge(h2); + maps.lambda0[e12] = 20.0; // l12 = exp(10) ≈ 22026 — hugely over-long + maps.lambda0[e23] = 0.0; + maps.lambda0[e31] = 0.0; + + int n = assign_euclidean_vertex_dof_indices(mesh, maps); + std::vector x(static_cast(n), 0.0); + + auto G = euclidean_gradient(mesh, x, maps); + + // v3 = source(h2). Its gradient component should include the π corner. + Vertex_index v3 = mesh.source(h2); + int iv3 = maps.v_idx[v3]; + ASSERT_GE(iv3, 0); + + // G_v3 = Θ_v3 − α3. The degenerate face gives α3 = π. + // Θ_v3 defaults to 2π, so G_v3 = 2π − π = π. + EXPECT_NEAR(G[static_cast(iv3)], PI, 1e-10) + << "Gradient at v3 must include the π limiting angle from the degenerate face"; + + // v1 and v2 get α = 0 from the degenerate face → G_vi = Θ − 0 = 2π. + Vertex_index v1 = mesh.source(h0); + Vertex_index v2 = mesh.source(h1); + int iv1 = maps.v_idx[v1], iv2 = maps.v_idx[v2]; + ASSERT_GE(iv1, 0); ASSERT_GE(iv2, 0); + EXPECT_NEAR(G[static_cast(iv1)], TWO_PI, 1e-10) + << "Gradient at v1 must be 2π (angle contribution = 0)"; + EXPECT_NEAR(G[static_cast(iv2)], TWO_PI, 1e-10) + << "Gradient at v2 must be 2π (angle contribution = 0)"; +} + // ════════════════════════════════════════════════════════════════════════════ // Gradient check: default right-isosceles triangle, vertex DOFs only // diff --git a/code/tests/cgal/test_euclidean_hessian.cpp b/code/tests/cgal/test_euclidean_hessian.cpp index 94829b7..80305c6 100644 --- a/code/tests/cgal/test_euclidean_hessian.cpp +++ b/code/tests/cgal/test_euclidean_hessian.cpp @@ -214,3 +214,42 @@ TEST(EuclideanHessian, FDCheck_MixedPinnedVertices) EXPECT_TRUE(hessian_check_euclidean(mesh, x, maps)) << "FD Hessian check failed for mixed pinned/variable vertices"; } + +// ════════════════════════════════════════════════════════════════════════════ +// Edge-DOF guard (Finding-G, java-port-audit item 2) +// +// euclidean_hessian() (vertex-only cotangent Laplacian) must throw +// std::logic_error when any edge DOF is active. Without this guard the +// function would silently return a Hessian with zero rows/cols for the +// edge DOFs, causing SimplicialLDLT to fail in a hard-to-diagnose way. +// ════════════════════════════════════════════════════════════════════════════ + +TEST(EuclideanHessian, EdgeDOFGuard_Throws) +{ + auto mesh = make_tetrahedron(); + auto maps = setup_euclidean_maps(mesh); + compute_euclidean_lambda0_from_mesh(mesh, maps); + assign_euclidean_all_dof_indices(mesh, maps); // assigns vertex + edge DOFs + + const int n = euclidean_dimension(mesh, maps); + std::vector x(static_cast(n), 0.0); + + EXPECT_THROW(euclidean_hessian(mesh, x, maps), std::logic_error) + << "euclidean_hessian must throw when edge DOFs are present"; +} + +TEST(EuclideanHessian, EdgeDOFGuard_VertexOnlyDoesNotThrow) +{ + // Vertex-only layout must NOT trigger the guard. + auto mesh = make_tetrahedron(); + auto maps = setup_euclidean_maps(mesh); + compute_euclidean_lambda0_from_mesh(mesh, maps); + auto gauge = *mesh.vertices().begin(); + assign_euclidean_vertex_dof_indices(mesh, maps, gauge); + + const int n = euclidean_dimension(mesh, maps); + std::vector x(static_cast(n), 0.0); + + EXPECT_NO_THROW(euclidean_hessian(mesh, x, maps)) + << "euclidean_hessian must not throw for vertex-only DOF layout"; +} diff --git a/code/tests/cgal/test_spherical_functional.cpp b/code/tests/cgal/test_spherical_functional.cpp index a8d745f..e82a6f4 100644 --- a/code/tests/cgal/test_spherical_functional.cpp +++ b/code/tests/cgal/test_spherical_functional.cpp @@ -652,3 +652,52 @@ TEST(SphericalGoldenJava, FullMeshEdgeDofGradient_Tetrahedron) EXPECT_NEAR(G[static_cast(maps.e_idx[eAB])], -0.35189517043413690, 1e-12); EXPECT_NEAR(G[static_cast(maps.e_idx[eCD])], -0.44101986058895950, 1e-12); } + +// ════════════════════════════════════════════════════════════════════════════ +// Degenerate spherical triangle — limiting angles (Finding-F, java-port-audit item 1) +// +// spherical_angles() must return the limiting angles (π opposite the +// over-long edge, 0/0 elsewhere) when the spherical triangle inequality +// is violated, with valid = false. This mirrors the Euclidean behaviour +// and is required for the convex C¹ BPS extension. +// ════════════════════════════════════════════════════════════════════════════ + +TEST(SphericalFunctional, DegenerateTriangle_LimitingAngles_S12TooLong) +{ + // s12 > s23 + s31: s12 = 2.5, s23 = s31 = 0.5 (all < π so valid arc lengths) + // s23 < 0 → actually use s-based check + // Easier: use s12 = π − ε (nearly degenerate hemisphere edge) + // and very short s23, s31 so s12 > s23 + s31. + const double s12 = 2.0, s23 = 0.4, s31 = 0.4; // s12 > s23+s31 = 0.8 + auto fa = spherical_angles(s12, s23, s31); + EXPECT_FALSE(fa.valid); + // s12 is the edge opposite v3 → α3 = π + EXPECT_NEAR(fa.alpha3, PI, 1e-12) + << "corner opposite over-long s12 must be π"; + EXPECT_NEAR(fa.alpha1, 0.0, 1e-12); + EXPECT_NEAR(fa.alpha2, 0.0, 1e-12); +} + +TEST(SphericalFunctional, DegenerateTriangle_LimitingAngles_S23TooLong) +{ + // s23 > s12 + s31 → α1 = π + const double s12 = 0.4, s23 = 2.0, s31 = 0.4; + auto fa = spherical_angles(s12, s23, s31); + EXPECT_FALSE(fa.valid); + EXPECT_NEAR(fa.alpha1, PI, 1e-12) + << "corner opposite over-long s23 must be π"; + EXPECT_NEAR(fa.alpha2, 0.0, 1e-12); + EXPECT_NEAR(fa.alpha3, 0.0, 1e-12); +} + +TEST(SphericalFunctional, DegenerateTriangle_LimitingAngles_S31TooLong) +{ + // s31 > s12 + s23 → α2 = π + const double s12 = 0.4, s23 = 0.4, s31 = 2.0; + auto fa = spherical_angles(s12, s23, s31); + EXPECT_FALSE(fa.valid); + EXPECT_NEAR(fa.alpha2, PI, 1e-12) + << "corner opposite over-long s31 must be π"; + EXPECT_NEAR(fa.alpha1, 0.0, 1e-12); + EXPECT_NEAR(fa.alpha3, 0.0, 1e-12); +} diff --git a/doc/reviewer/external-audit-2026-05-30.md b/doc/reviewer/external-audit-2026-05-30.md index e42688b..ba865c0 100644 --- a/doc/reviewer/external-audit-2026-05-30.md +++ b/doc/reviewer/external-audit-2026-05-30.md @@ -704,8 +704,8 @@ index. Add a comment: | C | `euclidean_hessian.hpp` | 26–27, 57–58 | Doc error | Medium | ✅ Fixed 2026-05-31 | | D | `euclidean_functional.hpp` + 2 others | 97–107 | Doc error | Medium | ✅ Fixed 2026-05-31 | | E | `cp_euclidean_functional.hpp` | 338–349, 373–387 | Inconsistency | Medium | ✅ Fixed 2026-05-31 | -| F | test files | — | Test gap | Medium | 🟠 Open | -| G | test files | — | Test gap | Medium | 🟠 Open | +| F | test files | — | Test gap | Medium | ✅ Fixed 2026-05-31 | +| G | test files | — | Test gap | Medium | ✅ Fixed 2026-05-31 | | H | test files | — | Test gap | Medium | 🟠 Open | | I | CI workflow | — | Arch risk | High | 🔵 Open | | MINOR-1 | `spherical_functional.hpp` | 404 | Doc error | Minor | 🟡 Open |