Triangular minimality for planar Coulomb renormalized energy

1 Introduction

The Coulomb interaction in two dimensions is logarithmic. An infinite configuration of equal point charges against a uniform neutralizing background has an infinite electrostatic energy, even in a bounded region containing a charge. Subtracting the logarithmic self-energy at each charge and then taking energy per unit area gives the renormalized energy studied here. Its proposed minimizer is the triangular lattice.

1.1 The energy and the main results

For R>1, let K_R=[-R,R]^2. Fix smooth cutoffs \chi_R\in C_c^\infty(\R^2) satisfying

Let \Lambda\subset\R^2 be locally finite and simple, so every point has multiplicity one, and put \nu_\Lambda=\sum_{p\in\Lambda}\delta_p. The admissible class \mathcal A_1 consists of locally integrable vector fields E such that

where the equations hold in distributions and B_R is the disk of radius R centered at the origin. These conditions imply E\in L^q_{\mathrm{loc}}(\R^2;\R^2) for every 1\leq q<2: locally subtract the finitely many Coulomb singularities and a smooth field with the prescribed constant divergence; the remaining field has harmonic components. Section 2 supplies the details. Consequently the same class results if local L^q integrability is required for any one, or for all, of these exponents.

For E\in\mathcal A_1 define

The local puncture limit exists and is finite. The outer limsup is initially allowed to take either infinite value. The puncture limit is taken first, and the cutoff layer has the fixed width in (1.1) throughout.

Identify \R^2 with \mathbb C when convenient. The triangular lattice of covolume one is

Let h_\triangle be the mean-zero periodic distributional solution of -\Delta h_\triangle=2\pi(\delta_0-1) on \R^2/\Lambda_\triangle, and set E_\triangle=-\nabla h_\triangle.

In particular W(E)=-\infty is impossible in this class. The inequality includes fields of energy +\infty. No separation, periodicity, or Voronoi-cell bound is imposed on an arbitrary competitor.

The finite-volume statement used to prove Theorem 1.1 has its own exact normalization. For an integer n\geq2, write T_n=\R^2/(\sqrt n\,\Z^2). Given distinct points a_1,\ldots,a_n\in T_n, let h have mean zero and satisfy -\Delta h=2\pi(\sum_{i=1}^n\delta_{a_i}-1) on T_n. Define

This is a lower bound on square tori. It does not assert equality on any particular square torus. Constant harmonic fields increase the energy by n|c|^2/2; the required orthogonality follows by integrating the periodic gradient, as proved in (2.12).

The planar minimum also determines the linear term in optimal logarithmic energy for finitely many points on the sphere.

This proves the d=2 case of the logarithmic-energy conjecture of Brauchart, Hardin, and Saff [4, Conjecture 4]. Bétermin and Sandier had already proved the existence of the linear coefficient and shown that it equals (1.8) precisely when the density-one triangular lattice minimizes their planar renormalized energy [1, Theorem 1.5]. Section 8 checks the normalization and applies that equivalence to Theorem 1.1. The conclusion is an asymptotic value, not a construction of near-optimal point sets or an algorithm for Smale’s seventh problem.

1.2 Context and prior work

The finite-vortex renormalized energy of Bethuel, Brezis, and Hélein [2] is an antecedent of the infinite-system energy considered here. Sandier and Serfaty introduced the planar functional in their analysis of microscopic vortex patterns in the Ginzburg–Landau model. They proved existence of a minimum, approximation of its value by square-periodic fields, and triangular optimality among Bravais lattices [17, Theorems 1 and 2]. Their Conjecture 1 asks whether the triangular value is also the minimum over the full admissible class. Theorem 1.1 resolves this conjecture positively. The current normalization in that paper has one point per area 2\pi/m; our density-one convention corresponds to its background m=2\pi.

The lattice comparison belongs to a longer history. The minimization of convergent Epstein zeta sums at fixed covolume was developed by Rankin [14] and Cassels [5], with an error noted in [6] and a repair supplied by Ennola [8] and simplified by Diananda [7]. Montgomery proved triangular minimality of the lattice theta function for every positive Gaussian parameter [11, Theorem 1]. These are comparisons within the lattice class; the issue here is the lower bound for arbitrary admissible configurations.

For the confining-potential class in their Theorem 2, Sandier and Serfaty show that microscopic limits of finite Coulomb minimizers concentrate on minimizers of the planar renormalized energy [18, Theorem 2]. Its minimum also determines a next-order coefficient for optimal logarithmic energy on the sphere [1, Theorem 1.5]. Petrache and Serfaty prove a conditional implication from universal optimality to minimality for Coulomb and Riesz renormalized energies in their configuration-energy normalization [13, Theorem 2 and Corollary 4]. Their periodic comparison uses tori compatible with the reference lattice, not arbitrary square tori, and their planar conclusion is conditional on the Cohn–Kumar conjecture. The companion universal-optimality theorem supplies the Gaussian comparison used in that implication [12, Theorem 1.1 and Corollary 1.3]. The present proof treats the logarithmic functional directly and uses no universal-optimality result as an input.

Two geometric precedents help explain the finite-torus argument. Screening a point charge by a uniform disk yields a separation argument for planar confined logarithmic minimizers; see Lieb, Rougerie, and Yngvason [10, Lemma 3 and Appendix A]. Voronoi moment inequalities combine sharp cell comparisons with Euler’s bound on the average number of sides, as in Gruber’s proof of a theorem of Fejes Tóth [9]. Affine bounds on cell energies, sharp at the regular hexagon, also appear in the Wasserstein model of Bourne, Peletier, and Theil for cells with at least three sides and mass v\ge1.5\cdot10^{-4} [3, Lemma 12]. Here the cell comparison must arise from a trial potential that is continuous off the sites and has the prescribed logarithmic singularities. Its nonradial correction and the two scalar estimates of Section 5 address that requirement, including sectors whose perpendicular foot lies outside the corresponding edge.

The whole-plane reduction below uses Sandier and Serfaty’s local energy estimates and periodic approximation. The finite-torus argument provides the comparison with the triangular value. Throughout, the conclusion concerns energy per unit area: compact changes to a configuration can preserve this value, so the minimum statement alone does not classify all minimizing configurations.

1.3 How the comparison works

The proof separates a finite-torus comparison from the passage to arbitrary fields. The finite comparison identifies the triangular value. Sandier and Serfaty’s approximation theorem supplies periodic fields approaching the unknown full-plane minimum; it does not identify that minimum.

Section 2 checks the conversion from that approximation theorem to our normalization and fixed-width square cutoffs. The main point is to control the boundary error for an arbitrary field whose energy limsup is bounded above, without assuming separation or density one in advance. Local energy estimates give the needed charge density and boundary-strip bound. Proposition 2.4 then transfers any common lower bound on square tori to every field in \mathcal A_1, including the possible infinite-energy cases.

We next prove the finite comparison. On each square torus the point energy attains a minimum. A screened logarithmic potential and the minimum principle separate the minimizing points, as shown in Section 3. Separation is used only for that minimum; once it is bounded below, every other configuration is bounded below by minimality. This is why the main theorems impose no separation condition on their competitors.

The lifted minimizer has convex Voronoi cells, the regions closer to one site than to any other. Each edge and its incident cell determine a sector joining that edge to the site. Counting incidences on the torus gives at most six sectors per point on average. The comparison is not an estimate of the electrostatic energy of an isolated cell. Instead, we construct a single trial potential H on the torus. Its sector formulas have matching traces and satisfy H=-\log r+O(r^2) near each site. Thus H and the Green potential h have matching logarithmic singularities, while the constant term is normalized to zero in H. Proposition 4.3 gives the exact identity

