Files
ConformalLabpp/doc/math/validation.md
Tarik Moussa 540f71a629
Some checks failed
C++ Tests / test-fast (pull_request) Successful in 2m25s
API Docs / doc-build (pull_request) Successful in 53s
C++ Tests / test-cgal (pull_request) Failing after 10m32s
release: v0.9.0 — finalise PR #11 with CHANGELOG, version bump, stub cleanup
Closes the v0.9.0 release loop on top of Phase 9a-Newton + Phase 8b-Lite:

* CHANGELOG.md (NEW) — Keep-A-Changelog format, with v0.9.0 entry
  detailing all Phase 9a / 9b / 8b-Lite contents and the doc-audit
  corrections that landed via PR #10.

* CITATION.cff — version 0.7.0 → 0.9.0, date 2026-05-18 → 2026-05-22.

* Stale HDS-port stubs removed (13 GTEST_SKIPs total):
  - code/tests/test_spherical_functional.cpp
  - code/tests/test_hyper_ideal_functional.cpp
  - code/tests/test_hyper_ideal_hyperelliptic_utility.cpp
  These referenced a "HDS port (Phase 4)" that never happened —
  CoHDS was intentionally replaced by CGAL::Surface_mesh, and the
  functional tests live in code/tests/cgal/test_*_functional.cpp.

* Test-count updates everywhere:
  - Non-CGAL  36 → 23  (drop = 13 deleted stubs)
  - CGAL      176 → 227
  - Total     212 → 250  (+38 net, 0 skipped)
  Files: README.md, CLAUDE.md, CHANGELOG.md, scripts/try_it.sh,
         doc/api/tests.md, doc/contributing.md, doc/getting-started.md,
         doc/math/validation.md.

Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
2026-05-22 04:27:24 +02:00

9.0 KiB
Raw Blame History

Mathematical Validation

This document lists analytically known results and explains how to verify them against conformallab++ output. It is the primary tool for an independent mathematician to check the correctness of the implementation.


How to run the examples

cmake -S code -B build -DWITH_CGAL=ON -DCMAKE_BUILD_TYPE=Release
cmake --build build --target conformallab_cgal_tests
ctest --test-dir build -R cgal --output-on-failure

All 227 tests pass, 0 skipped (see doc/api/tests.md).


1 — GaussBonnet (topology)

Theorem. For any closed triangulated surface M,

Σᵥ (2π  Θᵥ) = 2π · χ(M)

where χ(M) = 2 2g is the Euler characteristic.

Surface g χ Σ(2π Θᵥ)
Sphere (tetrahedron, cube, …) 0 2
Torus 1 0 0
Double torus 2 2

How to check:

#include "gauss_bonnet.hpp"
auto defect = gauss_bonnet_sum(mesh, maps);   // Σ(2π  Θᵥ)
auto chi    = mesh.euler_characteristic();
EXPECT_NEAR(defect, 2.0 * M_PI * chi, 1e-10);

Covered by: cgal.GaussBonnet.* tests in test_phase6.cpp.


2 — Period matrix: fundamental domain invariants

Theorem (SL(2,)-reduction). Every lattice τ ∈ has a unique representative in the standard fundamental domain

F = { τ ∈  :  |τ| ≥ 1,  |Re(τ)| ≤ 1/2,  Im(τ) > 0 }

After calling compute_period_matrix(hol), the returned τ must satisfy:

Condition Invariant
pd.tau_reduced.imag() > 0 τ lies in the upper half-plane
std::abs(pd.tau_reduced) >= 1.0 - 1e-10 τ outside unit disk
std::abs(pd.tau_reduced.real()) <= 0.5 + 1e-10 τ in vertical strip

These three conditions hold for any closed genus-1 triangulated surface processed through Euclidean uniformization — they are topology, not geometry.

Covered by: cgal.PeriodMatrix.TauInFundamentalDomain_* tests in test_phase7.cpp.


3 — Square-symmetric torus

Setup. Take a torus mesh with 4-fold rotational symmetry around the z-axis (e.g. code/data/off/torus_4x4.off, which has M=4 columns of vertices).

Expected. The symmetry group Z₄ acts conformally. Conformal automorphisms of the torus correspond to SL(2,) symmetries of τ. The unique fixed point of a rotation of order 4 in the modular group is τ = i. Therefore:

For a mesh with exact 4-fold symmetry and uniform edge lengths:
    Re(τ) = 0  (to machine precision, by symmetry)
    Im(τ) ≈ 1  (approaches 1 as mesh is refined)

The coarse 4×4 mesh (torus_4x4.off) gives Im(τ) in (0.7, 1.3) depending on the 3D embedding (R=2, r=1 torus of revolution has unequal inner/outer edge lengths). The uniformization algorithm finds the conformal class of the abstract metric encoded in the edge lengths.

Manual verification (run from the build directory after adding a small program or reading from the test output):

ConformalMesh mesh; load_mesh(mesh, "code/data/off/torus_4x4.off");
EuclideanMaps maps = setup_euclidean_maps(mesh);
compute_euclidean_lambda0_from_mesh(mesh, maps);
enforce_gauss_bonnet(mesh, maps);
auto res = newton_euclidean(mesh, maps);
CutGraph cg = compute_cut_graph(mesh);
HolonomyData hol;
euclidean_layout(mesh, res.x, maps, &cg, &hol, true);
PeriodData pd = compute_period_matrix(hol);
// pd.tau_reduced satisfies the fundamental domain invariants above

4 — Hexagonal-symmetric torus

