Centroidal Voronoi Tessellations

• CVTs are spatial partitions where each Voronoi cell’s generator coincides with its centroid, minimizing a quadratic energy functional.
• Lloyd's algorithm iteratively computes Voronoi partitions and updates centroids, ensuring energy reduction and convergence under nondegenerate conditions.
• Extensions include weighted, anisotropic, and geodesic CVTs, broadening applications in mesh generation, image segmentation, and resource allocation.

Centroidal Voronoi Tessellations (CVTs) are spatial partitions in which each generating point (generator) of a Voronoi cell coincides with the centroid (mass center) of that cell, with respect to a given density function. CVTs minimize a quadratic "quantization" or "centroidal" energy and arise as fixed points of Lloyd's algorithm, a prototypical iterative scheme. CVTs are central in computational geometry, numerical analysis, optimal quantization, image segmentation, mesh generation, and, in generalized forms, underpin modern resource allocation, physical modeling, and biological pattern formation.

1. Mathematical Foundations and Energy Minimization

Let $\Omega\subset \mathbb{R}^d$ be a compact domain with density $\rho:\Omega\to\mathbb{R}_+$. A Voronoi tessellation is defined by a set of $N$ generators $\{z_i\}_{i=1}^N\subset\Omega$: for each $i$,

$$V_i = \{x\in \Omega : \|x-z_i\| \le \|x-z_j\| \quad \forall j\ne i\} .$$

A tessellation $\{V_i\}$ is centroidal if $z_i = c_i$, where

$$c_i = \frac{\int_{V_i} x\, \rho(x)\, dx}{\int_{V_i} \rho(x)\, dx} .$$

CVTs are (local) minimizers of the centroidal energy functional

$$E(\{z_i\},\{V_i\}) = \sum_{i=1}^N \int_{V_i} \|x - z_i\|^2 \rho(x) dx.$$

The condition for stationarity is that each $z_i$ coincides with the centroid $c_i$ of $V_i$ (Mullaghy, 26 Mar 2025). The functional extends naturally to discrete (probability) measures, Riemannian manifolds, and weighted metrics.

2. Algorithmic Approaches, Optimization, and Convergence

The most widespread computational method for constructing CVTs is Lloyd's algorithm, which iterates:

• Voronoi partition construction for current generators.
• Centroid computation of each cell.
• Generator relocation to the centroid.

Per-iteration cost is $O(N\log N)$ in 2D, higher in higher dimensions. Under generic nondegeneracy conditions on $\rho$ and $\Omega$, Lloyd’s algorithm guarantees monotonic decrease of energy and converges (in the limit) to a fixed point—a (possibly local) minimizer (Mullaghy, 26 Mar 2025, Ye et al., 2019). When combined with quasi-Newton optimization (e.g., L-BFGS), convergence is accelerated, which is especially valuable for highly non-uniform and large-scale CVTs (Yang et al., 2017).

The energy landscape is highly nonconvex with many local minima. Deterministic global optimization, such as thermodynamic annealing, can be employed, where the electrostatic interpretation of the energy provides the basis for randomized sampling and annealing strategies to find (possibly) global minimizers (Mullaghy, 3 Apr 2025). Symmetry-preserving local minima can be escaped via controlled perturbations, such as small rotations or deterministic moves that guide the search out of symmetric energy plateaus (Mullaghy, 26 Mar 2025, Gonzalez et al., 2020).

3. Extensions: Weighted, Anisotropic, Geodesic, and Constrained CVTs

Weighted CVTs (power diagrams) generalize the unweighted case by incorporating weights on generators, leading to spatial partitions defined by power distance and allowing control over cell areas and locations (Bourne et al., 2014, 0912.3974). Anisotropic or Riemannian CVTs replace the Euclidean norm in the energy by a domain-dependent metric, frequently derived from physical fields such as stress, inducing non-Euclidean tessellations that adapt cell density and shape to anisotropy (Pietroni et al., 2014). On manifolds, geodesic CVTs use intrinsic geodesic distance on the underlying geometry; this is crucial for applications like mesh generation and data segmentation on surfaces or embedded manifolds (Ye et al., 2019).