Here J(S) is the trial potential’s sector contribution: its linear background-interaction term minus one half of its renormalized Dirichlet integral, as defined in (4.10). The identity cancels the matching logarithmic singularities and retains the nonnegative square exactly.

The construction of H uses radial functions common to every sector. At a shared Voronoi edge its trace depends only on the equal distances to the two neighboring sites; along a ray to a vertex it depends only on the vertex radius. These observations glue the local formulas continuously off the sites and give an H^1 function on every punctured torus. The identity is proved for any compatible flat torus whose area equals its number of sites, not just for square tori. A separate one-site triangular torus supplies six regular sectors; there the chosen potential has exactly the triangular Green gradient, so the square term vanishes.

Section 5 constructs the radial functions from the triangular Green expansion and subtracts fixed multiples of sector area and angle. The desired bound has the form

Areas and angles sum exactly to n and 2\pi n on a square torus. The incidence count then gives n(6K+\lambda+2\pi\mu), and exact calibration on the triangular torus identifies this with nW(E_\triangle).

The perpendicular from a site to a Voronoi edge need not lie on the edge, so one endpoint angle of a sector can be negative. A lower bound for a scalar primitive alone does not control the resulting difference of two primitive values. A second estimate bounds the total negative part of its derivative. Together these estimates handle both ordinary and signed sectors without imposing a regular-polygon hypothesis.

Section 6 proves summable analytic bounds for all omitted lattice modes. Section 7 first verifies the retained finite core and then uses the analytic error bounds to enclose the full functions and their derivatives. Midpoint and interpolation errors cover a compact rectangle between mesh points. At the exact triangular point, analytic stationarity and a positive Hessian bound establish the non-strict minimum; an exact formula handles the unbounded exterior. Section 8 combines these results with the geometry and the whole-plane reduction, then derives Corollary 1.3. Appendix A records two complementary methods: an independent coefficient recurrence and sharper exact tail bounds. They preserve additional quantitative information without adding premises to the main comparison.

2 From square tori to the full plane

This section proves the reduction needed for the full admissible class. A uniform lower bound on finite square tori gives the same lower bound for every field in \mathcal A_1. The argument uses the periodic approximation theorem of Sandier and Serfaty, after checking the normalization and the fixed-width cutoffs in the definition of W.

2.1 Local regularity and the puncture limit

The local integrability convention causes no restriction on the competitors: under the distributional equations, membership in L^1_{\mathrm{loc}} already implies membership in every L^q_{\mathrm{loc}}, 1\le q<2.

In particular the only possible infinite values of W(E) arise in the large-square limsup. No puncture limit is interchanged with that limsup below.

2.2 Approximation and local energy estimates

We distinguish the current normalization from the electric-field normalization of \mathcal A_1. For b>0, write \mathcal C_b for currents j\in L^1_{\mathrm{loc}}(\R^2;\R^2) with

where \Lambda is simple and locally finite. Denote their local puncture functional by \mathsf W(j,\xi); its formula is the same as that of W(E,\xi). Put Q_s=(-s,s)^2.

For such a current, the inverse component rotation F=(j_2,-j_1) satisfies \operatorname{div}F=\rho-b and \operatorname{curl}F=0. In the proof of Lemma 2.1, use the local model \nabla\log|x-p|-(b/2)(x-p) near a charge and -(b/2)x off the charges. This gives the radial singular expansion for F and its rotated expansion for j. Orthogonality preserves the squared norm, so the same finite local puncture conclusion holds for j and every smooth compactly supported cutoff. This local argument uses simplicity and local finiteness, but not the charge-growth bound.

For balls and for Q_s, choose cutoffs 0\le\zeta_{U_s}\le1, compactly supported in U_s, equal to one at distance at least one from its boundary, and with uniformly bounded gradients. Write

Theorem 1(1),(3),(4) of [17] makes the cutoff choice immaterial within each fixed shape family, gives a common finite minimum

and provides square-periodic j_k\in\mathcal C_1 with \mathsf W_{\mathrm{ball}}(j_k)\to M. Their square periods have side lengths \sqrt{2\pi n_k} for positive integers n_k. The common minimum does not assert equality of the ball and square values of an arbitrary fixed field. Their conditions (1.4)–(1.5) require increasing bounded open sets with \bigcap_{s>0}\overline U_s=\{0\}, left-continuous volumes, |U_s-U_s|\le C|U_s|, and, for every fixed a\in\R^2,

They further require constants C and \theta<2 such that

Balls and squares satisfy these conditions with \theta=1 and |U_s-U_s|=4|U_s|.

We also use the following nonsharp consequences of Lemma 4.7 and Proposition 4.9 of [17]. The latter is an application of the mass-displacement method developed for Ginzburg–Landau energies in [16]. In the proof of their Lemma 4.7, Sandier and Serfaty credit the dyadic-shell argument for the local L^r estimate to Struwe [19]; the precise weighted estimate used here is theirs. Let U be bounded and open, \widehat U=U+B_1, and suppose \operatorname{div}j=0 and \operatorname{curl}j=\rho-1 in \widehat U. The number n of simple charges there is finite, and j is locally square-integrable off those charges. For 1<r<2, 0\le\xi\le1 compactly supported in U, and \|\nabla\xi\|_\infty\le C_0, Lemma 4.7 implies

where a_+=\max\{a,0\}. For n\ge1, this follows by replacing the energy by its positive part, enlarging the logarithmic count, and absorbing the fixed cutoff bounds into the constant. The statement of Lemma 4.7 and the final display of its proof differ by a factor two on the energy term; this weakening covers either form. If n=0, there is no puncture term and Hölder’s inequality gives directly

which also implies (2.3). Proposition 4.9 provides one measure g_U, supported in \widehat U and satisfying g_U\ge-C\,\dd x, such that for every smooth compactly supported \xi in U,

where

The constants in (2.4) are universal. In particular the same g_U works for two different cutoffs. The local smoothness needed in these results follows as in Lemma 2.1, after rotating the components.

2.3 Scaling and comparison of fixed-width cutoffs

2.4 Periodic averaging and the harmonic component

For a periodic field F with one unit logarithmic singularity at each of N charges in a period cell P, define its cell energy by

The integral is interpreted on the flat torus, so disks meeting a chosen cell boundary are identified periodically.

On T_n=\R^2/(\sqrt n\,\Z^2), let h be the mean-zero solution of -\Delta h=2\pi(\sum_{i=1}^n\delta_{v_i}-1) for distinct points v_1,\ldots,v_n, with energy defined by (1.6). The compatibility condition holds because |T_n|=n; the mean-zero solution exists and is unique by the Fourier series for the Laplacian on the torus. Any periodic locally integrable curl-free field E with divergence -\Delta h differs from -\nabla h by a smooth field of zero divergence and curl. Each component is harmonic on the torus, hence constant. Therefore E=-\nabla h+c, with c\in\R^2. As a distributional derivative of a periodic function, \nabla h has integral zero; it is in L^1(T_n) by Lemma 2.1. The cross integral on the punctured torus thus tends to zero, giving the exact identity

Applying the same fixed orthogonal rotation to the fields preserves this energy identity.

3 Finite tori and separated minimizing configurations

We first isolate the finite-dimensional problem used in the periodic reduction. Separation will be proved only for its minimizing configurations. This suffices for a lower bound on every configuration.

For an integer n\geq2, put

The torus has area n. Given distinct points v_1,\ldots,v_n\in T_n, let h be the mean-zero distributional solution of

We use the energy (1.6). For this square-torus energy, \eta is smaller than the injectivity radius and half the minimum distance between the finitely many points.

3.1 Green representation and existence

By (2.12), every admissible periodic field with these charges differs from -\nabla h by a constant vector c, which adds n|c|^2/2 to the energy. It therefore suffices to minimize over the mean-zero potentials considered here.

3.2 Separation of minimizers

