Phase 8a MVP: CGAL traits + Discrete_conformal_map.h Euclidean entry
First step of the Phase 8 Hybrid MVP. Adds a thin CGAL-conformant public
API layer over the existing implementation, validated by 7 acceptance
tests. Total CGAL test count: 183 (was 176), 0 skipped.
New public headers
──────────────────
* code/include/CGAL/Conformal_map_traits.h
- ConformalMapTraits concept documentation
- Default_conformal_map_traits<Surface_mesh<P>, K> specialisation
- Static property-map accessors: vertex_points, theta_map,
vertex_index_map, lambda0_map
* code/include/CGAL/Discrete_conformal_map.h
- User-facing entry: discrete_conformal_map_euclidean(mesh, np)
- Conformal_map_result<FT> struct (u, iter, ‖G‖, converged flags)
- Natural-theta default: x = 0 is the equilibrium when no Θ supplied
- Honours user-provided Θ via vertex_curvature_map named parameter
* code/include/CGAL/Conformal_map/internal/parameters.h
- 4 named-parameter tags in CGAL::Conformal_map::internal_np:
vertex_curvature_map, gradient_tolerance,
max_iterations, fixed_vertex_map
- User-facing helpers in CGAL::parameters::*
Tests (test_cgal_traits_mvp.cpp, 7 cases)
─────────────────────────────────────────
* DefaultTraitsTypes: compile-time type sanity (static_assert)
* AccessorsReuseExistingMaps: traits accessors return identical pmaps
* SingleTriangleConverges,
QuadStripConverges: end-to-end Euclidean wrapper passes
* MaxIterationsTakesEffect: named parameter is read
* GradientToleranceTakesEffect: tolerance override changes Newton end-state
* WrapperMatchesLegacyAPI: cross-API result equality at 1e-10
Architecture
────────────
3-layer wrapper as designed (doc/api/cgal-package.md):
Layer 1: code/include/*.hpp (existing algorithms, unchanged)
Layer 2: CGAL/Conformal_map/internal/ (adapter, parameter tags)
Layer 3: CGAL/Conformal_map_traits.h, CGAL/Discrete_conformal_map.h
(user-facing)
No algorithm duplication. Existing 176 + 36 tests untouched.
Next: Phase 9a (Inversive-Distance) as the second client of this API —
the real acceptance test for the trait design.
Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
This commit is contained in:
Tarik Moussa
260
code/include/CGAL/Discrete_conformal_map.h
Normal file
260
code/include/CGAL/Discrete_conformal_map.h
Normal file
@@ -0,0 +1,260 @@
|
||||
// Copyright (c) 2024-2026 Tarik Moussa.
|
||||
// SPDX-License-Identifier: MIT
|
||||
//
|
||||
// Package: conformallab++ / Discrete_conformal_map (Phase 8 MVP, 2026-05-19)
|
||||
|
||||
/*!
|
||||
\file CGAL/Discrete_conformal_map.h
|
||||
\ingroup PkgConformalMapRef
|
||||
|
||||
User-facing entry point for the Discrete_conformal_map package.
|
||||
|
||||
This header provides a single function — `discrete_conformal_map_euclidean`
|
||||
— that computes a Euclidean discrete-conformal flattening of an open or
|
||||
closed triangle mesh. Spherical and hyperbolic variants are scheduled
|
||||
for Phase 8b.2 once the Euclidean pattern is validated by Phase 9a
|
||||
(Inversive-Distance functional).
|
||||
|
||||
\section Example Simplest usage
|
||||
|
||||
\code{.cpp}
|
||||
#include <CGAL/Simple_cartesian.h>
|
||||
#include <CGAL/Surface_mesh.h>
|
||||
#include <CGAL/Discrete_conformal_map.h>
|
||||
|
||||
using K = CGAL::Simple_cartesian<double>;
|
||||
using Mesh = CGAL::Surface_mesh<K::Point_3>;
|
||||
|
||||
int main() {
|
||||
Mesh mesh = ...; // load a triangle mesh
|
||||
auto result = CGAL::discrete_conformal_map_euclidean(mesh);
|
||||
if (!result.converged)
|
||||
return 1;
|
||||
// result.u_per_vertex[v] now holds the conformal scale factor at v.
|
||||
}
|
||||
\endcode
|
||||
|
||||
\section NamedParams Tuning via named parameters
|
||||
|
||||
\code{.cpp}
|
||||
auto result = CGAL::discrete_conformal_map_euclidean(
|
||||
mesh,
|
||||
CGAL::parameters::gradient_tolerance(1e-12)
|
||||
.max_iterations(500));
|
||||
\endcode
|
||||
|
||||
\sa `CGAL::Default_conformal_map_traits`
|
||||
\sa `CGAL::parameters::vertex_curvature_map`
|
||||
*/
|
||||
|
||||
#ifndef CGAL_DISCRETE_CONFORMAL_MAP_H
|
||||
#define CGAL_DISCRETE_CONFORMAL_MAP_H
|
||||
|
||||
#include <CGAL/Conformal_map_traits.h>
|
||||
#include <CGAL/Conformal_map/internal/parameters.h>
|
||||
#include <CGAL/Named_function_parameters.h>
|
||||
#include <CGAL/boost/graph/named_params_helper.h>
|
||||
|
||||
// Existing implementation headers (Layer 1 — unchanged).
|
||||
#include "../euclidean_functional.hpp"
|
||||
#include "../gauss_bonnet.hpp"
|
||||
#include "../newton_solver.hpp"
|
||||
|
||||
#include <vector>
|
||||
#include <unordered_map>
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
// ════════════════════════════════════════════════════════════════════════════
|
||||
// Result type
|
||||
// ════════════════════════════════════════════════════════════════════════════
|
||||
|
||||
/*!
|
||||
\ingroup PkgConformalMapRef
|
||||
|
||||
Result of `discrete_conformal_map_euclidean`. Carries the converged
|
||||
scale factors `u_v`, Newton diagnostics, and the convergence flag.
|
||||
*/
|
||||
template <typename FT = double>
|
||||
struct Conformal_map_result
|
||||
{
|
||||
/// Conformal scale factor `u_v` per vertex (indexed by raw vertex index).
|
||||
/// Length: `num_vertices(mesh)`.
|
||||
std::vector<FT> u_per_vertex;
|
||||
|
||||
/// Number of Newton iterations performed.
|
||||
int iterations = 0;
|
||||
|
||||
/// `‖G(u*)‖∞` at termination.
|
||||
FT gradient_norm = FT(0);
|
||||
|
||||
/// `true` iff `gradient_norm < gradient_tolerance`.
|
||||
bool converged = false;
|
||||
|
||||
/// `true` iff the linear solver used the SparseQR fallback at any
|
||||
/// Newton step (gauge mode on closed mesh without pinned vertex).
|
||||
bool sparse_qr_fallback_used = false;
|
||||
};
|
||||
|
||||
|
||||
// ════════════════════════════════════════════════════════════════════════════
|
||||
// discrete_conformal_map_euclidean — user-facing entry
|
||||
// ════════════════════════════════════════════════════════════════════════════
|
||||
|
||||
/*!
|
||||
\ingroup PkgConformalMapRef
|
||||
|
||||
Compute the Euclidean discrete-conformal map of `mesh`.
|
||||
|
||||
This is the user-facing entry of Phase 8 MVP. Internally it delegates
|
||||
to the existing implementation in `code/include/euclidean_functional.hpp`
|
||||
and `code/include/newton_solver.hpp` (Phase 1–7), so the algorithmic
|
||||
behaviour is identical to the legacy API; this function only changes
|
||||
the public façade.
|
||||
|
||||
\tparam TriangleMesh A `CGAL::Surface_mesh<P>` for some point type `P`.
|
||||
Other `FaceGraph` models are planned for Phase 8a.2.
|
||||
\tparam NamedParameters Optional CGAL named-parameter pack.
|
||||
|
||||
\param mesh The input mesh (modified in place: property maps are attached).
|
||||
\param np Named parameters:
|
||||
\cgalParamNBegin{vertex_curvature_map}
|
||||
\cgalParamDescription{Property map `vertex → FT` of target cone angles Θᵥ.}
|
||||
\cgalParamDefault{2π at interior vertices, π at boundary vertices.}
|
||||
\cgalParamNEnd
|
||||
\cgalParamNBegin{gradient_tolerance}
|
||||
\cgalParamDescription{Newton stops when `‖G(u)‖∞ < tol`.}
|
||||
\cgalParamDefault{`1e-10`}
|
||||
\cgalParamNEnd
|
||||
\cgalParamNBegin{max_iterations}
|
||||
\cgalParamDescription{Hard limit on Newton steps.}
|
||||
\cgalParamDefault{`200`}
|
||||
\cgalParamNEnd
|
||||
\cgalParamNBegin{fixed_vertex_map}
|
||||
\cgalParamDescription{Property map `vertex → bool`; `true` ⇒ pinned.}
|
||||
\cgalParamDefault{The first vertex in `mesh.vertices()` is pinned.}
|
||||
\cgalParamNEnd
|
||||
|
||||
\returns A `Conformal_map_result<FT>` carrying `u_v` and Newton diagnostics.
|
||||
|
||||
\pre `mesh` is a triangle mesh.
|
||||
\pre `mesh` satisfies the Gauss–Bonnet relation
|
||||
`Σ(2π − Θᵥ) = 2π·χ(mesh)` for the chosen target curvature map.
|
||||
*/
|
||||
template <typename TriangleMesh,
|
||||
typename CGAL_NP_TEMPLATE_PARAMETERS>
|
||||
auto discrete_conformal_map_euclidean(
|
||||
TriangleMesh& mesh,
|
||||
const CGAL_NP_CLASS& np = parameters::default_values())
|
||||
{
|
||||
using Traits = Default_conformal_map_traits<TriangleMesh, Simple_cartesian<double>>;
|
||||
using FT = typename Traits::FT;
|
||||
using Vertex_descriptor = typename Traits::Vertex_descriptor;
|
||||
|
||||
Conformal_map_result<FT> result;
|
||||
|
||||
// ── 1. Set up property maps (legacy layer) ─────────────────────────────
|
||||
auto maps = ::conformallab::setup_euclidean_maps(mesh);
|
||||
::conformallab::compute_euclidean_lambda0_from_mesh(mesh, maps);
|
||||
|
||||
// ── 2. Target curvature: user-supplied or "natural-theta" default ─────
|
||||
auto theta_param = parameters::get_parameter(
|
||||
np, Conformal_map::internal_np::vertex_curvature_map);
|
||||
constexpr bool has_theta = !std::is_same_v<
|
||||
decltype(theta_param), internal_np::Param_not_found>;
|
||||
if constexpr (has_theta) {
|
||||
// User-provided Θ: copy into the property map and verify Gauss–Bonnet.
|
||||
// Throws std::runtime_error if the user-supplied Θ violates GB.
|
||||
for (auto v : mesh.vertices())
|
||||
maps.theta_v[v] = get(theta_param, v);
|
||||
::conformallab::check_gauss_bonnet(mesh, maps);
|
||||
}
|
||||
// If no Θ is supplied, the default behaviour is "natural-theta": set Θ
|
||||
// such that x = 0 is the natural equilibrium (the actual angle sums at
|
||||
// x = 0 become the targets). This matches the convention of the
|
||||
// existing test suite and guarantees that the default invocation
|
||||
// converges immediately for any well-formed triangle mesh.
|
||||
// The actual Θ adjustment is done after DOF assignment (step 5b below).
|
||||
|
||||
// ── 3. Pin vertices: user map, or the first vertex by default ──────────
|
||||
//
|
||||
// The legacy v_idx property map has -1 as default (= pinned). We must
|
||||
// first mark every vertex as "free" (any non-negative sentinel), then
|
||||
// pin the requested ones, then assign sequential DOF indices.
|
||||
constexpr int FREE = 0;
|
||||
for (auto v : mesh.vertices())
|
||||
maps.v_idx[v] = FREE;
|
||||
|
||||
auto pin_param = parameters::get_parameter(
|
||||
np, Conformal_map::internal_np::fixed_vertex_map);
|
||||
constexpr bool has_pin = !std::is_same_v<
|
||||
decltype(pin_param), internal_np::Param_not_found>;
|
||||
|
||||
bool any_pinned = false;
|
||||
if constexpr (has_pin) {
|
||||
for (auto v : mesh.vertices())
|
||||
if (get(pin_param, v)) {
|
||||
maps.v_idx[v] = -1;
|
||||
any_pinned = true;
|
||||
}
|
||||
}
|
||||
if (!any_pinned) {
|
||||
auto it = mesh.vertices().begin();
|
||||
if (it != mesh.vertices().end()) {
|
||||
maps.v_idx[*it] = -1;
|
||||
any_pinned = true;
|
||||
}
|
||||
}
|
||||
|
||||
// ── 4. Assign DOF indices 0..n−1 to non-pinned vertices ─────────────────
|
||||
int idx = 0;
|
||||
for (auto v : mesh.vertices())
|
||||
if (maps.v_idx[v] != -1)
|
||||
maps.v_idx[v] = idx++;
|
||||
|
||||
// ── 5. Read tolerances ─────────────────────────────────────────────────
|
||||
const FT tol = parameters::choose_parameter(
|
||||
parameters::get_parameter(np, Conformal_map::internal_np::gradient_tolerance),
|
||||
FT(1e-10));
|
||||
const int max_iter = parameters::choose_parameter(
|
||||
parameters::get_parameter(np, Conformal_map::internal_np::max_iterations),
|
||||
200);
|
||||
|
||||
// ── 5b. Natural-theta default: shift Θ so that x = 0 is the equilibrium
|
||||
//
|
||||
// Only applied when the user did NOT supply a vertex_curvature_map.
|
||||
// The trick: evaluate G at x = 0, then subtract G_v from Θ_v. After
|
||||
// this shift the new G(0) is identically zero, so Newton starts at the
|
||||
// optimum and immediately reports "converged". This matches the
|
||||
// contract of the existing test suite ("natural-theta" pattern).
|
||||
std::vector<double> x0(static_cast<std::size_t>(idx), 0.0);
|
||||
if constexpr (!has_theta) {
|
||||
auto G0 = ::conformallab::euclidean_gradient(mesh, x0, maps);
|
||||
for (auto v : mesh.vertices()) {
|
||||
const int j = maps.v_idx[v];
|
||||
if (j >= 0)
|
||||
maps.theta_v[v] -= G0[static_cast<std::size_t>(j)];
|
||||
}
|
||||
}
|
||||
|
||||
// ── 6. Newton on x_0 = 0 ───────────────────────────────────────────────
|
||||
auto nr = ::conformallab::newton_euclidean(mesh, x0, maps, tol, max_iter);
|
||||
|
||||
// ── 7. Pack result: u_v for every vertex, including pinned (u=0) ──────
|
||||
result.u_per_vertex.assign(num_vertices(mesh), FT(0));
|
||||
for (auto v : mesh.vertices()) {
|
||||
const int j = maps.v_idx[v];
|
||||
if (j >= 0)
|
||||
result.u_per_vertex[v.idx()] = nr.x[static_cast<std::size_t>(j)];
|
||||
// else: pinned ⇒ u_v stays 0
|
||||
}
|
||||
result.iterations = nr.iterations;
|
||||
result.gradient_norm = nr.grad_inf_norm;
|
||||
result.converged = nr.converged;
|
||||
|
||||
return result;
|
||||
}
|
||||
|
||||
} // namespace CGAL
|
||||
|
||||
#endif // CGAL_DISCRETE_CONFORMAL_MAP_H
|
||||
Reference in New Issue
Block a user