CVTs can be formulated with geometric constraints, e.g., enforcing equal area, avoiding small edges, or target cell density, by incorporating penalty terms into the objective functional and adapting gradient computations via shape calculus (Birgin et al., 25 Aug 2025). In high-dimensional spaces, CVTs can be explicitly constructed as Cartesian products of one-dimensional CVTs if the density factorizes, producing efficient, albeit axis-aligned, grid-like tessellations (Telsang et al., 2022).

4. Applications: Sampling, Remeshing, Resource Allocation, and Physical Modeling

CVTs are foundational in adaptive sampling, mesh optimization, and quantization. In Blaschke–Santaló diagram approximation, they yield highly uniform samplings in low-dimensional images of high- or infinite-dimensional parameter spaces, massively outperforming Monte Carlo sampling for a given sample budget (Bogosel et al., 2023). In surface and volumetric remeshing, CVT-based schemes (Euclidean, geodesic, or restricted) produce uniform, high-quality meshes and support adaptive density via curvature- or feature-driven strategies (Fei et al., 20 May 2025, Thomas et al., 2024). Spherical and pseudometrically-constrained CVTs enable highly uniform and symmetric point sets for applications in global illumination, MRI sampling, and climate mesh generation (Koay, 2012, Yang et al., 2017).

CVT frameworks underpin decentralized resource allocation via local negotiation algorithms enforcing global sum constraints, e.g., for smart grid power assignments (Telsang et al., 2022). Theoretical results link CVTs to electrostatic equilibria, and to the minimization of mechanical strain energy in biophysical models, such as the emergence of CVT-like patterns in epithelia under mechanical rigidity (Mullaghy, 3 Apr 2025, Lim et al., 15 Dec 2025).

5. Theoretical Structure: Gersho's Conjecture, Optimality, and Geometric Complexity

In $d$ dimensions, optimal CVT cell shapes are conjectured to be periodic (Gersho's conjecture): regular hexagons in 2D, truncated octahedra (BCC lattice) in 3D. For $n$ large enough, the minimal energy is achieved by partitions into asymptotically congruent, space-tiling shapes—convex polytopes with bounded diameter and face complexity (Choksi et al., 2018). Theoretical bounds on cell size, volume, and face number have been obtained, reducing global structure conjectures to extremely high-dimensional convex programs.

On fractal and singular supports (e.g., Cantor sets), explicit CVTs and their quantization errors can be constructed by leveraging self-similarity, enabling exact computation of optimal quantizers and exposing rich nonuniqueness and bifurcation behaviors (Roychowdhury, 2015, Dettmann et al., 2015).

6. Practical Implementation: Acceleration Strategies, Complexity, and Quality Metrics

Numerical computation of large CVTs is expedited by local neighborhood tracking, KD-tree acceleration, domain decomposition for parallelism, and multi-grid refinement. In spherical and high-resolution contexts, Lloyd-preconditioned L-BFGS schemes and incremental neighbor update policies yield an order-of-magnitude speed-up over classical Lloyd's method (Yang et al., 2017, Koay, 2012). Quality metrics such as cell aspect ratio, minimal and maximal angle, area uniformity, and geometric fidelity (e.g., Hausdorff distance) guide mesh optimization and application-specific evaluation (Fei et al., 20 May 2025). Adaptive schemes balance runtime and quality by focusing computational efforts according to local geometric complexity.

7. Open Questions and Current Research Directions

Major open problems include the full resolution of the 3D Gersho conjecture, efficient CVT generation in very high dimensions with non-factorizing densities, and robust algorithms for CVTs on general curved manifolds and singular spaces. Incorporating additional geometric or application-specific constraints into the energy (e.g., angle bounds, anisotropy, hard area constraints) remains an area of active investigation, as does the development of fast and scalable algorithms for massive instances (e.g., global climate simulations, neural 3D modeling) (Birgin et al., 25 Aug 2025, Pietroni et al., 2014, Thomas et al., 2024). Connections to physical theory (e.g., electrostatics, elasticity) and their exploitation for better optimization and understanding of emergent patterns in biology and physics continue to drive both theory and application (Mullaghy, 3 Apr 2025, Lim et al., 15 Dec 2025).