The argument uses a screened logarithmic potential and the minimum principle. A closely related exclusion argument for minimizing configurations of a planar confined system, with the same single-charge screening radius, appears in Lieb, Rougerie, and Yngvason [10, Lemma 3 and Appendix A]. The following proof treats the torus Green function and its background directly.

4 Voronoi sectors and a glued dual potential

We first record the sector geometry of separated square-torus minimizers. We then prove gluing and duality for compatible periodic Voronoi sectors on any flat torus. This separates the general identity from its square-torus and triangular applications. The construction uses general radial functions; their particular choice and triangular calibration will be given later. Its purpose is to turn a sector lower bound, established by two scalar estimates, into an energy bound for a whole configuration.

The use of polygon inequalities together with an average of at most six sides is familiar from moment methods; see Gruber’s proof of Fejes Tóth’s theorem [9]. A related affine cell bound for a Wasserstein interaction is sharp at the regular hexagon of unit mass; the result of Bourne, Peletier, and Theil [3, arXiv v3, Lemma 12] assumes integer side counts at least three and cell mass v\geq1.5\cdot10^{-4}. Here we must construct a globally admissible singular potential before estimating its sector contributions. The arguments below supply that construction and the incidence count on a torus.

4.1 Counting sectors on the torus

For d>0 and real \alpha,\beta with |\alpha|,|\beta|<\pi/2 and \alpha+\beta>0, define the sector

Its polar axis is the normal from its site toward the supporting edge.

For reference, the area and angle of (4.1) are

In particular the sector has positive angle even when one endpoint parameter is negative.

4.2 Compatible flat-torus data and trace matching

Let \Gamma\subset\R^2 be a rank-two lattice and put \mathbb T=\R^2/\Gamma. Let N\geq1 be an integer, let v_1,\ldots,v_N be distinct points of \mathbb T, and let P be their full \Gamma-periodic lift. Suppose that the Voronoi cells of P are bounded convex polygons with finitely many face–edge incidences modulo \Gamma. Assume that each incidence has the sector form (4.1) and that, for one named constant d_0>0, every sector has d\geq d_0. We call these compatible periodic Voronoi data with lower edge distance d_0. This is conditional geometric data, not a separation assertion about an arbitrary competitor.

The sectors partition the N cells modulo \Gamma, up to their boundary rays. Their total area is |\mathbb T|, and their total angle is 2\pi N: each convex cell has its site in its interior because d_0>0, so its sectors have total angle 2\pi. These facts do not require an incidence bound.

Call a radius \eta>0 admissible if it is smaller than the injectivity radius of \mathbb T and, when N\geq2, smaller than half the minimum torus distance between distinct sites. Thus the torus disks B(v_i,\eta) are embedded and pairwise disjoint; for N=1 there is no pairwise-distance condition. Write

In puncture limits we may further take \eta<d_0, so every deleted disk lies in its Voronoi cell. This restriction does not change a limit as \eta\downarrow0.

Fix compatible periodic Voronoi data with lower edge distance d_0. Let I be a finite or countable index set, and choose integers k_i\geq6, i\in I; different indices may have the same integer. Suppose that

satisfy the following assumptions:

and, for every T_0>d_0,

These are the exact inputs needed in the gluing and duality argument. The later radial construction will verify them directly.

For a sector (4.1), define t(\varphi)=d\sec\varphi and

The same radial functions are used in every sector.

The common notation D_i=b_i-u_i, B=\sum_i u_i is one way to supply these data. In that notation the common edge trace is simply \sum_i b_i(r). The gluing argument itself only uses B and D_i.

4.3 The exact dual identity

Assume now that |\mathbb T|=N. Let h be the mean-zero distributional solution on \mathbb T of

The area condition is exactly the compatibility condition for this equation. With the convention

the Fourier mode e^{2\pi i k\cdot x} has eigenvalue 4\pi^2|k|^2 under -\Delta. The constant Fourier coefficient of the right-hand side vanishes. Dividing its coefficient at each nonzero k\in\Gamma^* by 4\pi^2|k|^2, and setting the coefficient at zero to zero, gives a distributional solution. It is unique because the kernel of the Laplacian on the torus consists of constants. Define

where the limit is taken through admissible radii. The local calculation in Lemma 2.1, applied in a torus chart, gives the finite puncture limit. When \mathbb T=T_n and N=n, this is exactly \mathcal W_n(h) in (1.6).

For a sector S=S(d,\alpha,\beta), define

Lemma 4.2 shows that this limit is finite. Its counterterm uses the positive angle \theta(S), also for a sector with a negative endpoint parameter. Define

The sectorwise estimates in the proof of Lemma 4.2 also show that these expressions are well defined on any individual sector with d\geq d_0 and the endpoint conditions in (4.1), without first choosing a periodic tiling.

If the gradient hypotheses were unchanged but the value expansion were instead H=-\log r+c_i+o(1) at v_i, the right-hand side of the last identity would contain the additional term 2\pi\sum_i c_i. This explains why the zero constant term in the radial construction is part of the normalization.

For a square-torus minimizer, Lemma 4.1 supplies the compatible data with \Gamma=\sqrt n\,\Z^2, N=n, and d_0=d_*=r_*/2=1/(2\sqrt\pi). Its torus has area n, so Proposition 4.3 applies with \mathcal W_{\mathbb T}=\mathcal W_n. The triangular application in the next section supplies its own period and sectors and does not use a square-torus minimizer.

4.4 From scalar estimates to configuration bounds

The two scalar bounds below control ordinary and signed sectors by different arguments. Once they yield an affine lower bound for each sector, the dual identity and the incidence count give the energy bound. We formulate this implication for compatible flat tori first; separation of a square-torus minimizer then removes the geometric hypotheses from the finite-volume comparison.

The next section chooses particular radial data and exact constants K<0,\lambda,\mu. For these data, the primitive identity (5.15) and the two scalar bounds yield

for every relevant sector. Proposition 5.4 applies the compatible-torus conclusion to this construction and identifies 6K+\lambda+2\pi\mu with the triangular energy. Section 8 uses the square-torus conclusion, whose minimizer argument has already removed the separation hypothesis.

5 Calibration on the triangular cell

We construct the radial functions used by the Voronoi test potential. The potential has the same gradient as the triangular Green function on a regular hexagon, while its energy on every other sector reduces to two scalar estimates. The continuation outside that hexagon is auxiliary; its explicit form makes those estimates accessible to interval arithmetic.

5.1 The local triangular Green function

Set

Thus p and R_0 are the inradius and circumradius of the covolume-one triangular Voronoi hexagon, and x_0=p\tan(\pi/6). Identify the plane with \mathbb C, and write \partial_z=(\partial_x-i\partial_y)/2, with a nearest lattice point on the positive real axis. For m\ge1 put

These sums converge absolutely and are real by reflection symmetry. We also introduce an auxiliary index 0, with k_0=6. It is distinct from index 1: both exponents equal 6, but they will carry different radial functions.

5.2 Matching the regular hexagon

Write \gamma(r)=-\log r+\pi r^2/2 for r>0. We first choose the radial data on the regular hexagon, where the test potential must equal h_\triangle-C_\triangle. In a sector whose outward normal is the positive real axis, the edge has distance p and the endpoint angles are \pm\pi/6. At an edge point of radius t=p\sec\varphi,

The boundary value of the mth Green-function mode is therefore c_mt^{k_m}\cos(k_m\varphi). Expressing it as a function of the radius gives the polynomial

On either vertex ray, the same mode equals (-1)^mc_mr^{k_m}, since k_m\pi/6=m\pi. This is the radial baseline we use for that mode inside the hexagon. The correction in (4.5) then recovers the angular dependence exactly:

