Repulsive Surface Algorithm in Ground-State Optimization

• Repulsive Surface Algorithm is a computational framework that models energy-minimizing configurations of repulsive particles and collision-avoiding surface meshes.
• It combines numerical Monte Carlo annealing with analytical helical lattice parameterization to efficiently determine ground states on cylindrical and other developable surfaces.
• The method leverages hierarchical acceleration and block-dynamics sampling to scale optimization tasks and estimate key statistical-mechanical observables like surface pressure.

The repulsive surface algorithm encompasses a set of computational and analytical methodologies for investigating ground-state configurations and constraint-aware optimization in systems governed by repulsive interactions on surfaces or within Euclidean space. Central settings include the determination of energy-minimizing phyllotactic states on cylindrical geometries, robust collision-avoiding geometric optimization for discretized surfaces, and the sampling and statistical analysis of repulsive point processes under Gibbsian measures. Frameworks in this domain target a range of mathematical physics, statistical mechanics, geometry processing, and combinatorial optimization objectives.

1. Ground-State Search for Repulsive Particles on Cylindrical Surfaces

The repulsive-surface algorithm for phyllotactic ground state search—introduced by Tomlinson & Wilkin (Tomlinson et al., 2020)—addresses the zero-temperature configuration problem for $N$ particles confined to a cylindrical surface, each subject to a repulsive pair potential. The methodology consists of a two-pronged approach, combining numerical zero-temperature Monte Carlo (MC) annealing with an analytical parameterization of minimal-energy helical lattices.

Numerical MC Annealing:

• The cylindrical surface is unwrapped into a planar rectangular domain with periodic boundary conditions in the angular ($y$) direction of width $c$ (circumference) and, optionally, the axial ($x$) direction for bulk properties. Particles are initialized randomly in an $L \times c$ box at linear density $\rho_1 = N/L$.
• Standard Metropolis moves at $T=0$ are used: one particle is chosen at random and displaced in $(x, y)$ by a small uniform offset. A move is accepted only if the total pairwise energy decreases ($\Delta E < 0$). Iteration continues until no further downhill moves are possible (energy plateau achieved).

Analytical Helical Lattice Parameterization:

• Ground states are analytically constructed as isosceles Bravais lattices parameterized by primitive vectors $a$ and $b$, constrained by: (i) $\hat{z}\cdot(b \times a) = c/\rho_1$ (cell area), (ii) $m a + n b = \hat{y} c$ (periodicity), and (iii) $|a| = |b|$ (isosceles).
• The integers $m$, $n$ serve as "parastichy numbers," with $p = m^2 - n^2$ and $\alpha = c\rho_1$. Closed-form expressions define the primitive vectors in terms of these parameters.
• Sites are indexed by $(k, l)$ and mapped back to cylindrical coordinates via $\theta_{k, l} = (2\pi / c) y_{k, l}$, $z_{k, l} = x_{k, l}$.

2. Universal Parabolic Energy Collapsing and Efficient Phase-Boundary Computation

A principal insight is the numerical collapse of the per-particle energy $E_{m,n}(\alpha)$, for a given $(m, n)$ lattice, as a universal quadratic function of the common vector length $L(\alpha; m, n) = |a| = |b|$: $$E_{m,n}(L) \approx E_{\text{min}} + \frac{1}{2} K (L - L_0)^2,$$ with constant stiffness $K$ and $L_0$ corresponding to the perfect triangular lattice. This allows substantial simplification:

• Rather than performing explicit pairwise summation or energy-minimization at each parameter value $(c, \rho_1)$, a catalog of candidate phyllotactic states (indexed by $(m, n)$) may be constructed.
• For each candidate, $L(\alpha; m, n)$ is computed, plugged into the parabola fit, and the state minimizing $E$ at each $\alpha$ located directly.

Phase transitions between states occur where the parabolic energy curves cross. The quadratic approximation yields algebraic transition conditions for $\alpha_T$, leading immediately to (critical) densities or circumference boundaries. This enables fast determination of the global ground state and phase boundaries across the parameter space without brute-force energy minimization.

3. Generalization to Other Developable Surfaces

The repulsive-surface approach is readily extended beyond cylinders:

• For any developable surface (such as cones or ribbons), the essential steps are: unwrapping to a representation with periodic boundary; constructing the primitive lattice vectors according to new geometric constraints; verifying (numerically) the universal parabola collapse; and locating transitions via curve crossings.
• The method's key components—analytical lattice construction and universal quadratic energy form—remain applicable, with adaptation of the periodicity and area constraints as dictated by the surface geometry.
• This provides a nearly exact, highly efficient algorithm for mapping sequences of phyllotactic or minimal-energy states for repulsive particles on a broad class of surfaces.

4. Repulsive Surfaces: Collision-Avoiding Optimization for Surface Meshes

The repulsive-surfaces algorithm (Yu et al., 2021) addresses the distinct, but conceptually related, problem of optimizing the geometry of surfaces in $\mathbb R^3$ while guaranteeing non-self-intersection and respecting various geometric constraints (area, volume, pins).

Key Components:

• The core is the tangent-point energy for a smooth immersion $f: M \rightarrow \mathbb{R}^3$, where $M$ is a 2-manifold. For any pair $(x, y) \in M \times M$, the kernel