Setup. Take a torus mesh with 6-fold rotational symmetry (code/data/off/torus_hex_6x6.off, M=6).

Expected. The unique τ fixed under a rotation of order 6 in SL(2,) is τ = e^{iπ/3} = ½ + i√3/2. So:

Re(τ) = 0.5  (to machine precision, by symmetry)
Im(τ) = √3/2 ≈ 0.8660

The coarse 6×6 torus of revolution approximates this: Re(τ) ≈ 0.5 by symmetry, Im(τ) approaches √3/2 as the mesh is refined toward a flat hexagonal lattice.


5 — Newton convergence rate

Theorem. Because the Euclidean and hyper-ideal energies are strictly convex (after gauge-fixing), Newton's method converges quadratically near the optimum.

Expected: for any mesh with up to a few hundred faces, Newton converges in fewer than 30 iterations starting from u = 0.

auto res = newton_euclidean(mesh, maps);
EXPECT_LT(res.iterations, 30);
EXPECT_LT(res.gradient_norm, 1e-10);

Covered by: cgal.EuclideanPipeline.ConvRates_* and similar tests.


6 — Gradient check (finite differences)

For each functional F(u), the gradient G = ∂F/∂u is verified by:

|G(u)ᵢ    (F(u + εeᵢ)  F(u  εeᵢ)) / (2ε)| < 1e-6

with ε = 1e-5. This check is run inside the test suite for all three geometries (Euclidean, Spherical, HyperIdeal) at u = 0 and at random u.

Relevant test suites:

cgal.EuclideanFunctional.GradientCheck_*
cgal.SphericalFunctional.GradientCheck_*
cgal.HyperIdealFunctional.GradientCheck_*

A failing gradient check means the energy and its derivative are inconsistent — the Newton solver will converge to the wrong point.


7 — Holonomy composition (Möbius maps)

For a closed surface, the composition of holonomies around any contractible cycle must be the identity. In genus 1 with a single handle:

T₁ · T₂ · T₁⁻¹ · T₂⁻¹ = Id    (commutator = Id for a torus)

because π₁(T²) = × is abelian.

For genus g ≥ 2, the fundamental group is non-abelian and this check does not hold, but the representation ρ: π₁(Σ_g) → SU(1,1) must still satisfy the relation

[T₁, T₂] · [T₃, T₄] · … = Id    (product of g commutators = Id)

These are the holonomy consistency checks implemented in test_phase7.cpp (cgal.HolonomyData.*).


9 — Cross-validation with geometry-central (optional / hypothetical)

Note: This section describes a possible external cross-validation that is not a prerequisite for the correctness of the implementation. It is of interest because geometry-central implements the same mathematical core (Gillespie, Springborn, Crane — SIGGRAPH 2021, building on Springborn 2020), but with a different algorithmic strategy (Ptolemaic flips + intrinsic triangulations instead of Newton on the original triangulation).

Which outputs are comparable?

Output conformallab++ geometry-central Comparable?
u-vector (scale parameters) res.x u after Yamabe flow ✓ after normalisation
UV coordinates layout.uv[v] conformal parameterisation ✓ up to Möbius transformation
Gauss-Bonnet deficit gauss_bonnet_sum() implicit via curvature flow ✓ (analytically identical)
Number of Newton iterations res.iterations Yamabe steps ~ (different algorithm)
Period matrix τ pd.tau_reduced not available
Möbius holonomy hol.T_a, T_b not available

Normalisation alignment

The u-vector in conformallab++ has one degree of freedom (global additive constant — gauge freedom after pin-fixing). geometry-central may use a different convention. Normalise before comparing:

// conformallab++: centre u
double mean_u = std::accumulate(x.begin(), x.end(), 0.0) / x.size();
std::vector<double> x_norm(x.size());
for (int i = 0; i < x.size(); ++i) x_norm[i] = x[i] - mean_u;

// Then compare with the geometry-central u-vector (also centred):
// max|x_norm[i] - gc_u[i]| < 1e-8  →  identical convergence point

When is the comparison useful?

Point in time What is possible
Now (Phase 7) Manual comparison using the same .off/.obj test meshes
After Phase 8 Automated comparison script (Python or separate C++ binary)
Phase 10 (research) Algorithm comparison: Newton vs. Ptolemaic flips on difficult meshes

Connection to the literature

The Springborn 2020 paper ("Ideal Hyperbolic Polyhedra and Discrete Uniformization") is already implemented in conformallab++ — it is the mathematical foundation for the HyperIdeal geometry mode (Phase 2/3). The geometry-central implementation is based on the extension by Gillespie, Springborn & Crane (2021), which uses the same variational principle of BobenkoSpringborn 2004 but additionally applies Ptolemaic flips to improve the triangulation during optimisation — an idea not yet implemented in conformallab++ (→ GC-2 in the phase roadmap).


8 — Checklist for an independent reviewer

Run these in order to validate the implementation:

  • ctest --test-dir build -R cgal --output-on-failure → 227 tests pass, 0 skipped
  • cgal.GaussBonnet.* all pass → topology is correctly read from mesh
  • cgal.EuclideanFunctional.GradientCheck_* pass → energy = integral of gradient
  • cgal.PeriodMatrix.TauInFundamentalDomain_* pass → SL(2,) reduction correct
  • cgal.MobiusMap.Compose_* and Inverse_* pass → Möbius arithmetic correct
  • cgal.HolonomyData.* pass → holonomy loops close up

All of the above are deterministic, analytic tests — no mesh loading, no file I/O, no floating-point non-determinism beyond standard IEEE-754.