Thus the edge data b_m=V_m and the vertex-ray baseline determine the desired potential on each regular sector. Its remaining part \gamma is already radial and needs no correction there.

5.3 Continuation to arbitrary sectors

Other Voronoi sectors can have edge radii below p or above R_0. We now extend the radial data to those radii, keeping two derivatives and making the correction vanish beyond a fixed radius. The assembled potential is required to agree with h_\triangle-C_\triangle only inside the regular hexagon; the continuation is chosen for the scalar estimates, not as an extension of that Green function.

Define the clipped polynomial

Let h_s=3/40 and c=3/50. For 0\le z\le h_s put

For a function V with two left derivatives at R_0, its continuation is defined for r\ge R_0 by

Indeed (P,P',P'')(0)=(0,1,0) and (Q,Q',Q'')(0)=(0,0,1); all six values at h_s vanish. Thus this operation preserves two derivatives at R_0 and joins a constant with two derivatives at R_0+h_s.

Define b_0=\gamma up to R_0, continued by (5.5). For m\ge1, set b_m=V_m on [p,R_0] and use \mathcal HV_m above R_0. Below p put

where

Direct differentiation of (5.4) gives T_{0m}=V_m(p), T_{1m}=V_m'(p), and 2T_{2m}=V_m''(p). Consequently every b_i is C^2 away from the origin.

Now define

and set

Here the auxiliary mode 0 lets the radial part \gamma join its constant continuation by the same correction mechanism as the other modes. Its exponent k_0=6 makes this correction vanish at the pole without changing the singular normalization. The coefficient estimates proved below justify the sums and their first three derivatives on every open spline interval. In particular, B,D_i are C^2 and piecewise C^3 away from zero. Near zero, B=\gamma+\sum_{m\ge1}(-1)^mc_mr^{k_m}; normal convergence of (5.3) therefore gives

The fixed radial data also satisfy (4.4) for every d_0>0. Indeed, there are only finitely many retained modes, while for the omitted modes (6.4) and the vanishing of D_i beyond l bound the summand on [d_0,T_0] by (1+k_i+T_0)w_{k_i}. The explicit geometric terms in (6.3) make this series summable. This conclusion uses the analytic mode bounds, not the later finite-core estimates.

For a sector with edge distance d and endpoint angles -\alpha<\beta, write t(\varphi)=d\sec\varphi and use the test potential

At the edge its trace is \sum_i b_i(t); at a ray ending at a vertex of radius r_v its trace is B(r)+\sum_i(r/r_v)^{k_i}D_i(r_v). These expressions explain both trace matching conditions in the glued Voronoi construction.

5.4 The sector integral

The following formulas separate radial integration from sector geometry. For t>0 define

Here H_k is indexed by the exponent, so the two modes with exponent six use the same H_6. They remain separate in every sum involving D_i. Put

All sums here and below include index 0. The angular coefficient is nonnegative, since

Let J denote the sector contribution to the dual energy,

where the renormalized integral adds (\alpha+\beta)\log\eta before letting the radial exclusion \eta tend to zero. Then

To verify this identity, differentiate (5.9):

The source integral gives S_1+t^2\sum_iD_i/(k_i+2). The radial square contributes \log t+S_2, its cross terms contribute -2D_i+2k_iD_iH_{k_i}, and its quadratic terms have coefficient k_ik_j/(k_i+k_j). Integrating the angular square gives \tan^2\varphi\sum Q_iQ_j/(k_i+k_j). These are precisely (5.11). The estimates below permit termwise integration.

5.5 Subtracting area and angle

Given real constants \lambda,\mu, define

The evenness of G in s makes g_d odd. If A is the sector area and \theta=\alpha+\beta, substitution s=d\tan\varphi in (5.13) gives

This identity holds for signed endpoint angles; only |\alpha|,|\beta|<\pi/2 and \alpha+\beta>0 are required.

We choose \lambda,\mu so that the adjusted primitive is stationary at the regular sector (d,x)=(p,x_0), and set K to twice its value there. The exact choices are

The first integrand extends continuously to s=0 with value f_0'(p)/(2p).

5.6 The two scalar estimates needed by the geometry

The remaining analytic task is to establish

The rational endpoint b=28209/100000 is below the geometric threshold d_*=1/(2\sqrt\pi). Thus the domain d\ge b covers all edge distances supplied by the separated-minimizer geometry; the elementary comparison is justified in Section 8.

The second estimate is essential when the perpendicular from a point to a Voronoi edge falls outside that edge. Its role is made precise in the following implication, which completes the geometric part of the argument once (5.18) has been certified.

6 Infinite-mode truncation and error propagation

This section supplies analytic estimates for the modes omitted after m=38. We use the radial functions b_m,U_m,u_m,B,D_m of the preceding construction, with

Thus k_0=6, k_m=6m for m\geq1, and

Tildes will always mean that the modes m\geq39 have been omitted; the index-zero term is retained. Differences between a full quantity and its truncation are denoted by \Delta. Our goal is an absolute error below 10^{-11} for f_0,e and their first two derivatives on [b,l], conditional on the retained-mode bounds verified in Section 7.

6.1 Coefficient estimates

Write a lattice point as z=a(n+q/2+i\sqrt3 q/2), with (n,q)\in\Z^2. Then

There are 8v integer pairs with \max(|n|,|q|)=v. For every d=6m\geq6 and every positive integer D, the sum outside the square |n|,|q|\leq D, before multiplication by a^{-d}, is therefore at most

The last inequality follows by comparison with the integral from D to infinity. The first shell has six points with squared norm 1 and two with squared norm 3. Moreover,

It follows that, for d=k_m\geq12,

Indeed, the bracket equals 6+50/729<7 at d=12, and its two nonconstant terms decrease with d. The finite coefficient enclosures may use (6.1) with D=256 for m=1 and D=48 for 2\leq m\leq38, including the factor 1/(k_m a^{k_m}). Reflection symmetry makes each coefficient real.

6.2 Uniform derivative bounds for omitted modes

For k=6m\geq234, define the rational number

For 0\leq j\leq3 and 0<r\leq l, the construction satisfies

At spline knots these inequalities mean both one-sided derivatives. Here are explicit checks of the constants.

The polynomial cutoff S(y)=10y^3-15y^4+6y^5 on [0,1], extended constantly outside that interval, obeys

The identities p^4=1/12 and R_0^4=4/27 give the convenient rational bounds

All terminating decimals in this section denote rational numbers.

For U_m(r)=(-1)^m c_mr^k(1-S((r-R_0)/.1)), Leibniz’ rule and (6.2) bound the third derivative by

Using l>.69 and k\geq234, the coefficient of k^2(13/20)^k is less than 25.445<30. The same calculation for orders 0,1,2 gives coefficients below 0.000001, 0.000197, and 0.070572, respectively. This proves the first term inside the parentheses in (6.3).

For p-c\leq r\leq p, write q_m(r)=T_{0m}+T_{1m}(r-p)+T_{2m}(r-p)^2, so that b_m(r)=S((r-p+c)/c)q_m(r). Since p/a=1/2 and p>1/2,

The last inequality uses ((k-2)(k-3)/6)-1/2\leq k^2/6. Consequently

In particular,

The lower derivative orders satisfy the same final bound, and below p-c the function is zero.

On [p,p+.01], put w=r^2-p^2 and X=p+\sqrt w. The polynomial formula for V_m gives, for 0\leq\ell\leq3,

For j\geq\ell, the coefficient comparison that proves this is

The bound p>.537 implies .02p+.0001<.04p^2, hence X\leq1.2p=(3/5)a on this interval. Also r<.55 and X>1/2. The chain rule gives

Thus its absolute value is at most

The lower derivatives satisfy smaller bounds. The common bound 5600k^5(3/5)^k therefore applies through order three.