$$k_{f,p}(x,y) = \frac{\|P_f(x)[f(x)-f(y)]\|^p}{\|f(x)-f(y)\|^{2p}},$$

with $P_f(x)$ the normal projector and typically $p=6$, yields a double-integral energy function that provides an infinite barrier to (self-)collision.

• The energy is discretized over triangle meshes, and its gradient is computed via midpoint quadrature and analytical differentiation.
• For optimization, a fractional Sobolev inner-product preconditioner is assembled, approximated by a combination of sparse Poisson solves and hierarchically-applied fractional operators, removing mesh-dependent timestep restrictions inherent to $L^2$-gradient flows.
• Hierarchical acceleration is implemented via a BVH (Barnes–Hut) for $O(|F| \log |F|)$ energy and gradient evaluations, and block-cluster hierarchical matrices for fast mat–vecs with the preconditioner.
• Constraints are handled through saddle-point solves and, for nonlinear cases, via corrector-projections.
• The resulting optimization loop (with remeshing, block-hierarchical evaluation, and inner iterative solves) achieves interactive rates even for $\sim 30,000$-triangle meshes.

Applications:

This framework supports minimization of genus-$g$ surfaces (knotted or unknotted), automated isotopy (e.g., torus eversion), variational geometric modeling, shrink-wrapping, nesting, collision-resolving flows, and generative modeling of surface ensembles.

5. Algorithmic Approximation of Surface Pressure in Repulsive Point Processes

In the context of repulsive Gibbs point processes with a pair potential $\phi$ (finite-range, symmetric, nonnegative), the "surface pressure" is a key statistical–mechanics observable describing finite-size corrections to pressure. The algorithmic framework in (Michelen et al., 2022) enables efficient approximation of both the partition function $Z_\Lambda(\lambda)$ and the surface pressure $sp_\Box(\lambda)$.

Main elements:

• The partition function and pressure are given by

$$p(\lambda) = \lim_{n\to\infty} \frac{1}{|\Lambda_n|} \log Z_{\Lambda_n}(\lambda),$$

and the surface pressure correction is

$$sp_\Box(\lambda) = \lim_{n\to\infty} \frac{\log Z_{\Lambda_n}(\lambda) - |\Lambda_n| p(\lambda)}{SA(\Lambda_n)}.$$

• Efficient randomized FPRAS algorithms are constructed via block-dynamics sampling. Each block step replaces the configuration inside a randomly chosen ball with a newly sampled configuration conditioned on the outside.
• For activities $\lambda < e/\Delta_\phi$ (with $\Delta_\phi$ the potential-weighted connective constant), strong spatial mixing (SSM) holds, enabling block-dynamics to mix in $O(N \log N)$ steps in a volume-$N$ box for approximate sampling, and allowing estimation of the partition function and surface pressure with $\varepsilon$-relative and $\varepsilon$-additive error in polynomial time.
• The approach utilizes the "box-identity" for surface pressure, which reduces its computation to integrals of one-point densities under Gibbs measures with modified activities; these densities are efficiently approximated by the block-dynamics procedure.

6. Computational Complexity and Implementation Considerations

Problem Domain	Algorithmic Complexity	Main Constraints/Assumptions
Phyllotactic ground states on cylinders	Quadratic in catalog size; analytical for parabolic ansatz	Zero temperature; pairwise repulsive potential; periodic surface constraints
Repulsive surfaces (meshes)	$O(|F|\log|F|)$ per iteration with hierarchical acceleration	Non-intersecting immersions; constraints via saddle-point systems; triangle meshes
Repulsive Gibbs point processes	$O(N^2 \log(N/\varepsilon) / \varepsilon^2)$ (partition), $O(\varepsilon^{-2}\log^{d+3}(1/\varepsilon))$ (surface pressure)	Finite-range, nonnegative pair potential; $\lambda < e/\Delta_\phi$ (SSM regime)

Principal implementation features include: the use of preconditioned iterative solvers and hierarchical data structures (BVH, block-cluster trees) for scalability; direct algebraic evaluation for phase boundaries in reduced models; and Markov chain Monte Carlo governed by strong spatial mixing bounds for high-dimensional sampling and partition function estimation. Notably, the combination of a few MC annealing runs with universal parabola collapse, or hierarchical energy evaluations with constraint-aware flows, consistently yields near-exact, scalable ground-state or surface-optimization algorithms.

7. Context, Significance, and Related Directions

The frameworks and algorithms described provide robust, efficient, and nearly universal tools for determining minimum-energy or collision-avoiding configurations of repulsive systems across both discrete (point processes, lattice packings) and continuum (surface immersions, mesh optimization) settings. The emergence of universal quadratic energy forms, strong spatial mixing, and combinatorial/analytic reductions underlies the computational tractability and rigor of these results. These approaches extend to a wide range of geometric and statistical mechanical models, including but not limited to phyllotaxis, surface geometry processing, spatial random processes, and variational collision-avoidance in computational geometry. All major results are grounded in the cited works (Tomlinson et al., 2020, Yu et al., 2021), and (Michelen et al., 2022).