On [p+.01,R_0], use z=p+i\sqrt{r^2-p^2}, so that V_m=\Re(c_mz^k) and |z|=r. Direct differentiation gives

For example, \sqrt{r^2-p^2}>.103, p<.54, and r<.63 suffice in the formulas |z'|=r/\sqrt{r^2-p^2}, |z''|=p^2/(r^2-p^2)^{3/2}, and |z'''|=3p^2r/(r^2-p^2)^{5/2}. Since 1/r<2, the third derivative of c_mz^k is bounded by

All derivative orders through three are therefore bounded by 3\cdot10^6 k^2(289/500)^k on this interval.

For the extension \mathcal HV_m=V_m(R_0)+P(r-R_0)V_m'(R_0)+Q(r-R_0)V_m''(R_0), direct bounds on the coefficients of P,Q show, for 0\leq j\leq3,

Thus the common bound for the extension is at most

It remains valid on the constant part beyond R_0+h_s.

Finally the cutoff \tau(r)=S((r-R_0-.001)/.075) satisfies

The maximal Leibniz factor through order three is 1+3\cdot27+3\cdot1067+142223=145506<10^6. Applying this to u_m=(1-\tau)U_m+\tau b_m and D_m=(1-\tau)(b_m-U_m) proves (6.4). The series and their indicated derivatives consequently converge absolutely and locally uniformly. For r<b, the omitted u_m derivatives additionally have their original O(r^{k-3}) behavior.

6.3 Norms and bounds for retained modes

For a function F with one-sided derivatives as needed, define

These norms are submultiplicative. For nonnegative weights \omega_i, a weighted mode norm means \sup_{t\in I}\sum_i\omega_i|F_i|_{j,t}: the sum is taken before the supremum.

The finite interval verification in Section 7 establishes the following bounds for the retained modes:

The finite coefficient enclosure also supplies |c_1|<1. These assertions are not conclusions of the analytic estimates in this section. Conditional on them, the remaining estimates below are independent of the finite calculation.

On [0,b],

has norm |\widetilde B_1|_{2,t}<10. For an explicit check, use t<b<2/7, a>1, \pi<22/7, and |c_1|<1 for the first mode. The modes k\geq12 contribute at most \sum_{k=12,18,\ldots}7k^2(2/7)^{k-3}<.014; the first term and its geometrically decreasing ratios prove this last rational inequality. Together with the linear term and the k=6 term this is less than 10.

6.4 Summing the omitted modes

Set

Here is a rational verification of the numerical inequality. For each term (k+1)^4k^sr^k, with

the ratio at the next allowed index is

It decreases with k, and at k=234 it is strictly below 1/10 for all four pairs. The first weighted summand obeys the exact rational comparison

Multiplying by 10/9 gives a bound less than 2\cdot10^{-21}, as asserted in (6.6).

Since \Delta B_1=\sum_{m\geq39}u_m', the derivative estimates give

and, on [b,l],

Indeed the divided-factorial sums use factors 5/2 and 8/3, respectively, both less than 4.

6.5 Propagation to the sector functions

Recall the quantities entering the sector formula:

Write N=2\cdot10^6 temporarily. Direct differentiation of \Delta S_1 gives

For example |\Delta B^{(j)}|\leq W_*/4 through order three, so the norm is at most [t^2/2+t+(1+t)/2]W_*/4 for t\leq l<1. For the other integral, use

The second-order norm of its integral is at most twice the supremum of the first-order norm of the integrand, because l<1. Submultiplicativity and |t|_{1,t}<2 yield

The identities and elementary bounds

imply, for b\leq t\leq l,

Indeed l/(k+1)+5/2+k/b<6(k+1). The bound for \Delta H_k follows by the same linear argument applied to \Delta B_1.

The radial source and energy functions are

The preceding bounds give

For completeness, |t^2|_{2,t}=(t+1)^2<4 gives the 16\mathfrak a source term. The sum of -2D_i contributes 8\mathfrak a. Split the mixed term as

The finite and omitted weighted norms give 12NW_* and 48(6000+W_*)\mathfrak a after including its factor 2. Bounding k_ik_j/(k_i+k_j) by k_ik_j, the quadratic difference is at most 8N\mathfrak a+16\mathfrak a^2. These are precisely the terms in (6.8). Finally W_*\leq4\mathfrak a, \pi<22/7, and \mathfrak a<2\cdot10^{-21} make its total coefficient less than 56,264,181<6\cdot10^7.

Since k_i\geq6 and t<1, direct differentiation of Q_i shows

For example, its coefficients on |D_i|,|D_i'|,|D_i''|,|D_i'''| are bounded respectively by k_i, t+k_i-1, t+(k_i-2)/2, and t/2. The finite and omitted sums of these norms are bounded by 2N and 8\mathfrak a. Dropping the denominator k_i+k_j\geq1 in the energy difference therefore gives

The norm contains half the absolute second derivative. Accounting for this factor, (6.8) and (6.9) imply absolute truncation errors less than 10^{-11} for every member of

In fact the respective second derivative bounds are below 2.4\cdot10^{-13} and 1.3\cdot10^{-13}.

6.6 The constant region

For r\geq l, all D_i vanish and B is constant. The omitted constant satisfies

Since B_1(t)=1/t there, the derivatives of \log t and S_2(t) cancel. The full infinite-mode functions consequently obey the exact formula

Thus outside [b,l] one propagates the enclosures of f_0(l) and B(l) through this exact formula; no uniform absolute error for arbitrarily large t is asserted.

7 The finite scalar verification

The geometric argument requires the two inequalities (5.18) for the exact functions G and g_d. We now specify a finite calculation that implies them. Its inputs are the radial data (5.7), the formulas (5.11), and the constants defined by the integrals (5.16). Its outputs are interval bounds for these constants, two strict bounds on a compact rectangle, and a positive Hessian bound near the exact stationary point. The estimates of Section 6 connect the finite radial calculation to the infinite series.

Every terminating decimal in this section denotes an exact rational number. An interval enclosure always contains the exact mathematical quantity in question. In particular, the constants \lambda,\mu,K are not replaced by rounded numerical choices. Their enclosures are used in all subsequent operations, while Lemma 5.2 supplies stationarity for their exact values.

7.1 Outward arithmetic and derivative intervals

The calculation follows the interval inclusion principle: every operation encloses its exact result, and analytic remainders are added explicitly. Revol and Rouillier [15] discuss this principle and the dependence of enclosure width on expressions and subdivision. We specify the integer operations here, so the proof’s arithmetic assumptions and rounding rules are explicit.

Put q=2^{144}. A scalar interval is represented by integers u\le v, meaning [u/q,v/q]. Rational inputs are rounded down at the lower endpoint and up at the upper endpoint. Addition and negation are exact at this denominator. Multiplication takes the least and greatest of the four endpoint products, then divides by q with outward rounding. If 0<u\le v, inversion gives

negative intervals are treated by negation, and division by an interval containing zero is forbidden. Integer powers use repeated multiplication. For a positive interval, integer square roots give the enclosure

The extra unit in the upper endpoint is harmless when the upper square root is an integer.

For logarithms, first extract an integer power of two and reduce the argument to [1/4,1]. The arguments used below are positive table intervals or enclosures of fixed positive constants; their widths allow this reduction. With z=(v-1)/(v+1), evaluate

The same series at z=1/3 encloses \log2, used to restore the extracted power. All summation errors are already enclosed by outward arithmetic; the displayed remainder is added separately. Similarly, use

The alternating-series bound gives |\mathcal R_v|\le v^{-301}/301<q^{-1} for both values of v. Adding [-2q^{-1},2q^{-1}] to each arctangent enclosure before the linear combination is sufficient. These procedures also enclose p,R_0,x_0,l using their exact defining formulas.

To bound derivatives we use interval jets

Addition, multiplication, inversion, logarithms, and square roots are performed as formal power series through degree three. For example, the product coefficient is \sum_{r=0}^jF_rG_{j-r}, and the inverse coefficients follow recursively from the product being one. A variable on a radial interval T is represented by (T,1,0,0); the resulting coefficients enclose its derivatives at every point of T on the chosen spline piece. At a single point, the leading interval encloses that point and the same derivative entries are used. After differentiation only the coefficients supported by the available derivative orders are used. The calculation below never requires a fourth derivative of a spline: third derivatives of B,D_i suffice for the second derivatives of f_0,e and for the radial quadrature errors.

For a differential equation Y'=F(t,Y) with an enclosed value at the evaluation point, its jet coefficients are reconstructed in increasing order by

Thus the enclosed value and the known differential equation supply the derivatives without numerical differentiation.

Spline polynomials are evaluated by interval Horner arithmetic. One may first translate a polynomial P(t)=\sum_{\ell=0}^Na_\ell t^\ell to the midpoint m of its leading interval:

This is an exact polynomial identity, also for jets. It reduces interval widths without changing the spline or introducing an approximation; direct Horner evaluation remains a valid enclosure.

7.2 Coefficients and the radial table

Retain the auxiliary mode 0 and modes 1\le m\le38. The two modes with exponent six remain separate in sums involving D_i. Enclose each c_m by summing over |n|,|r|\le D, where D=256 for m=1 and D=48 otherwise, and then adding the symmetric error in (6.1), multiplied by a^{-k_m}/k_m. The summands can be evaluated using only integers before that last factor: if k=k_m, A=2n+r, and N=n^2+nr+r^2, then

The pair (0,0) is omitted. The resulting enclosure also tests |c_1|<1.

The radial breakpoints are

Each intervening interval is subdivided by affine interpolation between its exact endpoints. For the first interval use n=\lfloor8000(u-v)\rfloor+1 subintervals, and for all later intervals use n=\lfloor24000(u-v)\rfloor+1, where v,u are its endpoints. In execution an outward upper bound on the product is used to choose n, which can only add subdivisions. All true node positions are enclosed by the corresponding interval arithmetic. The length of every actual subinterval is bounded above in the calculation; these upper bounds are used in all error estimates.

For the fixed table, write [u_i/q,v_i/q] for the valid outward enclosure of its ith node, 0\le i\le m. The interval constructor enforces u_i\le v_i, and the node arithmetic preserves this invariant. After appending the last node and before propagating any integral state, the program requires

it throws an exception if a comparison fails. Thus, on a completed table construction, u_i\le v_i<u_{i+1}\le v_{i+1}. Both endpoint arrays are strictly increasing, and the stored step interval is

with positive lower endpoint. This also justifies every division by H_i below. These facts concern the initialized table and its unchanged aligned node, point, and cell arrays in the fixed calculation.

On 0<t\le p-.06, the truncated functions have the monomial form

For compact formulas introduce the finite list

This list indexes monomials of \widetilde B_1, rather than the sector modes; its exponent e_0=2 is unrelated to k_0=6. The starting integrals, including their derivatives, are obtained from the exact identities

In particular, the table starts propagation at the exact node p-.06 with these values; preceding nodes and cells use (7.8) directly.

There is a useful exact simplification for the sector functions on this first interval. Since b_m(t)=0, D_m(t)=-a_mt^{k_m} and D_0(t)=0, substitution in (5.9) cancels every nonconstant mode and gives H(r,\varphi)=\gamma(r) for r\le t. Equivalently, every Q_i vanishes. Direct radial integration gives

Indeed,

Combining these with the source and energy factors in the sector formula proves (7.9). The identity holds for the full infinite series as well as for every finite truncation. It may therefore replace the mode sums when evaluating f_0,e and their first two derivatives on the first interval, producing tighter enclosures. The full finite-mode values (7.8) at p-.06 are still required to start the continuation of S_1,S_2,H_k.

For a subsequent radial step [v,u], propagate the integrals by

For each integrand Y, the midpoint rule has the rigorous enclosure

Midpoint values use point jets; second-derivative bounds use a jet whose leading interval contains [v,u]. The formula applies to weak second derivatives bounded almost everywhere as well.

For values inside a step, use the integral equations rather than the midpoint approximation. If T=[v,u] and h=u-v, then

where Y_1=r\widetilde B and Y_2=-2\widetilde B_1+r\widetilde B_1^2. For the last inclusion, (r/t)^k\le1 whenever v\le r\le t\le u. Recover the positive jet orders from

Substitution into (5.11) produces both node enclosures and cell enclosures of the ordered six-tuple

No numerical differentiation of table values is used.

Separately, subdivide each interval in (7.6) into steps of length at most .00022. Direct interval jets for \widetilde B_1 and \widetilde D_i test the two finite-core inequalities (6.5). Together with |c_1|<1, these are the finite premises of the truncation analysis; their verification does not use any infinite-mode error bound. Only after these tests pass is (6.10) invoked.

The spline data are C^2 and have bounded one-sided third derivatives. Accordingly f_0,e are C^1 with locally bounded weak second derivatives on (b,\infty). A radial cell uses its own spline formula at its endpoints; either one-sided second derivative is valid for its interior. When a query meets a knot, the hull includes both adjacent cells. Thus no jump of a second derivative is mistaken for global classical differentiability.

7.3 Queries for the infinite-mode functions

Let \epsilon=[-10^{-11},10^{-11}], and let T be a valid outward enclosure of the queried true radii. In the fixed production calls every such radius is at least b, although T may extend below b. For derivatives, take the componentwise hull of (7.14) over every radial cell meeting T and add \epsilon to each component. By (6.10), this encloses the full infinite-mode tuple (f_0,f_0',f_0'',e,e',e'') at every queried true radius in T\cap[b,l]. The exact continuation below covers queried true radii at least l. Here is the cell-selection argument for the implementation. In integer endpoint units write the valid query as [a,b'], let j be the first index with v_j\ge a, and let k be the first index with u_k>b', using m+1 if either index is absent. The closed envelope of cell i is [u_i,v_{i+1}]. It meets [a,b'] exactly when v_{i+1}\ge a and u_i\le b'. By the ordering proved above, these conditions are equivalent to i\ge j-1 and i<k. Therefore the binary searches on the upper and lower endpoint arrays select exactly

The non-strict intersection inequalities include both cells at a shared endpoint. Every true cell meeting the query lies in its envelope and is therefore selected; any extra cells only widen the hull. This is a statement about the fixed initialized table, not a general-purpose lookup interface.

For a query requiring values alone, improve the enclosure by linear interpolation of the endpoint values. On a radial cell [v,u] put

An empty intersection is discarded. The interpolation error is

The cell enclosure supplies the derivative bound. Add this symmetric error and \epsilon, then take the hull over all relevant cells. Shared coefficient uncertainties need not be independent: interval operations enclose all their possible values and hence also their correlated exact values. In particular, the positive stored interval H_i contains the true step u-v, so outward evaluation of the quotient defining X contains the normalized coordinate of every true queried point in that cell. Intersecting with [0,1] cannot discard such a point.

Queries that meet [l,\infty) also include the exact continuation (6.11). Enclose its two constants by the last node value plus \epsilon and by

Thus the continuation contributes the tuple

Using the entire query interval in this expression is conservative even if it extends below l, since the relevant exterior subset is contained in it. This rule propagates the constant uncertainties through the exact formula and asserts no radius-independent error on the unbounded tail.

7.4 Enclosing the calibrated constants

We evaluate the exact integrals defining \lambda in (5.16). Put s_0=.01 and t(s)=\sqrt{p^2+s^2}. Taylor’s formula at p removes the apparent division by zero in the first integral. Specifically,

Indeed, t-p=s^2/(t+p), while

The quadratic Taylor remainder of f_0 supplies the second line of (7.18). The derivative range comes from the radial cells with the infinite-mode error added.

The implementation obtains this range by a direct cell scan over [p,t(s_0)], and obtains the derivative bounds below by a second direct scan over [p,R_0]. Both use the same closed-cell selection rule, with outward endpoint bounds, but neither adds the exterior formula. This is sufficient here: p^4=1/12 gives p>1/2>b, x_0=R_0/2>s_0, and

Hence both true scan domains are contained in the finite table [b,l]. A nonempty scan alone would not justify a query extending beyond that table.

On [s_0,x_0] integrate p(f_0(t(s))-f_0(p))/s^2 by a midpoint sum with 45000 equal steps. Integrate e(t(s))/p on [0,x_0] by a separate midpoint sum with 45000 equal steps. Obtain integer upper bounds

from the cell enclosures plus \epsilon. The two integrands have second derivatives bounded respectively by

For the first bound, set V=f_0(t(s))-f_0(p); then

and differentiate pV/s^2 twice. The resulting terms are p(V''/s^2-4V'/s^3+6V/s^4), bounded by the first expression in (7.19). The second follows by twice differentiating e(t(s))/p. For an interval of length L and equal step length h, add the symmetric error Lh^2\sup|Y''|/24 to each complete midpoint sum. Together with (7.18), the endpoint terms and the factors in (5.16), this encloses \lambda. The same formula then encloses \mu.

To enclose K=2\int_0^{x_0}G(p,s)\,\mathrm ds, choose an integer n and equal steps of length h=x_0/n\le.0001, and use their midpoint sum. Independently enclose G_{ss}(p,s) on a partition with step length at most .0002. If h is the midpoint step length and L_K bounds this second derivative, the enclosure is

The factor two multiplies both the sum and its error.

The first sign tests are

Each strict inequality is accepted only when the corresponding interval endpoint proves it. The bound on \lambda also serves as a check on the normalization. The last two signs and (6.11) give

Here e(t)=0 and d/t^2>0.

7.5 The compact rectangle and its derivative bounds

Because l<.697 and \sqrt{b^2+.65^2}>l, the signs (7.22) reduce the remaining work to

For d\ge.697, G(d,s)>0 for every s\ge0 and g_d(x)\ge0. For d\ge b and x\ge.65, the primitive is increasing in x. The negative part of G vanishes for s\ge.65.

We obtain derivatives of G by differentiating its defining formula, using the full radial query enclosures. The following identities specify all needed factors. Put

and write f'=f_0'-\lambda t, f''=f_0''-\lambda. Then

For v=d or v=s let r_v have the meaning above. The derivative formulas are

These formulas hold on each smooth piece and almost everywhere across the radial knots. No mixed second derivative of G is needed.

Partition both sides of (7.23) into intervals of length at most .0008. Interval substitution in (7.25) on every closed rectangle supplies global upper bounds

The value-only interpolation is not used to obtain these derivative bounds.

7.6 The primitive away from its stationary point

Let

Its interior contains (p,x_0). To verify the primitive bound on the complement, partition the d interval at the two vertical sides of \mathcal D_0, and partition the s interval at its two horizontal sides. Subdivide these pieces into steps with maximal lengths

Every cell therefore lies either in \mathcal D_0 or outside its interior. For each d node, form successive midpoint sums from 0 to every x node, using the actual step lengths:

The sum at x=0 is zero and is tested as well. The value queries in Section 7.3 may be used here. At every node outside the interior of \mathcal D_0, enclose 2Q(d,x)-K and let m_* be a lower bound for all these enclosures.

The quadrature error at grid nodes is bounded by E_0, and the additional interpolation error between nodes by E_I, where

Indeed, the bilinear interpolation estimate follows from

Apply the one-dimensional interpolation error in each variable successively; the resulting bounds add and require no mixed derivative. Thus the finite assertion

proves 2g_d(x)>K throughout \mathcal D\setminus\operatorname{int}\mathcal D_0. All the derivative estimates are weak derivative estimates, so the same conclusion holds when a quadrature segment or rectangle crosses a radial knot.

7.7 The negative part on the whole rectangle

At every d node, including those between the vertical sides of \mathcal D_0, use the same s partition to enclose

Let N_* be an upper bound for all these sums. The function y\mapsto\max\{0,-y\} is 1-Lipschitz. On a subinterval of length h, its midpoint error is therefore at most L_2h^2/4: integrate L_2|s-s_{\rm mid}|. Summing the errors gives at most .65S_hL_2/4. Every d lies within D_h/2 of a node; changing d to that node changes the integral by at most .65D_hL_1/2. Consequently the assertion

proves

Together with (7.22), this gives the entire negative-part inequality in (5.18).

7.8 Convexity in the stationary rectangle

It remains to prove the primitive inequality on \mathcal D_0. Partition its d interval into steps of length at most .00004. Partition [0,.32621] into steps of length at most .00008, including .29420 as a node. The accumulation in s must begin at zero even though the Hessian is tested only above .29420.

For a d interval D and an s interval S_j, use the full derivative queries to enclose

For a cell D\times S_j meeting \mathcal D_0, enclosures of the three Hessian entries of (d,x)\mapsto g_d(x) are

The first line is direct interval integration, not midpoint quadrature: it bounds both the completed subintervals and the partial final subinterval. Test, on every such cell, that

Interval multiplication encloses the determinant even when its factors are correlated. These inequalities imply positive definiteness, and in particular positive semidefiniteness, of the Hessian almost everywhere on \mathcal D_0.

The function (d,x)\mapsto g_d(x) belongs to W^{2,\infty} on a neighbourhood of this rectangle. To justify the passage from Hessian bounds to convexity, mollify on each smaller interior rectangle: the mollified Hessian is positive semidefinite, so the mollified functions are convex. Local uniform convergence proves convexity in the interior, and continuity includes the boundary. By Lemma 5.2, (p,x_0) is an interior stationary point. Hence

This step uses the exact stationary point and does not infer convexity or stationarity from sampled function values.

7.9 Finite assertions and the scalar conclusion

For clarity, the finite procedure has the following fixed parameters and obligations.

Finite bounds. The program certificate.cpp implements this procedure using core.hpp and grids.hpp. All three files and an execution wrapper are included in the verification directory. One completed guarded run is documented in verification/README.md, which binds the code and formula hashes, run metadata, saved integer pairs, and their interpretation below. The wrapper compiles and runs the computation; it does not read the TeX manuscript. The source-to-code and analytic correspondences are separate arguments. The interval denominator is 2^{144}. The table gives outward-rounded projections, while interval assertions use uncoarsened endpoints.

Here the determinant is enclosed on each individual rectangle before the bounds are combined; it is not computed from separately combined entry ranges. The global derivative bounds, rounded outwards for display, are

These derivative numbers are projections of the computed hulls, not separate executable threshold assertions. For the identified sources, under the documented runtime assumptions, the completed run verifies the listed finite assertions using uncoarsened interval endpoints. Proposition 7.1 explains why these assertions imply the required bounds for the exact infinite-mode functions: the interval formulas and range-coverage arguments identify the quantities enclosed, while the analytic truncation estimates supply the omitted modes. Neither the rounded table nor the execution record alone supplies those mathematical correspondences.

8 The finite-torus and whole-plane bounds

We now apply the geometric comparison to the certified scalar bounds, then pass from square tori to the plane. The square-torus conclusion of Proposition 4.5 already incorporates the separation of minimizers and the passage to arbitrary configurations.

The comparison establishes the minimum of this logarithmic renormalized energy. It does not classify all configurations having that energy per unit area. In particular it supplies no uniqueness or defect-exclusion assertion.

8.1 The spherical logarithmic-energy constant

A Supplementary coefficient and tail estimates

The proof uses the coefficient enclosures and omitted-mode estimates in Sections 6 and 7. This appendix records an independent integer recurrence for the same finite coefficient sums and alternative exact bounds for the same omitted-mode series. These arguments are supplementary: they do not add a hypothesis to the scalar theorem or replace the analytic derivative bounds and finite retained-mode premises used there.

A.1 An independent recurrence for the lattice coefficients

Fix integers (n,q)\ne(0,0) and write

Thus w=n+q/2+i\sqrt3q/2, N=|w|^2>0, and |\zeta|^2=4N. A point of the covolume-one triangular lattice is aw, where a=\sqrt{2/\sqrt3}. For m\geq0, put

The two numbers u=\zeta^6 and v=\overline{\zeta}^{6} have

and hence satisfy X^2-2r_1X+(4N)^6=0. Multiplying this equation by u^{m-1} and v^{m-1}, respectively, and adding gives the integer recurrence

Here all the quantities on the right are integers. This proves the recurrence for every individual lattice point, without a numerical approximation to a complex power.

The binomial theorem identifies its numerator exactly. Since only even powers of i\sqrt3q have nonzero real part,

This is the numerator in (7.5), after writing q for its second integer coordinate. Indeed, w^{-6m}=\overline{w}^{6m}/N^{6m}, and therefore

Thus the recurrence and the binomial expression agree for every point and every mode, not only for a selected finite list of comparisons.

Grouping by the norm. For a positive integer D, let

After forming the numerator for each point, one may combine the equal denominators in (A.3):

Only finitely many positive norms occur. It is important that the recurrence precede this grouping. For fixed D, summing (A.1) over a norm shell gives

In general the first sum cannot be replaced by one value of r_1 times R_{m,N}(D). For example, (7,0) and (5,3) both have N=49, but their respective values of r_1 are 7\,529\,536 and -4\,912\,064. One may apply the recurrence either pointwise or within groups having the same pair (N,r_1), and only then combine all terms of a common norm.

For the finite coefficient enclosures, take D=256 when m=1 and D=48 when 2\leq m\leq38. The omitted dimensionless sum has absolute value at most

by (6.1). The coefficient c_m is consequently enclosed by adding [-E_{m,D},E_{m,D}] to the finite sum in (A.4) and then multiplying by the positive factor a^{-6m}/(6m). The square and its omitted tail are unchanged by grouping.

Complex pairs and interval containment. The finite calculation may instead keep both real components. Write

The multiplication rule

shows that, starting from (x_0,y_0)=(1,0),

represents \zeta^{6(m+1)} exactly. In particular x_m=r_m and its real component obeys (A.1). This is the complex-pair recursion used in the finite calculation, again applied per point before norm grouping. The half square n>0, or n=0 and q>0, contains one point of each pair (n,q),(-n,-q). Because 6m is even, the two real numerators agree, which justifies doubling each half-square contribution.

For completeness, grouping also explains an interval containment relation. At fixed denominator Q=2^{144}, let

For any finite rational list x_1,\ldots,x_h, the elementary inequalities

give

Apply this to the exact rationals in each norm shell. Converting the shell sum once and then summing shells is contained in the enclosure obtained by converting every point separately. Adding the same outward tail interval and multiplying by the same positive scaling enclosure preserve containment, because interval addition, multiplication, and outward rounding are monotone under inclusion. This statement assumes the same primitive rounding rule, tail enclosure, and scaling enclosure on both sides; exact algebra alone does not assert endpoint identity or containment between different choices of those enclosures.

A.2 Alternative exact bounds for the omitted modes

For the allowed indices k=234,240,\ldots, split the weighted summand (k+1)^4w_k of (6.3) into four positive rationals:

Then

This is the same \mathfrak a as in (6.6), distinct from the lattice scale a. If

in the order i=1,2,3,4, the exact consecutive ratios are

Both rational factors decrease for k>0, so each \rho_i(k) decreases. The four initial ratios, with no rounding, are

Here is an exact simultaneous check that all four are less than 1/10. The integer comparisons

give, since p_i\leq5 and 0<\beta_i\leq13/20,

Three geometric bounds. Set T_0=\sum_{i=1}^4t_i(234). One convenient entirely integer form of its exact comparison is obtained by putting

The common denominator 500^{234} gives

The displayed integer inequality is exactly T_0<156/10^{23}. In particular,

where the last comparison reduces to 1560<1800.

Let \rho_{\max}=\max_{1\leq i\leq4}\rho_i(234), and define three exact rational upper bounds

Decreasing ratios imply t_i(234+6j)\leq t_i(234)\rho_i(234)^j. Summing these geometric majorants and using (A.7) proves

The second inequality holds because 1-\rho_i(234)\geq1-\rho_{\max}>0. The common-ratio estimate keeps only the largest initial ratio, whereas the componentwise estimate retains each one; both refine the coarse factor 10/9.

An exact partial sum and a positive remainder bound. There are exactly 62 allowed indices from 234 through 600, because 600=234+6\cdot61. Define the finite rational sum

and the positive rational number

Every numerator is positive and 0<\rho_i(606)<\rho_i(234)<1/10, so R_{606}>0. Starting the geometric estimate at the next mode 606, rather than at the last included mode 600, gives

The first inequality is strict because all subsequent summands are positive. The number R_{606} is an upper bound for the remainder, not an assertion that the remainder equals it.

This refinement is no larger than the componentwise bound. For the i-th component put r=\rho_i(234). Then t_i(606)\leq t_i(234)r^{62} and \rho_i(606)\leq r, while each term of its partial sum is at most t_i(234)r^j. Consequently

Summing over i places P_{600}+R_{606} below A_{\rm component} in (A.9). These finite bounds are exact rationals; a rounded display that makes P_{600} and P_{600}+R_{606} look equal does not change (A.12).

Rational propagation to the radial errors. The refined bounds can be substituted directly into the propagation estimates of Section 6.5. They do not replace its analytic mode bounds or the finite premises (6.5): the retained second-order B_1 norm is still bounded by 6000, and the weighted third-order D_i norm by 2\cdot10^6. Recall also that

so replacing W_* by 4\mathfrak a is an upper bound, not an identity. Collecting the integral and mode-norm contributions in (6.8) and using \pi<22/7 gives the rational expression, for A,W\geq0,

The terms are respectively the source integral, the S_2 error, the linear D_i error, the retained and omitted mixed D_iH_{k_i} errors, and the quadratic D_iD_j error. The bound \pi<22/7 and (6.8) give |f_0-\widetilde f_0|_{2,t}\leq\mathcal F(\mathfrak a,W_*). Every coefficient of \mathcal F is nonnegative. Thus any one of the rational upper bounds A above, for which \mathfrak a\leq A, satisfies

Direct collection of the rational terms gives the exact identity

For example, the source term contributes 176A; the remaining linear terms contribute 56\,264\,004A, and the quadratic terms contribute 104A^2.

Put A_0=2/10^{21}. Each rational upper bound for \mathfrak a above is less than A_0, and the coefficient in (A.15) increases with A. The exact comparison

therefore supplies the claimed coefficient bound. More explicitly,

The order-two norm is |F|_{2,t}=|F(t)|+|F'(t)|+|F''(t)|/2. Hence the second derivative, unlike the value and first derivative, requires the factor 2:

For the angular function, (6.9) similarly yields the entirely rational comparison

The value and first derivative are bounded by their norms without this factor. Thus these rational alternatives preserve the error threshold 10^{-11} for all six quantities in (6.10) on [b,l]. For t\geq l, the exact continuation (6.11), rather than a uniform absolute error at unbounded radius, remains the applicable statement.