# The universal optimality of the triangular lattice

*An expository walkthrough, written to be readable with high-school mathematics, while exactly stating every theorem and every step of the proof.*

---

## How to read this

- For the big picture: Sections 0, 1, 2, 3.
- For the construction itself: Sections 4, 5, 6, 7.
- For turning the construction into energy: Sections 8, 9, 10.
- For the certification: Section 11.
- For a one-page summary: Section 12.6.
- For a map of every fact invoked: Appendix C.

## Contents

- **0. Prologue** -- the paper in one sentence.
- **1. Cast of characters** -- lattice, configurations, kernels, energy.
- **2. Main theorem** -- universal energy minimality; history; quantifiers; Riesz and log corollary.
- **3. Strategy** -- turn energy into Fourier signs; sharp pairs; roadmap.
- **4. Stage A: candidate functions** -- sine product, cardinal functions, Fourier transform, jet operator.
- **5. Stage B: finite data** -- truncation, Fejer quadrature, mass matrix, power sums, tables.
- **6. Stage C: exact interpolation** -- operator equation, Schur complement, correction bound.
- **7. Signs** -- deleted-slack quotient, Bernstein rows, quadratic barrier.
- **8. Gaussian to energy** -- density LP, sharp pairs, Poisson summation.
- **9. Gaussian to universal** -- positive mixtures, singular kernels, shifted potentials.
- **10. Renormalized Riesz and jellium** -- field energy, heat comparison, screening, thermodynamic limit.
- **11. Arithmetic verification** -- rational enclosures, scalar budgets, floor-rounding, Arb scripts.
- **12. Epilogue** -- summary, quantifier subtlety, open questions.
- **Appendix A** -- where the constants come from.
- **Appendix B** -- floor-rounded algorithm reference.
- **Appendix C** -- statement of every hypothesis.

---

## 0. Prologue: the paper in one sentence

Among all locally finite point sets in the plane with average density one, the triangular (hexagonal) lattice is the minimizer of average pairwise energy per point for every nonnegative completely monotone function of squared distance -- including when the energy is infinite. The same lattice also minimizes the logarithmic and Riesz renormalized energies with $0 \le s < 2$, and the corresponding jellium problems.

If that was too dense, here is a friendlier version. Pick a favorite pairwise interaction. Some examples:

- **Inverse power**: $g(t) = t^{-p}$, $p > 0$. Points close together repel strongly.
- **Gaussian**: $g(t) = e^{-\pi \alpha t}$, $\alpha > 0$. Repulsion dies off quickly.
- **Logarithm**: $g(t) = -\frac{1}{2} \log t$. The planar Coulomb potential, singular at $t = 0$ and slowly growing at infinity.
- **Riesz**: $g(t) = t^{-s/2}$, $0 < s < 2$. A power law interpolating between Coulomb-like and slower decay.

Now consider all ways to scatter identical points in the plane at one point per unit area. The paper proves that the triangular lattice minimizes the average pairwise energy per point for every smooth completely monotone kernel (all inverse powers, all Gaussians, and any positive mixture of these). It also proves the analogous statement for the log and Riesz kernels, which are too singular to compare directly and require a standard neutralization ("renormalization") step.

The word "universal" means: one configuration minimizes the whole family at once, not one kernel at a time. The Cohn-Kumar conjecture (2007) predicted this for the triangular lattice, the $E_8$ lattice in dimension eight, and the Leech lattice in dimension twenty-four. This paper settles the planar case.

---

## 1. The cast of characters

### 1.1 The triangular lattice

Let $b = \sqrt 3 / 2$. The triangular lattice $A$ of unit covolume is

$$A = b^{-1/2}\bigl\{ j\,(1,0) + k\,(1/2, b) : j, k \in \mathbb Z \bigr\}.$$

The prefactor $b^{-1/2}$ normalizes area: the parallelogram spanned by $b^{-1/2}(1,0)$ and $b^{-1/2}(1/2, b)$ has area exactly one. So $A$ is the triangular lattice at density one.

Two constants appear everywhere:

$$b = \frac{\sqrt 3}{2} = 0.8660254\ldots, \qquad B = \frac{1}{b^{2}} = \frac{4}{3}.$$

$B$ is the "area scale" that will reappear in every Gaussians Fourier-transform computation.

The **dual lattice** $A^{*}$ is the set of all $w \in \mathbb R^{2}$ such that $w \cdot a \in \mathbb Z$ for every $a \in A$. Since the basis of $A$ has determinant one, a quarter-turn $J$ gives $A^{*} = JA$. Because $A$ is rotation-invariant, $A$ and $A^{*}$ have the same set of $\lvert a \rvert^{2}$ values -- this will let us make contacts on the direct and dual side simultaneously.

### 1.2 Configurations, density, and disks

A **configuration** is a locally finite subset $C \subset \mathbb R^{2}$: infinite, but with only finitely many points in any bounded region.

For $R > 0$, let $B_{R}$ be the closed disk of radius $R$ centered at the origin, $C_{R} = C \cap B_{R}$, and $N_{R} = \# C_{R}$.

The configuration $C$ has **centered disk density one** if

$$\frac{N_{R}}{\pi R^{2}} \longrightarrow 1 \qquad (R \to \infty).$$

This is *weaker* than periodicity. Points can cluster locally, so long as the global average count per unit area is one, and there is no separation assumption at all.

### 1.3 Completely monotone kernels

Let $g : (0, \infty) \to [0, \infty)$ be smooth. Say that $g$ is **completely monotone** if

$$(-1)^{r} g^{(r)}(t) \ge 0 \quad\text{ for every integer } r \ge 0 \text{ and every } t > 0.$$

The derivatives alternate in sign, and every $(-1)^{r} g^{(r)}$ is nonnegative.

Examples:

- **Inverse powers**: $g(t) = t^{-p}$. Its $r$th derivative is $(-1)^{r} p(p+1) \cdots (p+r-1)\, t^{-p-r}$, and $(-1)^{r} g^{(r)}(t) = p(p+1) \cdots (p+r-1)\, t^{-p-r} > 0$.
- **Gaussians**: $g(t) = e^{-\alpha t}$. Its $r$th derivative is $(-\alpha)^{r} e^{-\alpha t}$, so $(-1)^{r} g^{(r)}(t) = \alpha^{r} e^{-\alpha t} \ge 0$.
- **Products**: $t^{-p} e^{-t}$ has the same sign pattern, by the Leibniz rule applied to nonnegative-magnitude derivatives.
- **Positive mixtures**: if each $g_{\alpha}$ is completely monotone and $\mu$ is a positive measure, then $\int g_{\alpha}(t) \, d\mu(\alpha)$ is completely monotone. This closure under mixtures is the engine of Section 9.

The classical **Hausdorff-Bernstein-Widder theorem** (proved in shifted form in Section 9.1) says the above is *equivalent* to the representation

$$g(t) = \int_{[0,1]} v^{t} \, d\rho(v)$$

for some positive Borel measure $\rho$ on $[0, 1]$. The change of variables $v = e^{-\alpha}$ turns this into a Laplace transform:

$$g(t) = \int_{0}^{\infty} e^{-\alpha t} \, d\tilde\rho(\alpha).$$

For $g = t^{-p}$, the representing measure has infinite mass at $\alpha = 0$: the decomposition is valid, but the total mass is not finite. This is exactly why the paper will *shift* the potential before mixing (Section 9.2).

### 1.4 The energy of a configuration

Fix a nonnegative kernel $g : (0, \infty) \to [0, \infty)$. The **average ordered-pair energy** at radius $R$ is

$$E_{g}(C_{R}) \;=\; \frac{1}{N_{R}} \sum_{\substack{x, y \in C_{R} \\ x \ne y}} g\bigl(\lvert x - y \rvert^{2}\bigr).$$

The sum counts *ordered* pairs: $(x, y)$ and $(y, x)$ are separate terms; $x = y$ is excluded.

The full energy is

$$E_{g}(C) \;=\; \liminf_{R \to \infty} E_{g}(C_{R}) \in [0, \infty].$$

We use $\liminf$ rather than $\lim$ because $E_{g}(C_{R})$ may oscillate with $R$; the lower limit captures the asymptotic "worst-case-lower" energy. Nonnegativity of $g$ ensures $E_{g}(C)$ is well-defined in the extended reals.

The lattice has a clean value:

$$E_{g}(A) = \sum_{a \in A \setminus \{0\}} g\bigl(\lvert a \rvert^{2}\bigr).$$

For $g = t^{-p}$ with $p \ge 1$ this sum is infinite, which is fine: the theorem's conclusion is that every competitor is also infinite.

---

## 2. The main theorem: universal energy minimality

**Theorem (Universal energy minimality).** *Let $g : (0, \infty) \to [0, \infty)$ be smooth and completely monotone. For every locally finite configuration $C \subset \mathbb R^{2}$ of centered disk density one,*

$$E_{g}(C) \;\ge\; E_{g}(A) \;=\; \sum_{a \in A \setminus \{0\}} g\bigl(\lvert a \rvert^{2}\bigr),$$

*where both sides are interpreted in the extended nonnegative reals $[0, \infty]$.*

Three things stand out.

- **Universal.** One configuration $A$ witnesses the inequality for the whole family of completely monotone $g$ at once.
- **Against everything.** The inequality is not restricted to lattices, or to periodic configurations, or to bounded density fluctuations; it applies to any locally finite $C$ with centered disk density one.
- **Infinite energies allowed.** When the right-hand side is infinite (e.g. $g = t^{-1}$), the inequality says $E_{g}(C) = \infty$ for every competitor. The formulation is meaningful throughout.

The theorem immediately covers every inverse power $g(t) = t^{-p}$, $p>0$, and every Gaussian $g(t) = e^{-\alpha t}$, $\alpha > 0$. It covers each positive mixture of these by linearity of the representing measure.

### 2.1 What is new, and what is not

Lattice-versus-lattice comparisons go back decades:

- **Rankin (1953)** and **Cassels (1959)** proved that $A$ minimizes the Epstein zeta sum $\sum_{a \ne 0} \lvert a \rvert^{-2p}$ over all unit-covolume planar lattices for $p > 1$. Cassels's proof had a gap that Ennola (1964) and Diananda (1964) closed.
- **Montgomery (1988)** proved that $A$ minimizes $\sum_{a \ne 0} e^{-\pi \alpha \lvert a \rvert^{2}}$ over all unit-covolume lattices, for every $\alpha > 0$.

The **Cohn-Kumar conjecture (2007)** asked for something strictly stronger: comparison not just against other lattices, but against every configuration, and one universal auxiliary proof technique that works for the whole completely monotone family.

**Cohn, Kumar, Miller, Radchenko, and Viazovska (2022)** proved the Cohn-Kumar conjecture for the $E_{8}$ and Leech lattices, constructing sharp Fourier interpolation formulas. Their methods inspired the present paper.

Previous planar work covered restricted competitor classes:

- Faulhuber-Shafkulovska-Zlotnikov (2023): $A$ vs. the honeycomb configuration, and vs. products / unions of translates of rectangular or triangular lattices.
- Hardin-Tenpas (2025): local optimality among configurations that are four or six points per fixed triangular lattice.
- Leblé (2025): $A$ cannot be lowered by sufficiently small bounded displacements, with the allowed displacement size depending on the kernel.

None of these reaches the full "all locally finite configurations of centered disk density one" conclusion proved in Theorem 2.

### 2.2 The subtlety about quantifiers

The Cohn-Kumar conjecture asked for a single auxiliary function $f$, allowed to depend on the kernel $g$, with certain properties (see Section 3.3). This is strictly stronger than what the present paper establishes.

The present paper constructs a sharp auxiliary function $f_{\alpha}$ for *each Gaussian scale* $\alpha > 0$. For a general mixed potential (e.g. $t^{-1/2}$), it does not construct a single auxiliary function; it only proves the *energy* bound by integrating the Gaussian-scale bounds over the positive-mixture representation. The paper states this honestly: the stronger per-potential quantifier of Cohn-Kumar's original conjectured form is not claimed.

The energy conclusion (the theorem) is fully established regardless.

### 2.3 Corollary: Riesz and log

For $0 < s < 2$, the raw sum $\sum_{a \ne 0} \lvert a \rvert^{-s}$ diverges. So the "energy comparison" above is trivially $\infty \ge \infty$ and carries no information. Physically this is because an infinite configuration of all-same-sign charges has infinite energy: standard Coulomb physics neutralizes the charges against a uniform oppositely charged background.

Two renormalizations are standard.

- **Periodic.** Fix a torus $T = \mathbb R^{2} / L$. Neutralize the points against a uniform background of opposite total charge on $T$. Remove each point's self-energy. The result is a finite periodic energy $W_{s, L}$.
- **Field.** Build the infinite "jellium" configuration: point charges plus uniform background. Introduce a self-energy cutoff $\eta > 0$, remove it, and take the field-energy per unit area. Minimize over compatible gradient fields. The result is the infinite-system energy $W_{s}$.

The paper proves the following corollary (Section 10.4).

**Corollary (planar renormalized and jellium energies).** *Let $0 \le s < 2$, with interaction $g(x) = -\log \lvert x \rvert$ at $s = 0$ and $g(x) = \lvert x \rvert^{-s}$ for $0 < s < 2$.*

*(i) For each $n \ge 1$, the configuration $A/nA$ (i.e. the $n^{2}$ distinct classes of the triangular lattice in the torus $\mathbb R^{2}/(nA)$) minimizes $W_{s, nA}$ among all configurations of $n^{2}$ distinct points in that torus.*

*(ii) The infinite-system energy $W_{s}(A)$ minimizes $W_{s}$ over all density-one planar configurations. In addition, the ordered-pair thermodynamic jellium minimum equals $W_{s}(A)/\kappa_{s}$, where $\kappa_{0} = 2\pi$ and $\kappa_{s} = 4\pi$ for $0 < s < 2$.*

The proof routes through the Gaussian comparison proved in the main theorem, by way of a "heat-kernel representation" (Petrache-Serfaty, 2020) that turns the renormalized energy into an integral of Gaussian energies over a heat-time parameter (Section 10.2).

---

## 3. The road: from energies to Fourier signs

The theorem has a simple shape after Fourier transformation. We describe the shape now; every piece gets its own section afterward.

### 3.1 Fourier transform and nonnegativity

For a Schwartz function $f : \mathbb R^{2} \to \mathbb R$ (smooth, rapidly decaying, and so is every derivative), define

$$\widehat f(\xi) \;=\; \int_{\mathbb R^{2}} f(x)\, e^{-2\pi i\, x \cdot \xi}\, dx.$$

For finite signed measures $\mu, \nu$ on $\mathbb R^{2}$ set

$$Q_{f}(\mu, \nu) = \iint f(x - y)\, d\mu(x)\, d\nu(y).$$

Fourier inversion and Fubini give

$$Q_{f}(\mu, \nu) = \int_{\mathbb R^{2}} \widehat f(\xi)\, \overline{\widehat\mu(\xi)}\, \widehat\nu(\xi)\, d\xi.$$

If $\widehat f \ge 0$ pointwise, then $Q_{f}$ is positive semi-definite: $Q_{f}(\mu, \mu) \ge 0$ for every finite measure $\mu$. So Cauchy-Schwarz applies to $Q_{f}$.

For a configuration $C$ let $\mu_{R} = \sum_{x \in C_{R}} \delta_{x}$. Then $Q_{f}(\mu_{R}, \mu_{R}) = \sum_{x, y \in C_{R}} f(x - y)$. The diagonal terms contribute $N_{R} \cdot f(0)$, and the off-diagonal terms are the pairwise energy. Hence

$$\sum_{\substack{x, y \in C_{R} \\ x \ne y}} f(x - y) = Q_{f}(\mu_{R}, \mu_{R}) - N_{R}\, f(0).$$

When $\widehat f \ge 0$ we get $Q_{f}(\mu_{R}, \mu_{R}) \ge 0$ and therefore

$$\sum_{\substack{x, y \in C_{R} \\ x \ne y}} f(x - y) \ge -N_{R}\, f(0).$$

Dividing by $N_{R}$ and taking the lower limit gives the trivial lower bound

$$\liminf_{R \to \infty} \frac{1}{N_{R}} \sum_{\substack{x, y \in C_{R} \\ x \ne y}} f(x - y) \ge -f(0).$$

This is non-vacuous when $\widehat f(0) - f(0) > 0$, but it neglects density. Density one sharpens the bound below.

### 3.2 The Cohn-Zhao linear program

Density one enters by comparing $\mu_{R}$ to the area measure $\nu_{R} = \mathbf 1_{B_{(1+\epsilon)R}}$ (uniform area measure with total $V_{R} = \pi (1+\epsilon)^{2} R^{2}$). This trick is due to Cohn-Zhao (2014) and was made explicit for packings by Cohn-de Courcy-Ireland (2018).

Fix $\epsilon > 0$. For every $x \in B_{R}$ the disk $B_{(1+\epsilon)R} - x$ contains $B_{\epsilon R}$. So

$$Q_{f}(\mu_{R}, \nu_{R}) = N_{R}\, \widehat f(0) + O\!\bigl(N_{R} \cdot \int_{\lvert z \rvert \ge \epsilon R} \lvert f(z) \rvert\, dz\bigr),$$

and the error term is *uniform* in the positions of the points of $C_{R}$: clustering does not help or hurt. Simultaneously $Q_{f}(\nu_{R}, \nu_{R}) / V_{R} \to \widehat f(0)$.

Cauchy-Schwarz gives $Q_{f}(\mu_{R}, \mu_{R})\, Q_{f}(\nu_{R}, \nu_{R}) \ge \lvert Q_{f}(\mu_{R}, \nu_{R}) \rvert^{2}$. Divide by $N_{R}^{2}$:

$$\frac{Q_{f}(\mu_{R}, \mu_{R})}{N_{R}} \cdot \frac{V_{R}}{N_{R}} \cdot \frac{Q_{f}(\nu_{R}, \nu_{R})}{V_{R}} \ge \left( \frac{Q_{f}(\mu_{R}, \nu_{R})}{N_{R}} \right)^{2}.$$

As $R \to \infty$: $V_{R} / N_{R} \to (1+\epsilon)^{2}$, and the two ends converge to $\widehat f(0)^{2}$ each. Taking lower limits in $R$ then $\epsilon \downarrow 0$ yields

$$\liminf_{R \to \infty} \frac{Q_{f}(\mu_{R}, \mu_{R})}{N_{R}} \ge \widehat f(0).$$

Removing the $N_{R}$ diagonal terms gives the **density LP inequality**:

$$\liminf_{R \to \infty} \frac{1}{N_{R}} \sum_{\substack{x, y \in C_{R} \\ x \ne y}} f(x - y) \;\ge\; \widehat f(0) - f(0).$$

If $f \le g$ pointwise off the origin and $\widehat f \ge 0$, then the same argument runs with $f$ in place of $g$ inside the sum (termwise domination), giving

$$E_{g}(C) \;\ge\; \widehat f(0) - f(0).$$

**Key message.** A *minorant* $f \le g$ with nonnegative Fourier transform produces a lower bound on the average $g$-energy.

### 3.3 The Gaussian "sharp pair"

The main theorem will be proved in the Gaussian case, $g(t) = e^{-\pi \alpha t}$, by constructing a suitable minorant at each $\alpha$. Define $G_{\alpha}(x) = e^{-\pi \alpha \lvert x \rvert^{2}}$. We want a real radial Schwartz function $f_{\alpha}$ with

$$f_{\alpha} \le G_{\alpha}, \qquad \widehat f_{\alpha} \ge 0, \qquad (\text{two sign conditions})$$

and two families of *contacts*:

$$f_{\alpha}(a) = G_{\alpha}(a) \text{ for every } a \in A \setminus \{0\}, \qquad (\text{direct contacts})$$

$$\widehat f_{\alpha}(w) = 0 \text{ for every } w \in A^{*} \setminus \{0\}. \qquad (\text{dual contacts})$$

**Why the contacts produce the sharp estimate.**

Poisson summation for the lattice $A$ says, for any decaying $f$,

$$\sum_{a \in A} f(a) = \sum_{w \in A^{*}} \widehat f(w).$$

Split off the origin of each side:

$$f(0) + \sum_{a \ne 0} f(a) = \widehat f(0) + \sum_{w \ne 0} \widehat f(w).$$

With $f = f_{\alpha}$, the dual contacts give $\sum_{w \ne 0} \widehat f_{\alpha}(w) = 0$, so

$$\sum_{a \ne 0} f_{\alpha}(a) = \widehat f_{\alpha}(0) - f_{\alpha}(0).$$

The direct contacts set $f_{\alpha}(a) = G_{\alpha}(a)$ at every nonzero lattice point, so

$$\sum_{a \ne 0} G_{\alpha}(a) = \widehat f_{\alpha}(0) - f_{\alpha}(0).$$

Combining this identity with the density LP inequality for $f_{\alpha}$ against $G_{\alpha}$ yields

$$E_{e^{-\pi\alpha t}}(C) \ge \widehat{f_{\alpha}}(0) - f_{\alpha}(0) = \sum_{a \ne 0} G_{\alpha}(a) = E_{e^{-\pi\alpha t}}(A).$$

That is the Gaussian case of the main theorem. The remaining sections construct the pair $(f_{\alpha}, \widehat f_{\alpha})$ for every $\alpha > 0$.

### 3.4 What is hard about the construction

The two sign conditions and the two contact conditions must hold *simultaneously*, *at every scale* $\alpha$, and at every radius.

- Contact conditions are *interpolation* conditions. There are infinitely many of them, one at every nonzero lattice shell.
- Sign conditions are *continuous* everywhere.

Neither condition alone is hard. Simultaneously, they multiply the difficulty.

A further subtlety: values and first derivatives of $f$ and $\widehat f$ on the actual shells of $A$ do *not* suffice to determine a radial Schwartz function. Sardari (2021) proved that infinitely many linearly independent radial Schwartz functions have all those data zero. So the interpolation problem must be enlarged: we work with a superset of shells (the periodic node set $\mathcal N$ of Section 4.1), and impose contacts at every node of the superset. In a summable coefficient space this *does* determine the function uniquely (Section 6).

### 3.5 The three-stage strategy

The construction splits into three stages, each handled in its own section.

**Stage A (Section 4).** Set up the *machinery*: express every admissible interpolation datum (value and first derivative at every node) as an entire function of a complex variable, Fourier-transform it, and prove the resulting Schwartz function pair has explicit, uniform bounds. In this stage there is no "interpolation" yet; we just write down the *form* of every candidate.

**Stage B (Sections 5-6).** Solve the *interpolation problem* exactly. Turn the coupled set of "value-and-slope at every node" equations into a linear operator equation in a summable space. Approximate it by an explicit finite matrix $D$ (Stage B1, Section 5), invert it there, and solve the remaining infinite tail by a Schur-complement argument (Stage B2, Section 6). The crucial small quantity is that the operator's outputs at far nodes are *tiny*: interpolation with distant inputs never significantly reinforces a given node.

**Stage C (Section 7).** Prove the *two signs*. This turns the contact data from Stage B into a genuine pair of upper-lower bounds. The strategy is:

1. Subtract value and first derivative at every node, so the residues are quotient functions with bounded behavior near nodes.
2. Bound the quotient from below on each *half-gap* between adjacent nodes using *Bernstein polynomial* coefficient estimates (Section 5.3 provides the input, Section 7.2 applies it).
3. Bound the quotient from below on the *unbounded tail* using a *quadratic barrier* for the sine product (Section 7.4).

The final sections (8-10) convert the Gaussian pair into the universal mixture inequality and then into the renormalized/jellium corollary. Two appendices supply the numerical inequalities and the alternative rational verification procedure used in Stages B and C.

The next sections build Stage A concretely.

---

## 4. Stage A: building entire candidate functions

We now construct the entire infinite family of "candidate" functions that the interpolation problem of Sections 5-6 will pin down to a specific pair. These candidates are of the shape "sine product times a rational sum." Their point is: each candidate has control over value and first derivative at every node of the periodic set $\mathcal N$, and the map from "dataset" to "function" is *linear* and *summable*.

### 4.1 Shells, nodes, and the periodic set $\mathcal N$

The construction uses four numerical parameters:

$$b = \frac{\sqrt 3}{2}, \qquad h = \frac{2}{5}, \qquad B = \frac{4}{3} = \frac{1}{b^{2}}, \qquad s = b \lvert x \rvert^{2}.$$

$s$ is the "rescaled squared radius" used throughout Sections 4-7. $h = 2/5$ is a Gaussian damping exponent chosen to make numerical bounds tight; other values would work but produce worse constants.

Take the basis

$$u = \sqrt{2/\sqrt{3}}\;(1, 0), \qquad v = \sqrt{2/\sqrt{3}}\;(1/2, \,b).$$

The determinant of $(u, v)$ is one, so this is a unit-covolume lattice -- exactly the triangular lattice $A$ from Section 1.1.

For any $j, \ell \in \mathbb Z$:

$$b \lvert j u + \ell v \rvert^{2} = j^{2} + j \ell + \ell^{2}.$$

This quadratic form takes the values $0, 1, 3, 4, 7, 9, 12, 13, 16, 21, \ldots$ -- the *shells* of $A$.

Two local computations pin down these shells modulo $12$.

- **Modulo 3.** $j^{2} + j \ell + \ell^{2} \equiv (j - \ell)^{2} \pmod 3$, so its residue is $0$ or $1$ in $\{0, 1, 2\}$.
- **Modulo 2.** $j^{2} + j \ell + \ell^{2}$ is even iff both $j$ and $\ell$ are even, in which case it is divisible by 4.

The two restrictions, combined, force

$$j^{2} + j \ell + \ell^{2} \in \{\,0, 1, 3, 4, 7, 9\,\} \pmod{12}.$$

Define

$$L = \{0, 1, 3, 4, 7, 9\}, \qquad \mathcal Z = 12\mathbb Z + L, \qquad \mathcal N = \mathcal Z \cap (0, \infty), \qquad \mathcal N_{0} = \{0\} \cup \mathcal N.$$

**Claim.** Every nonzero shell of $A$ lies in $\mathcal N$.

**Claim (strictness).** The inclusion is proper: for example $24 \in \mathcal N$, but $24$ is not a shell of $A$. Indeed, if $j^{2} + j\ell + \ell^{2} = 24$, parity forces $j, \ell$ even, so $(j/2)^{2} + (j/2)(\ell/2) + (\ell/2)^{2} = 6$; parities force $j/2, \ell/2$ even too, so the original value would be divisible by 16, contradicting $16 \nmid 24$.

The extra nodes (like 24) are *artificial*: they are not shell radii of $A$. But they will also carry interpolation data. This is exactly the enlargement that the Sardari non-uniqueness theorem requires (Section 3.4).

Finally, recall the dual lattice $A^{*} = JA$ from Section 1.1 --- the same shells. So contacts at the nodes of $\mathcal N$ serve *both* the direct ($A$) side and the dual ($A^{*}$) side simultaneously.

### 4.2 The normalized Gaussian pair

The construction will produce *two* entire functions $H_{1}, H_{2}$ and a single radial Schwartz $f$. Everything is normalized to the parameter $k > 0$.

For each $k \ge 2.36$ define

$$T_{k}(s) = k^{-1} e^{-k(s - 1)}.$$

$T_{k}$ is a stretched exponential "reference profile," scaled so $T_{k}(1) = 1/k$ and decaying away from $s = 1$.

**Theorem (Normalized Gaussian Fourier pair).** *For every finite real $k \ge 2.36$, there exist entire functions $H_{1}, H_{2}$ real on the real axis, and a real radial Schwartz function $f$ on $\mathbb R^{2}$, with*

$$f(x) = e^{-\pi h b \lvert x \rvert^{2}}\, H_{1}(b \lvert x \rvert^{2}), \qquad \widehat f(\xi) = e^{-\pi h b \lvert \xi \rvert^{2}}\, H_{2}(b \lvert \xi \rvert^{2}),$$

*such that for every real $s \ge 0$:*

$$H_{1}(s) \le T_{k}(s), \qquad H_{2}(s) \ge 0, \qquad (\text{two signs})$$

*and at every node $n \in \mathcal N$:*

$$H_{1}(n) = T_{k}(n), \qquad H_{1}'(n) = T_{k}'(n) = -e^{-k(n-1)}, \qquad H_{2}(n) = H_{2}'(n) = 0. \qquad (\text{four contacts per node})$$

This is the technical heart of the paper. Everything else unfolds from it.

**Why a pair, and why this normalization.** The two functions $H_{1}, H_{2}$ will eventually form the two halves of $f$ and $\widehat f$. The Gaussian damping $e^{-\pi h b \lvert x \rvert^{2}}$ ensures $f, \widehat f$ are Schwartz; the $H_{i}$ are *entire* and control the jets at nodes. The two constants $b$ and $B = 1/b^{2}$ are dictated by the Gaussian-transform algebra: the transform of a Gaussian with parameter $\pi h b$ is another Gaussian with the same parameter, and the coupling constants $B, h$ inside $H_{i}$ are chosen to make this work out.

**Converting to the Gaussian of the main theorem.** If we take $k = \pi(\alpha / b - h)$, a straightforward algebra computation gives

$$k e^{-k} e^{-\pi h s} T_{k}(s) = e^{-(\pi h + k) s} = e^{-\pi \alpha \lvert x \rvert^{2}} = G_{\alpha}(x).$$

So $k e^{-k} \cdot f$ is a minorant of $G_{\alpha}$ (with $f \le e^{-\pi h s} T_{k}$ from the theorem), its Fourier is nonnegative (multiplied by $k e^{-k} > 0$, retains nonnegativity), and its direct and dual contacts at $\mathcal N \supseteq A \setminus \{0\}$ match $G_{\alpha}$ and zero respectively. This is the Gaussian pair of Section 3.3, for every $\alpha \ge 1$.

For $0 < \alpha < 1$ (i.e. $k$ too small), a Fourier-complement trick handles it (Section 8.2). We focus on $\alpha \ge 1$ for the construction.

**What "existence" means here.** The theorem asserts *there exists* a pair $(H_{1}, H_{2}, f)$. It does not claim uniqueness --- in fact the space of such pairs, in a fixed coefficient system, is one-dimensional (Section 6 will pin it down by the *specific* structure we build). The construction is explicit; the numerical inequalities are certified separately.

### 4.3 The sine product $P(s)$

The first ingredient is a *periodic* function with a double zero at every node and nowhere else. Define

$$P(s) = \prod_{a \in L} \left( 2 \sin\frac{\pi(s - a)}{12} \right)^{2}.$$

Six squared sine factors, one for each residue class of $\mathcal N$. Each factor is $2$-periodic in the appropriate variable (period 12 in $s$); the product is nonnegative, periodic with period 12, has a double zero at every point of the periodic set $\mathcal Z$, and no other real zeros.

$P$ has a short complex-exponential expansion:

$$P(s) = \sum_{\ell = -6}^{6} P_{\ell}\, e^{i \pi t_{\ell} s}, \qquad t_{\ell} = \frac{\ell}{6}, \qquad P_{-\ell} = \overline{P_{\ell}}.$$

**Coefficient values.** Setting $w = e^{i \pi s / 6}$ converts the product to $w^{-6} \prod_{a \in L} (w - e^{i \pi a / 6})^{2}$. The scalar factor is $1$ because $\sum_{a \in L} a = 24$ and the leading coefficient picks up exactly $1$. Multiplying the degree-12 polynomial gives

$$(P_{0}, P_{1}, P_{2}, P_{3}, P_{4}, P_{5}, P_{6}) = \left( 5, \; -1 + 2 i b, \; 1 + 2 i b, \; -2, \; -\tfrac{1}{2} + i b, \; -1 - 2 i b, \; 1 \right).$$

The absolute values that matter later:

$$\left( \lvert P_{1} \rvert, \ldots, \lvert P_{6} \rvert \right) = (2, 2, 2, 1, 2, 1), \qquad \sum_{\ell = 1}^{6} \lvert P_{\ell} \rvert = 10, \qquad \sum_{\ell = 1}^{6} t_{\ell} \lvert P_{\ell} \rvert = \frac{16}{3}.$$

**Node constants $Q_{n}, D_{n}$.** At any node $n \in \mathcal N_{0}$, write

$$P(n + u) = Q_{n} u^{2} \bigl( 1 + D_{n} u + O(u^{2}) \bigr).$$

The periodicity of $P$ makes $Q_{n}, D_{n}$ depend only on the residue of $n$ mod 12:

| $n \bmod 12$ | $Q_{n}$ | $D_{n}$ |
|---|---|---|
| 0 | $\pi^{2}/3$ | $-7 \pi \sqrt 3 / 18$ |
| 1 | $(2\pi^{2}/3)(2 - \sqrt 3)$ | $\pi(3 + \sqrt 3)/18$ |
| 3 | $(2\pi^{2}/3)(2 - \sqrt 3)$ | $-\pi(3 + \sqrt 3)/18$ |
| 4 | $\pi^{2}/3$ | $7 \pi \sqrt 3 / 18$ |
| 7 | $(2\pi^{2}/3)(2 + \sqrt 3)$ | $-\pi(3 - \sqrt 3)/18$ |
| 9 | $(2\pi^{2}/3)(2 + \sqrt 3)$ | $\pi(3 - \sqrt 3)/18$ |

To get the table: pull out one vanishing $\sin^{2}$ factor, take the value and log-derivative of the remaining product. Concretely,

$$Q_{a} = \left( \frac{\pi}{6} \right)^{2} \prod_{d \in L \setminus \{a\}} \left( 2 \sin\frac{\pi(a - d)}{12} \right)^{2}, \qquad D_{a} = \frac{\pi}{6} \sum_{d \in L \setminus \{a\}} \cot\frac{\pi(a - d)}{12}.$$

Two bounds (proved by rational comparisons of $\pi, \sqrt 3$; see Appendix A):

$$Q_{n} > \frac{7}{4} = 1.75, \qquad \lvert D_{n} \rvert < \frac{53}{25} = 2.12.$$

These constants show up again in every jet bound and every sign estimate. The specific rational witnesses are $\frac{157}{50} < \pi < \frac{22}{7}$ and $\sqrt 3 < \frac{1733}{1000}$.

### 4.4 The entire cardinal functions

Now we turn "specify a value and a slope at every node" into a single entire function. The trick: multiply the sine product by a sum of *double poles* (to control values) and *simple poles* (to control slopes).

**Definition (coefficient list).** A coefficient list is data of the form

$$\mathbf{a} = \left( C ; \; \{(c_{n}, d_{n})\}_{n \in \mathcal N_{0}} \right),$$

where $C \in \mathbb R$ is a constant, and $c_{n}, d_{n} \in \mathbb R$ for each node $n$. Its *unweighted* norm is

$$\lVert \mathbf{a} \rVert_{1} = \lvert C \rvert + \sum_{n \in \mathcal N_{0}} \left( \lvert c_{n} \rvert + \lvert d_{n} \rvert \right).$$

Summability ($\lVert \mathbf{a} \rVert_{1} < \infty$) is the fixed admissibility criterion.

Define

$$p_{\mathbf{a}}(s) = P(s) \left( C + \sum_{n \in \mathcal N_{0}} \left( \frac{c_{n}}{(s - n)^{2}} + \frac{d_{n}}{s - n} \right) \right), \quad s \in \mathbb C \setminus \mathcal N_{0}.$$

The double poles $c_{n}/(s-n)^{2}$ interact with the double zero of $P$ at $n$ to give a *finite* value $Q_{n} c_{n}$ at $n$. The simple poles give a finite derivative. Away from nodes, $p_{\mathbf{a}}$ is analytic.

**Lemma (Cardinal measure).** *For every summable $\mathbf{a}$:*

*(a) $p_{\mathbf{a}}$ extends to an entire function, real on the real axis, with jets*

$$p_{\mathbf{a}}(n) = Q_{n} c_{n}, \qquad p_{\mathbf{a}}'(n) = Q_{n}\left( d_{n} + D_{n} c_{n} \right) \qquad (n \in \mathcal N_{0}).$$

*(b) There is a finite complex measure $\mu_{\mathbf{a}}$ on $[-1, 1]$ such that*

$$p_{\mathbf{a}}(s) = \int_{-1}^{1} e^{i \pi t s} \, d\mu_{\mathbf{a}}(t) \qquad (s \in \mathbb C),$$

*and $\mu_{\mathbf{a}}(- E) = \overline{\mu_{\mathbf{a}}(E)}$ for every Borel set $E$, plus the TV bound*

$$\lVert \mu_{\mathbf{a}} \rVert_{\mathrm{TV}} \le 25 \lvert C \rvert + 36 \sum_{n \in \mathcal N_{0}} \left( \lvert c_{n} \rvert + \lvert d_{n} \rvert \right).$$

*(c) The same $p_{\mathbf{a}}$ has a "folded" evaluation on real $s$:*

$$p_{\mathbf{a}}(s) = \operatorname{Re} \int_{0}^{1} e^{i \pi t s} \, d\nu_{\mathbf{a}}(t) \qquad (s \in \mathbb R),$$

*where $\nu_{\mathbf{a}}$ is a finite complex measure supported on $[0, 1]$, with atoms $C e_{j} P_{j}$ at $t_{j} = j/6$ (with $e_{0} = 1$, $e_{j} = 2$ for $1 \le j \le 6$) and a density on each open interval $((j-1)/6, j/6)$ given by*

$$d_{j}(t) = 2\pi \sum_{n \in \mathcal N_{0}} e^{-i \pi t n} \sum_{\ell = j}^{6} \left[ -\pi c_{n}(t_{\ell} - t) + i d_{n} \right] P_{\ell}\, e^{i \pi t_{\ell} n}.$$

**Proof sketch.** For one node $n$, use $u = s - n$ and the elementary integral identities

$$\frac{e^{i \pi t_{\ell} u} - 1}{u} = i \pi \int_{0}^{t_{\ell}} e^{i \pi t u}\, dt, \qquad \frac{e^{i \pi t_{\ell} u} - 1 - i \pi t_{\ell} u}{u^{2}} = -\pi^{2} \int_{0}^{t_{\ell}} (t_{\ell} - t)\, e^{i \pi t u}\, dt.$$

Multiply the first by $P_{\ell} e^{i \pi t_{\ell} n}$ and sum over $\ell$; the subtracted constant and linear terms cancel because $P(n) = P'(n) = 0$. The result is an integral representation of $c_{n}/(s-n)^{2}$ and $d_{n}/(s-n)$ as spectral measures. The negative-frequency parts are conjugate reflections, so the "folding" gives item (c). Summability over $n$ gives (b). $\square$

**Why these particular constants 25 and 36.** The constant column $C P$ has variation $\lvert P_{0} \rvert + 2 \sum_{\ell = 1}^{6} \lvert P_{\ell} \rvert = 5 + 2 \cdot 10 = 25$. The $c_{n}$ column has variation $\pi^{2} \sum_{\ell=1}^{6} t_{\ell}^{2} \lvert P_{\ell} \rvert = 65 \pi^{2}/18 < 36$. The $d_{n}$ column has variation $2\pi \sum_{\ell = 1}^{6} t_{\ell} \lvert P_{\ell} \rvert = 32\pi/3 < 36$. So $25$ and $36$ are just convenient upper bounds.

### 4.5 The Gaussian transform of a cardinal function

The whole point of the sine product coefficients $P_{\ell}$ is that they pair beautifully with complex Gaussian transforms. Define, for $t \in [-1, 1]$,
$$z(t) = -\frac{B}{t + i h} - i h, \qquad \lambda(t) = \frac{i}{b(t + i h)}.$$
Both are holomorphic in $t$ except at $t = -ih$ (which is not on $[-1, 1]$; $h = 2/5$, so the pole is below the contour).

**Lemma (Transformed cardinal).** For every summable $\mathbf{a}$, the *transformed cardinal function*
$$K_{p_{\mathbf{a}}}(s) = \int_{-1}^{1} \lambda(t)\, e^{i \pi z(t) s}\, d\mu_{\mathbf{a}}(t)$$
is entire, and satisfies the reality condition $K_{p_{\mathbf{a}}}(\overline{s}) = \overline{K_{p_{\mathbf{a}}}(s)}$ for all complex $s$. In particular it is real on the real axis.

**Why the formula works.** The complex Gaussian transform in dimension two, with our convention $\widehat{f}(\xi) = \int f(x) e^{-2\pi i x \cdot \xi} dx$, is
$$\widehat{e^{-\pi \gamma \lvert x \rvert^{2}}}(\xi) = \gamma^{-1}\, e^{-\pi \lvert \xi \rvert^{2}/\gamma} \qquad (\operatorname{Re} \gamma > 0).$$
Apply it to the Gaussian $e^{-\pi h b \lvert x \rvert^{2}} \cdot e^{i \pi t s}$ with $s = b \lvert x \rvert^{2}$. The result is another Gaussian with parameter $\gamma = b(h - it)$, whose transform is another Gaussian with parameter $1/\gamma$. Writing $\gamma = b(h - it) = -ib(t+ih)$ and simplifying, the transform of $e^{-\pi h b \lvert x \rvert^{2}} e^{i \pi t s}$ equals
$$e^{-\pi h s}\, \lambda(t)\, e^{i \pi z(t) s}.$$
Integrating against $d\mu_{\mathbf{a}}(t)$ and swapping integrals gives the claimed $K_{p_{\mathbf{a}}}$.

**Pointwise frequency consistency.** $z$ is an involution on the upper half-plane: $z(z(t)) = t$. And $\lambda(t) \lambda(z(t)) = 1$. These identities hold pointwise; they do not define a composition of the folded integral because $z$ does not preserve the integration contour $[-1, 1]$. They are also not needed for anything in the argument except as sanity checks; the Fourier-pair conclusion below is proved from the actual complex Gaussian formula, not from these.

**Lemma (Fourier pairs).** For every summable $\mathbf{a}$, the radial functions $e^{-\pi h s} p_{\mathbf{a}}(s)$ and $e^{-\pi h s} K_{p_{\mathbf{a}}}(s)$, viewed as functions on $\mathbb R^{2}$ via $s = b \lvert x \rvert^{2}$, are real Schwartz functions, and are Fourier transforms of one another:
$$\widehat{\bigl(e^{-\pi h s}\, p_{\mathbf{a}}\bigr)} = e^{-\pi h s}\, K_{p_{\mathbf{a}}}, \qquad \widehat{\bigl(e^{-\pi h s}\, K_{p_{\mathbf{a}}}\bigr)} = e^{-\pi h s}\, p_{\mathbf{a}}.$$

**Lemma (Coupled pair).** Suppose $\mathbf{a}_{1}, \mathbf{a}_{2}$ are summable. Define
$$H_{1}(s) = p_{\mathbf{a}_{1}}(s) + K_{p_{\mathbf{a}_{2}}}(s), \qquad H_{2}(s) = p_{\mathbf{a}_{2}}(s) + K_{p_{\mathbf{a}_{1}}}(s).$$
Then the radial function $f(x) = e^{-\pi h b \lvert x \rvert^{2}} H_{1}(b \lvert x \rvert^{2})$ has Fourier transform $\widehat{f}(\xi) = e^{-\pi h b \lvert \xi \rvert^{2}} H_{2}(b \lvert \xi \rvert^{2})$.

Both lemmas are proved in the paper by direct use of the complex Gaussian transform and the exchange of Fourier integration with the finite measure integration. The details: measure is finite, frequency is on a compact interval, every derivative of a Gaussian is bounded by a polynomial times a decaying Gaussian. So the exchanges are legitimate.

**Why two lists.** The point of the coupled form is symmetry on the direct/dual sides. A single list would produce $f$ and $\widehat f$ tied together by one function; two lists decouple them and let us impose independent interpolation conditions on $H_{1}$ and $H_{2}$ separately. This is the difference between the Cohn-Kumar scalar setting and the present coupled setting.

### 4.6 The jet operator and the distant-output bound

Everything we care about is *local data at nodes*: value and first derivative. Define the **jet projection**
$$(\mathcal{J} u)_{c, n} = \frac{u(n)}{Q_{n}}, \qquad (\mathcal{J} u)_{d, n} = \frac{u'(n) - D_{n} u(n)}{Q_{n}} \qquad (n \in \mathcal N).$$
The normalizations $Q_{n}, D_{n}$ undo the leading behavior of $P$ at that node. With these factors,
$$\mathcal{J} p_{\mathbf{a}} = \mathbf{a}^{+}$$
where $\mathbf{a}^{+}$ denotes the positive-node part of $\mathbf{a}$, i.e. $\{c_{n}, d_{n}\}_{n \in \mathcal N}$ with the "extra coordinates" $(c_{0}, d_{0}, C)$ suppressed.

Define the **transformed jet operator**
$$S \mathbf{a} = -\mathcal{J} K_{p_{\mathbf{a}}}.$$
Then the coupled pair $H_{1}, H_{2}$ satisfies
$$\mathcal{J} H_{1} = \mathbf{a}_{1}^{+} - S \mathbf{a}_{2}, \qquad \mathcal{J} H_{2} = \mathbf{a}_{2}^{+} - S \mathbf{a}_{1}.$$

So the interpolation problem --- prescribe $\mathcal{J} H_{1} = \mathcal{J} T_{k}$ and $\mathcal{J} H_{2} = 0$ at every positive node --- becomes the *linear system*
$$\mathbf{a}_{1}^{+} - S \mathbf{a}_{2} = \mathcal{J} T_{k}, \qquad \mathbf{a}_{2}^{+} - S \mathbf{a}_{1} = 0.$$

Whether this system has a solution depends on how small the off-diagonal operator $S$ is. The construction is designed so $S$ is not small everywhere, but is exponentially small *at distant nodes*, which is enough to fit into a fixed-point framework.

**Lemma (Kernel bounds).** For every $t \in [0, 1]$:
$$\lvert \lambda(t) \rvert < 3, \qquad \lvert z(t) \rvert \le \frac{44}{15} < 3, \qquad \pi \operatorname{Im} z(t) > 0.187.$$
On $0 \le t \le 5/6$ the last improves to $\pi \operatorname{Im} z(t) > 0.70$. On $5/6 \le t \le 1$ the bounds are sharper:
$$\lvert \lambda(t) \rvert < 1.25, \qquad \lvert z(t) \rvert < 1.32, \qquad \pi \operatorname{Im} z(t) > 0.187 + 2.48(1 - t).$$
In particular, $\lambda$ and $z$ are uniformly bounded, and the damping exponent $\pi \operatorname{Im} z(t)$ is uniformly positive on $[0, 1]$, with $\operatorname{Im} z(t) \ge 0.187 / \pi$ in the worst case.

**Verification.** Direct calculation gives
$$\lvert \lambda(t) \rvert = \frac{1}{b \sqrt{t^{2} + h^{2}}}, \qquad \lvert z(t) \rvert^{2} = h^{2} + \frac{B^{2} - 2B h^{2}}{t^{2} + h^{2}}, \qquad \pi \operatorname{Im} z(t) = \pi h \left( \frac{B}{t^{2} + h^{2}} - 1 \right).$$
Since $B > 2 h^{2}$, all three decrease with $t$. Their largest / smallest values on $[0, 1]$ occur at endpoints, giving the displayed numerics. The sharper endpoint $5/6$ comes from plugging $t = 5/6$ into these formulas and using the rational bounds on $\pi$ and $\sqrt 3$.

**Tail bound (distant outputs of $S$).** For every positive integer $r$, every pair (sum of absolute node coefficients $\ge w$, absolute constant coordinate $\ge c$) of inputs satisfies
$$\sum_{n \in \mathcal N,\; n \ge r} \left( \lvert (S\mathbf{a})_{c, n} \rvert + \lvert (S\mathbf{a})_{d, n} \rvert \right) \le U_{r}(w, c) := \frac{800(w+c) e^{-0.70 r}}{1 - e^{-0.70}} + \frac{5.3(2.54 w / r + 2c) e^{-0.187 r}}{1 - e^{-0.187}}.$$

**Why it works.** For $n \ge r \ge 1$, the integrand in $K_{p_{\mathbf{a}}}$ carries an explicit exponential $e^{i \pi z(t) n}$; its real part is $e^{-\pi \operatorname{Im} z(t) n} \le e^{-0.187 n}$. The frequency derivative from $d/ds$ brings down a polynomial in $\lvert z(t) \rvert < 3$, controlled by the $2.54 w / r$ factor in the second term. On $[0, 5/6]$ the stronger $e^{-0.70 n}$ decay dominates, giving the first term; on $[5/6, 1]$ the weaker but polynomial-decay factor dominates, giving the second. The constants $800$ and $5.3$ are the total-variation product $36 \cdot \pi \cdot \lvert \text{dimension factor} \rvert$ and a uniform bound on the growth of $\lvert z \rvert$ on that segment.

**Operator norm.** Collecting the two bounds over all node columns (summing the $\le U_{1}(1, 1)$ bound, which is about $2500$) gives
$$\lVert S \rVert_{\text{op (finite lists)}} < 4000.$$
The specific bound $4000$ is what shows up everywhere later. It is not sharp; it is what the elementary estimates give.

**Uniform derivative bounds.** For every integer $j \ge 0$ and real $s \ge 0$:
$$\lvert p_{\mathbf{a}}^{(j)}(s) \rvert \le \pi^{j} V, \qquad \lvert K_{p_{\mathbf{a}}}^{(j)}(s) \rvert \le 3(3\pi)^{j} e^{-0.187 s} V,$$
where $V = 25 \lvert C \rvert + 36 \sum_{n} (\lvert c_{n} \rvert + \lvert d_{n} \rvert)$. The refined two-range bound separates the two frequency bands:
$$\lvert K_{p_{\mathbf{a}}}^{(j)}(s) \rvert \le 3(3\pi)^{j}(36w + 25c) e^{-0.70 s} + 1.25 (1.32 \pi)^{j} \left( \frac{2.54 w}{s} + 2 c \right) e^{-0.187 s}.$$

These bounds are uniform in the node locations and do not require finite support.

### 4.7 The coupled interpolation problem

The task of Stages A-C is:

**Find summable coefficient lists $\mathbf{a}_{1}, \mathbf{a}_{2}$ with prescribed extra coordinates $e_{i} = E_{i} q(k)$ (see Section 5) such that**
$$\mathcal{J} H_{1} = \mathcal{J} T_{k}, \qquad \mathcal{J} H_{2} = 0 \text{ on all positive nodes } n \in \mathcal N.$$

Once found, the pair $(H_{1}, H_{2})$ satisfies all four contacts of Theorem 4.2. The next two sections
- split the infinite problem into a finite block ($84$ nodes at $n \le 168$) and a tail;
- solve the finite block numerically, and prove the tail is contractive;
- bound the combined correction by $2 \cdot 10^{-8}$ in $\ell^{1}$ norm.

The Gaussian pair's existence follows by this fixed-point argument. The signs $H_{1} \le T_{k}$ and $H_{2} \ge 0$ (which are also part of Theorem 4.2) come later, in Section 7, from a Bernstein-polynomial certification.

---

## 5. Stage B: cutting to finite data

The interpolation problem of Section 4.7 lives in the infinite summable space. In practice it is solved in two steps:

1. **Truncate.** Restrict the positive nodes to those with $n \le 168$ and the spectral variable $s = b \lvert x \rvert^{2}$ to a finite quadrature grid. The truncated system is finite-dimensional; it is a $168 \times 171$ matrix equation.
2. **Correct for truncation.** Prove that the full infinite system differs from the finite one by a factor small enough that the finite solution, extended by zero, is within $2 \cdot 10^{-8}$ of the full solution.

Step 2 is the harder of the two analytically, and uses the tail bound $U_{r}(w, c)$ of Section 4.6.

### 5.1 Finite data

The finite truncation parameters are

$$M = 168, \qquad \mathcal F = \mathcal N \cap [1, M], \qquad N = 384.$$

There are exactly $84$ positive nodes in $\mathcal F$: each successive block of $12$ integers contains $6$ permitted residues $\{0, 1, 3, 4, 7, 9\}$, and $14 \times 6 = 84$.

A finite positive-node list has $168$ coordinates: first all $c_{n}$ in increasing $n$, then all $d_{n}$ in increasing $n$. Three extra coordinates are appended: $e = (c_{0}, d_{0}, C)$ in that order. So the complete finite input space has $171$ coordinates.

### 5.2 The Fejér quadrature rule

Recall the folded spectral measure of Lemma 4.4(c): atoms $C e_{j} P_{j}$ at $t_{j} = j/6$, and a density $d_{j}(t)$ on each open interval $((j-1)/6, j/6)$ for $1 \le j \le 6$. The atoms are kept exact (no quadrature error). The density on each interval is replaced by a finite quadrature sum.

The rule used is the **first Fejér rule**, split by interval. For each interval index $j \in \{1, \dots, 6\}$ and each quadrature index $v \in \{0, 1, \dots, N - 1\}$:

$$\alpha_{v} = \frac{v + 1/2}{N}, \qquad \tau_{j v} = \frac{2j - 1 + \cos(\pi \alpha_{v})}{12}, \qquad o_{v} = \frac{1}{N} \left( 1 - 2 \sum_{a = 1}^{N/2 - 1} \frac{\cos(2 a \pi \alpha_{v})}{4 a^{2} - 1} \right).$$

The points $\tau_{j v}$ lie in the *interior* of the interval: $0 < \pi \alpha_{v} < \pi$ implies $\cos(\pi \alpha_{v}) \in (-1, 1)$, hence $(2j - 2)/12 < \tau_{j v} < 2j / 12$. The weights $o_{v}$ are positive and sum to one on $[0, 1]$-normalized interval mass.

**What the quadrature replaces.** The density contribution to $p_{\mathbf{a}}(s)$ is
$$\int_{(j-1)/6}^{j/6} d_{j}(t)\, e^{i \pi t s}\, dt.$$
The quadrature replacement is
$$\sum_{v = 0}^{N-1} \frac{o_{v}}{6}\, d_{j}(\tau_{j v})\, e^{i \pi \tau_{j v} s}.$$
The factor of $1/6$ is the physical length of the interval.

**Quadrature error bound.** For every integrand of the form $d_{j}(t) h(t)$ with $h$ entire and $\lvert h^{(m)}(t) / m! \rvert \le R^{m}$ for $\lvert t - u \rvert \le 1/12$ around each midpoint $u = (2j - 1)/12$, the difference between the integral and the quadrature sum is bounded by $10^{-31}$.

**Why this scale.** Fejér's first rule is exact for polynomials of degree $\le N - 1 = 383$ on the standard interval. The Taylor series of $h$ on a disc of radius $1/12$ (smaller than the physical gap, so it stays inside the permitted region) converges fast enough that the tail after degree $383$ is negligible. Cauchy's coefficient estimate gives $\lvert h^{(m)}(u)/m! \rvert \le \text{bound}$, and the geometric tail is $\sum_{m \ge 384} 2^{-m} < 2^{-383}$. The result is dominated by $10^{83} \cdot 2^{-384}$, which is below $10^{-31}$. The exact numerical allowance in the paper is $4 \cdot 10^{83} \cdot 2^{-384} < 10^{-31}$; the integer inequality $4 \cdot 10^{114} < 2^{384}$ is what makes it work.

$10^{-31}$ is much smaller than every subsequent error budget ($10^{-8}$), so it will effectively be zero in the correction calculus.

### 5.3 The mass matrix $W$

The finite system is expressed as a matrix multiplication. Define the **mass matrix** $W$, with

$$\text{rows} = 6N + 7 = 2311, \qquad \text{columns} = 171.$$

Its rows are indexed by $(\text{j}, v) \in \{1, \dots, 6\} \times \{0, \dots, N-1\}$ (the density quadrature points, in lexicographic order) followed by $j \in \{0, \dots, 6\}$ (the seven atoms, in increasing $j$ order). Its columns use the finite-coordinate order specified in Section 5.1.

The entries are, for $t = \tau_{j v}$ and each finite node $n \in \mathcal F \cup \{0\}$:

$$W_{(j, v), c_{n}} = -\frac{2 \pi^{2} o_{v}}{6}\, e^{-i \pi t n} \sum_{\ell = j}^{6} (t_{\ell} - t)\, P_{\ell}\, e^{i \pi t_{\ell} n},$$
$$W_{(j, v), d_{n}} = \frac{2 \pi i o_{v}}{6}\, e^{-i \pi t n} \sum_{\ell = j}^{6} P_{\ell}\, e^{i \pi t_{\ell} n},$$
$$W_{(j, v), C} = 0,$$

and for atoms at $j \in \{0, \dots, 6\}$:

$$W_{\text{atom}(j), C} = e_{j} P_{j}, \qquad W_{\text{atom}(j), c_{n}} = W_{\text{atom}(j), d_{n}} = 0,$$

where $e_{0} = 1$ and $e_{j} = 2$ for $j > 0$.

**What $W\mathbf{a}$ computes.** For any input column $\mathbf{a}$ (a $171$-vector), the vector $W\mathbf{a} \in \mathbb C^{2311}$ has entries equal to the finite quadrature values of $p_{\mathbf{a}}$ (or more precisely, of the folded density and atomic contributions) at the mass points. Multiplying by the moment factors (Section 5.4) recovers the derivatives of $p_{\mathbf{a}}$ at any real center.

Two notes.

1. The column ordering matters: the $c_{n}$ and $d_{n}$ blocks are grouped so that $W$ can be split cleanly into the "positive coordinates" block (columns $1, \dots, 168$) and the "extras" block (columns $169, \dots, 171$).
2. All entries are exact rationals times exponentials $e^{i \pi t n}$; later we will floor-round to rationals with a floor $\epsilon = 10^{-80}$ (Appendix B).

### 5.4 Moment factors and Taylor jets on the finite side

The mass matrix $W$ produces the *zeroth* finite value of $p_{\mathbf{a}}$ at each mass point. To get higher derivatives, multiply by the **moment factors**. Define, for integers $a \ge 0$, real numbers $m, y$, and each mass point $t$:

$$F_{a}(m, y; t) = e^{i \pi t m} \frac{(i \pi y t)^{a}}{a!}, \qquad G_{a}(m, y; t) = \lambda(t)\, e^{i \pi z(t) m} \frac{(i \pi y z(t))^{a}}{a!}.$$

For a real coefficient column $A \in \mathbb R^{171}$, the Taylor coefficient of order $a$ at $m$ (scaled by $y^{a}/a!$) of $p_{\mathbf{a}}$ is
$$\frac{y^{a}}{a!} p_{\mathbf{a}}^{(a)}(m) = \operatorname{Re}\bigl( F_{a}(m, y; t)^{\top} W A \bigr).$$
Similarly, for the transformed cardinal function,
$$\frac{y^{a}}{a!} K_{p_{\mathbf{a}}}^{(a)}(m) = \operatorname{Re}\bigl( G_{a}(m, y; t)^{\top} W A \bigr).$$

**No conjugation.** Note that $F_{a}$ and $G_{a}$ are *vectors* (over the $2311$ mass points), and the multiplication is an ordinary dot product (adjoint $)^{\top}$, without taking the complex conjugate. This is important: the construction uses the linearity of the finite sum and differentiates it termwise; complex conjugation would convert the derivative to a conjugate derivative and give the wrong signs.

**Why this works on the finite side only.** The formulas use $W A$, which is the *finite* quadrature evaluation of $p_{\mathbf{a}}$. This is not the same as the exact $p_{\mathbf{a}}(m)$; the two differ by the quadrature error, which is $10^{-31}$. So the finite moment products approximate the exact Taylor coefficients to that precision.

### 5.5 The finite operator, the target, and the reference columns

We now assemble the finite interpolation problem.

**The finite operator $S^{\#}$.** With $W$ the mass matrix as above, the $168 \times 171$ matrix representing "the transformed jet operator on finite data" is:

$$S^{\#}_{c(m), :} = -\frac{1}{Q_{m}} \operatorname{Re}\bigl( G_{0}(m, 1; t)^{\top} W \bigr),$$
$$S^{\#}_{d(m), :} = \frac{1}{Q_{m}} \operatorname{Re}\bigl( D_{m} G_{0}(m, 1; t)^{\top} W - G_{1}(m, 1; t)^{\top} W \bigr),$$

for each finite node $m \in \mathcal F$. The rows are indexed by the same $168$ coordinates as before (all $c_{m}$ in increasing $m$, then all $d_{m}$).

**Interpretation.** For a finite input $\mathbf{a} \in \mathbb R^{171}$, the vector $S^{\#} \mathbf{a} \in \mathbb R^{168}$ is the finite quadrature evaluation of $\mathcal{J} K_{p_{\mathbf{a}}}$ at the finite nodes. The sign convention: $S^{\#}$ is the negative of $\mathcal{J} K_{p}$ (matching the definition of $S$ in Section 4.6).

**The target column $Y$.** The right-hand side at five "fundamental nodes" $n_{1}, \dots, n_{5} = 1, 3, 4, 7, 9$ is the target $T_{k}$. Because $T_{k}$ is an exponential, its jet projection has a compact form. Let $q(k) \in \mathbb R^{10}$ be the vector
$$q(k) = \left( T_{k}(1), T_{k}(3), T_{k}(4), T_{k}(7), T_{k}(9), \; -T_{k}'(1), \dots, -T_{k}'(9) \right).$$

Define $Y$ by
$$Y_{d(n_{j}), j} = -\frac{1}{Q_{n_{j}}}, \quad Y_{c(n_{j}), j + 5} = \frac{1}{Q_{n_{j}}}, \quad Y_{d(n_{j}), j + 5} = -\frac{D_{n_{j}}}{Q_{n_{j}}}, \quad j = 1, \dots, 5,$$
with all other entries zero.

**What $Y q(k)$ is.** At the nodes $n \in \{1, 3, 4, 7, 9\}$, $Y q(k)$ gives the target jet $\mathcal{J} T_{k}(n)$. At every other finite node, $Y q(k) = 0$. This is the finite target data, exact, no quadrature error, because $T_{k}$ evaluated at integers is a finite expression.

The target at the remaining finite nodes (all those that are not in $\{1, 3, 4, 7, 9\}$) is small -- bounded by a geometric tail, $< 10^{-11}$ in $\ell^{1}$ norm -- so it will be treated as an error term in the correction argument. This is the "5 nodes carry the data, the rest are handled as perturbations" structure.

**The extra coordinates.** The three extra coordinates $e_{i} = (c_{0}, d_{0}, C)$ for list $i$ are chosen explicitly, independent of $k$. Collect them into vectors:

- $E_{1} \in \mathbb R^{3 \times 10}$: row of $e$ coordinates for $H_{1}$'s list.
- $E_{2} \in \mathbb R^{3 \times 10}$: row of $e$ coordinates for $H_{2}$'s list.

Concretely, the first column of $E$ (viewed as $6 \times 10$) is
$$(0.567, \; 0.250, \; -0.0074, \; 0, \; -0.209, \; 0.0096)^{\top}$$
and the sixth column is
$$(-0.976, \; 0.408, \; 0.0115, \; 1, \; 0.938, \; -0.017)^{\top}.$$
The first three entries of each column give $E_{1}$, the last three give $E_{2}$. So $E_{i} q(k) = $ (extras of list $i$ for target data $q(k)$).

**Why this choice.** The extras matter because the coupled equations for $H_{1}, H_{2}$ include cross-contributions $S$ from the *other* list's extras. Without an $E_{i}$ source, the constant-coordinate column would be zero and the finite system would be degenerate. The specific numeric values (.567, .250, ...) were found by solving for the extras that put the full infinite interpolation system into a well-posed regime; they are documented and verifiable, but the "why these specific numbers" question is really "why any numbers with these properties work" -- the answer is an interval of valid choices, of which these are a convenient point.

**The coupled finite system.** Splitting $S^{\#}$ into its finite-positive columns $S^{\#}_{F F}$ ($168 \times 168$), finite extras columns $S^{\#}_{F e}$ ($168 \times 3$), and writing $D = -S^{\#}_{F F}$ (so that $D$ is the finite off-diagonal in the block equation), the finite version of the coupled interpolation system reads

$$A_{1} - D A_{2} = Y + S^{\#}_{F e} E_{2}, \qquad A_{2} - D A_{1} = S^{\#}_{F e} E_{1}.$$

Here $A_{i} \in \mathbb R^{168 \times 10}$ are the finite-positive parts of the two lists, one column per coordinate of $q(k)$. The sign of $D$ follows from the sign convention in $S$: $-S$ on the diagonal becomes $+D$ when the equation is written in this "shifted" form.

**The finite geometric-sum solution.** Summing and differencing the two equations,
$$(I - D)(A_{1} + A_{2}) = g_{1} + g_{2}, \qquad (I + D)(A_{1} - A_{2}) = g_{1} - g_{2},$$
where
$$g_{1} = Y + S^{\#}_{F e} E_{2}, \qquad g_{2} = S^{\#}_{F e} E_{1}.$$

For the finite system, we replace the exact inverses $(I - D)^{-1}$ and $(I + D)^{-1}$ by truncated geometric sums:
$$V_{+} = \sum_{j = 0}^{511} D^{j}, \qquad V_{-} = \sum_{j = 0}^{511} (-D)^{j}.$$

Then
$$X_{1 \mathcal F} = \frac{1}{2}\left( V_{+}(g_{1} + g_{2}) + V_{-}(g_{1} - g_{2}) \right), \qquad X_{2 \mathcal F} = \frac{1}{2}\left( V_{+}(g_{1} + g_{2}) - V_{-}(g_{1} - g_{2}) \right),$$
and
$$X_{i} = \begin{pmatrix} X_{i \mathcal F} \\ E_{i} \end{pmatrix} \in \mathbb R^{171 \times 10}, \qquad i = 1, 2.$$

The reference list for a given $k$ is then $X_{i} q(k)$, extended by zero on the positive nodes $n > 168$ (the "tail").

**Error from the truncation of the geometric sum.** This is where the tables in Section 5.6 come in. The identities
$$(I - D) V_{+} = I - D^{512}, \qquad (I + D) V_{-} = I - D^{512}$$
show that the finite geometric sums solve the equations up to $D^{512}$. If $\lVert D^{512} \rVert$ is small enough, the difference between $X_{i} q(k)$ and the true solution is small. Section 5.6 verifies $D^{2^9} = D^{512}$ has norm below $10^{-12}$ (row 8 of the power table).

**What is left.** Given the reference lists $X_{i} q(k)$, the full infinite solution differs from them by a correction that lives in the difference between the finite quadrature and the exact spectral integrals (bounded by $10^{-29}$), the difference between the finite nodes and all nodes (bounded by $10^{-14}$ at $168$ from the tail bound $U_{168}$), and by the difference between $D^{512}$ and zero. Section 6 assembles these bounds into a $2 \cdot 10^{-8}$ correction.

### 5.6 The power-sum computation and the two tables

The argument of Section 5.5 requires explicit bounds on the powers of $D$ that enter the geometric sums. These are computed by **repeated squaring**, and the bounds are strict rationals at each step.

**Algorithm.** Start with $U = D$ and $V_{+} = V_{-} = I$. Then for $r = 0, \dots, 8$:

- Update $V_{+} \leftarrow V_{+}(I + U)$.
- Update $V_{-} \leftarrow V_{-}(I + \varepsilon_{r} U)$ where $\varepsilon_{0} = -1$ and $\varepsilon_{r} = 1$ for $r > 0$.
- Record the column-sum norms of $U$, $V_{+}$, $V_{-}$.
- Square $U \leftarrow U^{2}$ (before the next step).

After step 8 has been executed, $U = D^{512}$ is discarded, and $V_{+} = \sum_{j=0}^{511} D^{j}$, $V_{-} = \sum_{j=0}^{511} (-D)^{j}$, exactly the sums used in Section 5.5.

**Why this works.** Before step $r$ (i.e. at the top of the loop with the current $U$), $U = D^{2^{r}}$. Then $V_{+}$ is being multiplied by $(I + D^{2^{r}})$, and
$$\prod_{r=0}^{9} (I + D^{2^{r}}) = \sum_{j=0}^{511} D^{j},$$
the standard binary identity. The $V_{-}$ product is analogous: it uses $(I - D)(I + D)(I + D^{2}) \cdots$, which is the alternating sum $\sum_{j=0}^{511} (-D)^{j}$.

**Power table.** The bounds are verified numerically; the entries below are strict upper bounds on the column-sum norm at each stage.

| Step $r$ | $\lVert U \rVert$ | $\lVert V_{+} \rVert$ | $\lVert V_{-} \rVert$ |
|---|---|---|---|
| 0 | 1.74 | 2.74 | 2.1 |
| 1 | 2.09 | 6.94 | 2.2 |
| 2 | 2.27 | 15.32 | 2.13 |
| 3 | 1.75 | 25.3 | 1.97 |
| 4 | 0.75 | 31.1 | 1.88 |
| 5 | 0.12 | 32.09 | 1.85 |
| 6 | 0.0031 | 32.12 | 1.85 |
| 7 | 0.000002 | 32.12 | 1.85 |
| 8 | 0.000000000001 | 32.12 | 1.85 |

**Reading the table.**

- Row 6 gives $\lVert D^{64} \rVert \le 0.0031$, a serious power contraction.
- Row 8 gives $\lVert D^{256} \rVert \le 10^{-12}$.
- Therefore $\lVert D^{512} \rVert \le (10^{-12})^{2} = 10^{-24}$: the geometric-sum truncation error is essentially zero compared to every other error budget.

This is the entire mechanism by which the finite geometric sums $V_{\pm}$ compute the reference lists. The bounds are discovered by running the algorithm and checking each column-sum in rational arithmetic; they are not derived from general principles. The construction is empirical at this step, and the resulting numerical table is what makes the proof go through.

**Column-norm table.** Each reference column has a bounded $\ell^{1}$ norm; the bounds are:

| Quantity | Columns 1-5 | Columns 6-10 |
|---|---|---|
| $\max_{j} \lVert (X_{1})_{j} \rVert_{1}$ | 1.29 | 2.65 |
| $\max_{j} \lVert (X_{2})_{j} \rVert_{1}$ | 1.01 | 3.03 |
| $\max_{j} \sum_{i=1}^{2} \lVert (X_{i})_{j} \rVert_{1}$ | 2.20 | 5.40 |

The two right-hand columns correspond to the $Y$-target columns (indices $1, \dots, 5$) and the derivative-target columns (indices $6, \dots, 10$). The combined column norms are what appear in the reference-list norm bound:
$$\sum_{i=1}^{2} \lVert X_{i} q(k) \rVert_{1} < 4.8 \quad \text{(for every } k \ge 2.36\text{)}.$$

### 5.7 Verification

Every entry of both tables is a strict rational inequality that can be checked by exact arithmetic on the finite matrix $D$. The two accompanying Python programs do this:

- `verification/arithmetic_bounds.py`: checks the scalar budgets and the numerical error-propagation majorants, using 512-bit Arb balls and exact rational arithmetic.
- `verification/numeric_balls.py`: encloses the exact quadrature arrays at 256-bit precision and tests every strict comparison in the power and column tables, plus the sign tables in Section 7.

Both are documented in `verification/README.md`. Together they certify the finite data. The passage back to the exact infinite spectral integrals is supplied by the quadrature-error lemma of Section 5.2.

---

## 6. Stage C: exact interpolation

We now prove that the finite reference lists of Section 5 can be corrected to *exact* solutions of the full infinite interpolation system, with a correction bounded by $2 \cdot 10^{-8}$ in $\ell^{1}$.

### 6.1 The proposition

**Proposition (Exact interpolation).** *For every real $k \ge 2.36$, there exists a unique pair of summable positive-node coefficient lists $\mathbf{a}_{1}, \mathbf{a}_{2} \in \mathcal E$ (where $\mathcal E$ is the space of lists with prescribed extra coordinates $e_{i} = E_{i} q(k)$), whose full lists $\ell_{1}, \ell_{2}$ satisfy the exact jets*
$$H_{1}(n) = T_{k}(n), \quad H_{1}'(n) = T_{k}'(n), \quad H_{2}(n) = 0, \quad H_{2}'(n) = 0 \qquad (n \in \mathcal N),$$
*where $H_{1} = p_{\ell_{1}} + K_{p_{\ell_{2}}}$, $H_{2} = p_{\ell_{2}} + K_{p_{\ell_{1}}}$.*

*Let $\ell_{i}^{\mathrm{ref}} = X_{i} q(k)$, extended by zero on the tail $\mathcal T$ (all positive nodes $n > 168$). Then*
$$\sum_{i=1}^{2} \lVert \ell_{i}^{\mathrm{ref}} \rVert_{1} < 4.8, \qquad \sum_{i=1}^{2} \lVert \ell_{i} - \ell_{i}^{\mathrm{ref}} \rVert_{1} < 2 \cdot 10^{-8}, \qquad \sum_{i=1}^{2} \lVert \ell_{i} \rVert_{1} < 5.$$

*The functions $H_{1}, H_{2}$ are entire and real on the real axis. In particular $H_{1} - T_{k}$ and $H_{2}$ have double zeros at every positive node, including the artificial nodes that are not actual lattice radii.*

The third inequality (total $\ell^{1}$ norm $< 5$) is the sum of the reference bound $4.8$ and the correction bound $2 \cdot 10^{-8}$.

**Uniqueness.** The bounded inverse below gives uniqueness, within the prescribed coefficient system, for summable lists with fixed extras and jets prescribed on the enlarged periodic node set $\mathcal N$.

This is worth contrasting with what is *not* claimed. Talebizadeh Sardari (2021, Theorem 1.7) proved that infinitely many linearly independent radial Schwartz functions have all their values and first derivatives (and those of their Fourier transforms) *vanish on the actual triangular shells* alone. So the interpolation problem that prescribes jets only on $A \setminus \{0\}$ is *not* well-posed. The present paper enlarges the node set to $\mathcal N$, which contains all the actual shells and many artificial extras. The well-posedness is on the larger set, within the fixed summable coefficient system.

### 6.2 The nonlinear block equation

For summable lists with fixed extras, the coupled interpolation system of Section 4.7 becomes a nonlinear equation in the unknown tail $\mathbf{a}_{\mathcal T} = (\mathbf{a}_{1 \mathcal T}, \mathbf{a}_{2 \mathcal T})$:

$$\mathbf{a}_{\mathcal F} = \mathcal L_{\mathcal F \mathcal F}^{-1} (r_{\mathcal F} - \mathcal L_{\mathcal F \mathcal T} \mathbf{a}_{\mathcal T}),$$

$$\mathbf{a}_{\mathcal T} = \mathcal L_{\mathcal T \mathcal T}^{-1} (r_{\mathcal T} - \mathcal L_{\mathcal T \mathcal F} \mathbf{a}_{\mathcal F}),$$

where

- $\mathbf{a}_{\mathcal F}$ is the finite-positive part (the first $168$ coordinates of each list);
- $\mathbf{a}_{\mathcal T}$ is the tail (nodes $n > 168$);
- $\mathcal L$ is the $2 \times 2$ block system satisfying the coupled equation from Section 4.7, split into finite/finite, finite/tail, tail/finite, tail/tail blocks;
- $r_{\mathcal F}, r_{\mathcal T}$ are the residual vectors (target minus the current approximation), bounded below.

The Schur complement of this system, with $\mathbf{a}_{\mathcal F}$ eliminated, is a tail-space equation
$$\mathcal Q \mathbf{a}_{\mathcal T} = \psi_{\mathcal T}, \qquad \mathcal Q = \mathcal L_{\mathcal T \mathcal T} - \mathcal L_{\mathcal T \mathcal F} \, \mathcal L_{\mathcal F \mathcal F}^{-1} \mathcal L_{\mathcal F \mathcal T},$$
where $\psi_{\mathcal T}$ is the effective right side.

### 6.3 The operator bounds

Four bounds establish that $\mathcal Q$ is a small perturbation of the identity.

**Bound on the finite inverse.** $\lVert \mathcal L_{\mathcal F \mathcal F}^{-1} \rVert < 40$. This follows from the finite column-sum bound: the finite operator, viewed as a $336 \times 336$ block matrix (two lists combined), has column sums bounded by the reference-list column norms.

**Bound on the tail/tail block.** $\lVert \mathcal L_{\mathcal T \mathcal T} - I \rVert < \delta$, where $\delta = 1.1 \cdot 10^{-14}$. The tail decay is controlled by $U_{168}(w, c)$ from Section 4.6, with the specific values $w = 5$ and $c = 0.02$ from the reference-list bounds.

**Bound on the tail/finite block.** $\lVert \mathcal L_{\mathcal T \mathcal F} \rVert < \delta$. Same mechanism; the transformed contribution from a finite input vanishes rapidly at distance $n > 168$.

**Bound on the finite/tail block.** $\lVert \mathcal L_{\mathcal F \mathcal T} \rVert < 4000$. This uses the operator norm of $S$ (Section 4.6). It is not small; this is why one cannot just invert the finite block naively. But the Schur complement absorbs it.

### 6.4 The Schur complement is invertible

Combining the tail/finite, finite/finite, finite/tail, and tail/tail bounds gives
$$\lVert \mathcal Q - I \rVert < \delta(1 + 40 \cdot 4000) = 1.760011 \cdot 10^{-9} < 1.$$

Since this is $< 1$, $\mathcal Q$ has a convergent Neumann series $\mathcal Q^{-1} = \sum_{k \ge 0} (I - \mathcal Q)^{k}$, and $\mathcal Q$ is invertible with
$$\lVert \mathcal Q^{-1} \rVert < \frac{1}{1 - 1.76 \cdot 10^{-9}} \approx 1 + 1.76 \cdot 10^{-9}.$$

Explicitly, for any right side $\psi = (\psi_{\mathcal F}, \psi_{\mathcal T})$, the solution is reconstructed as
$$\mathbf{a}_{\mathcal T} = \mathcal Q^{-1} \left( \psi_{\mathcal T} - \mathcal L_{\mathcal T \mathcal F} \mathcal L_{\mathcal F \mathcal F}^{-1} \psi_{\mathcal F} \right),$$
$$\mathbf{a}_{\mathcal F} = \mathcal L_{\mathcal F \mathcal F}^{-1} \left( \psi_{\mathcal F} - \mathcal L_{\mathcal F \mathcal T} \mathbf{a}_{\mathcal T} \right).$$

This is the exact full solution for any given right side $\psi$. The question is now: how far is our reference list $X q(k)$ from this exact solution?

### 6.5 The three sources of residual

The reference list $\ell^{\mathrm{ref}} = (X_{1} q(k), X_{2} q(k))$ satisfies the interpolation system approximately. Three sources of discrepancy:

1. **The target at the ten "carrier" nodes.** The five $q(k)$-target nodes $\{1, 3, 4, 7, 9\}$ are hit exactly. But the target $T_{k}$ is nonzero at the other $79$ finite nodes; the reference list gives $Y q(k) = 0$ there. The total target loss from omitted finite nodes is bounded by
$$\frac{1 + (1 + 2.12)/2.36}{1.75} \sum_{n = 12}^{\infty} e^{-2.36(n - 1)} < 8 \cdot 10^{-12}.$$

2. **The quadrature vs. exact measure.** Section 5.2's quadrature error is $10^{-31}$ per quadrature point; multiplied by $2311$ entries it stays well below $10^{-24}$.

3. **The geometric-sum truncation.** As computed in Section 5.6, $\lVert D^{512} \rVert < 10^{-24}$. Three-times the residual from omitted $D$ powers, times the zero-tail norm, is negligible.

Adding these three bounds:
$$\lVert r_{\mathcal F} \rVert_{1} < 8 \cdot 10^{-12} + 10^{-24} + 1.4 \cdot 10^{-15} < 10^{-11}.$$
(The middle number is a target-tail calculation; the last is the geometric-sum bound.)

**Tail residual.** At tail nodes $n > 168$, the reference list has zero positive-node coefficients, but the transformed operator still acts. Its total norm there is dominated by
$$\lVert r_{\mathcal T} \rVert_{1} \le U_{168}(5, 0.02) + 10^{-170} < 8.3 \cdot 10^{-14}.$$
The $.02$ is the combined extras constant bound (Section 6.6 below); the $10^{-170}$ is the target $T_{k}$ decayed to node $n = 169$ and beyond, via the local geometric majorant of Section 8.

### 6.6 The extras constant and reference corner

The three extras coordinates satisfy, for $E_{i} q(k)$ with $k \ge 2.36$,
$$\lvert E_{i} q(k) \rvert < .02 \quad \text{(each coordinate)}.$$
This is verified numerically: the specific values in Section 5.5 give $\lvert -.0074 + .0115/k \rvert + \lvert .0096 - .017/k \rvert = .017 - .0285/k < .02$ for $k \ge 2.36$.

Remember that the sum of the three extra coordinates, plus the reference list's other coordinates, is $< 4.8$ in $\ell^{1}$ norm (Section 5.5). Adding the tail $10^{-170}$ and extras $.02$ gives $\lVert \ell^{\mathrm{ref}} \rVert < 4.8 + 10^{-170}$ still $< 5$.

### 6.7 The correction size

Put $u_{\mathcal F} = \lVert (\mathbf{a} - \mathbf{a}^{\mathrm{ref}})_{\mathcal F} \rVert_{1}$ and $u_{\mathcal T} = \lVert (\mathbf{a} - \mathbf{a}^{\mathrm{ref}})_{\mathcal T} \rVert_{1}$. The block equations for the correction $(\mathbf{a} - \mathbf{a}^{\mathrm{ref}})$, together with the Schur bound, give the system

$$u_{\mathcal F} \le 40 \left( 10^{-11} + 4000 \, u_{\mathcal T} \right),$$
$$u_{\mathcal T} \le 8.3 \cdot 10^{-14} + 1.1 \cdot 10^{-14} (u_{\mathcal F} + u_{\mathcal T}).$$

Substituting the first inequality into the second:
$$u_{\mathcal T} \le \frac{8.3 \cdot 10^{-14} + 4.4 \cdot 10^{-24}}{1 - 1.760011 \cdot 10^{-9}} < 8.31 \cdot 10^{-14},$$
and from that
$$u_{\mathcal F} < 1.37 \cdot 10^{-8}.$$

Their sum is $< 2 \cdot 10^{-8}$, which is the claimed correction bound.

Adding $u_{\mathcal F} + u_{\mathcal T}$ to the reference norm $\lVert \ell^{\mathrm{ref}} \rVert < 4.8$ gives $\lVert \ell \rVert < 5$ for the exact list.

The extras did not change in the correction (they are pinned by the prescribed $E_{i} q(k)$), so this is the full-list bound, not just the positive-node part.

### 6.8 Finishing the proof of Proposition 6.1

The exact solution $\mathbf{a}$ of the full interpolation system of Section 6.2 exists and is unique (the Schur complement is invertible). Its lists $\ell_{1}, \ell_{2}$ satisfy
$$\mathcal J H_{1} = \mathcal J T_{k}, \qquad \mathcal J H_{2} = 0$$
at every positive node, and hence by definition of $\mathcal J$, the four contacts
$$H_{1}(n) = T_{k}(n), \qquad H_{1}'(n) = T_{k}'(n), \qquad H_{2}(n) = 0, \qquad H_{2}'(n) = 0$$
at every $n \in \mathcal N$.

The $H_{i}$ are entire and real on the real axis by construction (Section 4.5). Lemma 4.5 gives them as the two halves of a real radial Schwartz Fourier pair, which is $f$ and $\widehat f$ of the normalized Gaussian theorem.

This completes the existence half of Theorem 4.2. The signs $H_{1} \le T_{k}$ and $H_{2} \ge 0$ are established in Section 7.

---

## 7. The signs at every radius and every parameter

We now prove the two sign conditions
$$H_{1}(s) \le T_{k}(s), \qquad H_{2}(s) \ge 0 \qquad (s \ge 0, \; k \ge 2.36).$$
The contacts of Proposition 6.1 give equality at every node. So the sign conditions only need to be proved *strictly between* nodes.

### 7.1 The problem, stated

Recall:
$$H_{i}(s) = p_{\ell_{i}}(s) + K_{p_{\ell_{3-i}}}(s), \qquad T_{k}(s) = k^{-1} e^{-k(s-1)}, \qquad i = 1, 2.$$
We want $H_{1} - T_{k} \le 0$ and $H_{2} \ge 0$. The signs $-\sigma_{i}$ with $\sigma_{1} = -1$, $\sigma_{2} = 1$ unify the two statements:
$$\sigma_{1}(H_{1} - T_{k}) \ge 0, \qquad \sigma_{2} H_{2} \ge 0, \qquad \text{i.e.} \quad \sigma_{i}(H_{i} - \delta_{i} T_{k}) \ge 0$$
with $\delta_{1} = 1$, $\delta_{2} = 0$.

Since we have *exact* contacts at nodes, both sides vanish at every $n \in \mathcal N$. So the difference is *twice differentiable through zero* there. That is what makes the proof work.

### 7.2 The deleted-slack reduction

Define, for a twice-continuously-differentiable function $u$ and a node $m$,
$$\mathcal D_{m, 2} u(s) = \frac{u(s) - u(m) - (s - m) u'(m)}{(s - m)^{2}} \qquad (s \ne m),$$
with removable value $\mathcal D_{m, 2} u(m) = u''(m)/2$. At the origin, we set $\mathcal D_{0, 0} u = u$. Define also the "slack" function
$$F_{i, \mathrm{slack}}(s) = \sigma_{i}(H_{i}(s) - \delta_{i} T_{k}(s)),$$
so that the sign we want is $F_{i, \mathrm{slack}}(s) \ge 0$.

**Key identity.** At every positive node $m \in \mathcal N$, both $H_{i}$ and $\delta_{i} T_{k}$ have the target jets (by the contacts), so:
$$F_{i, \mathrm{slack}}(s) = (s - m)^{2}\, \mathcal D_{m, 2}\bigl(\sigma_{i}(H_{i} - \delta_{i} T_{k})\bigr)(s).$$
The sign $F_{i, \mathrm{slack}}(s) \ge 0$ is therefore equivalent to
$$\sigma_{i}\, \mathcal D_{m, 2}\bigl( H_{i} - \delta_{i} T_{k} \bigr)(s) \ge 0.$$
(The sign of $(s - m)^{2}$ is positive, so it drops out.)

The strategy is now clear:
- For each half-gap $[m, m']$ covered by $s \in [0, s_{0}]$ with $s_{0} = 89.5$, prove the quotient inequality by finite polynomial computation.
- For $s \ge s_{0}$, prove the same inequality using a quadratic barrier for the sine product $P(s)$.

### 7.3 Bernstein rows and the half-gap reduction

We construct polynomial rows to which we will apply a Bernstein-polynomial argument. Recall that on $[0, 1]$, the *Bernstein basis functions* $\binom{40}{d} x^{d}(1 - x)^{40 - d}$ for $d = 0, \dots, 40$ are nonnegative and sum to one. So if a polynomial has all Bernstein coefficient $> c$, then it is $> c$ on the whole interval $[0, 1]$.

For a fixed center $m \in \mathcal N_{0} \cap [0, 88]$ and a signed half-gap length $y$, the half-gap consists of $s = m + yx$, $0 \le x \le 1$. The Taylor coefficients of the "deleted quotient" are
$$\mathcal D_{m, \nu} u(m + yx) = \sum_{r = 0}^{\infty} \frac{y^{r}}{(r + \nu)!} u^{(r + \nu)}(m) x^{r},$$
where $\nu = 0$ if $m = 0$ and $\nu = 2$ if $m > 0$.

Define the **Bernstein Taylor row** at center $m$, signed half-gap $y$:
$$l_{i, r} = \sigma_{i} \frac{y^{r}}{(r + \nu)!} \left( \widetilde p_{X_{i}}^{(r + \nu)}(m) + \widetilde K_{X_{3 - i}}^{(r + \nu)}(m) \right), \qquad 0 \le r \le 40,$$
where the tilde denotes the finite quadrature evaluation (an approximation to the exact $p_{X_{i}}^{(r + \nu)}(m)$ and $K_{X_{3 - i}}^{(r + \nu)}(m)$).

And define the **Bernstein row**
$$L_{i, d} = \sum_{r = 0}^{d} \omega_{d, r}\, l_{i, r}, \qquad \omega_{d, r} = \frac{\binom{d}{r}}{\binom{40}{r}}.$$

By design, $L_{i, d}$ is the $d$-th Bernstein coefficient of the partial Taylor series $\sum_{r=0}^{40} l_{i, r} x^{r}$ of the deleted quotient, with the **partial Taylor series being of the reference lists** $X_{i} q(k)$ rather than the exact lists $\ell_{i}$.

**The correction from reference to exact.** By (5.4), the difference between reference and exact is controlled by $\lVert \ell_{i} - X_{i} q(k) \rVert_{1} < 2 \cdot 10^{-8}$. The derivative of a cardinal function is bounded by $\pi^{j} V$, and the derivative of the transformed by $3(3\pi)^{j} e^{-.187 s} V$, with $V \le \min(36 \cdot 5, 180)$. The correct difference at every step is bounded by
$$\sup_{s \text{ in half-gap}} \left| \mathcal D_{m, \nu} (H_{i} - H_{i}^{\mathrm{ref}})(s) \right| < .00011.$$
This uses $\nu = 0$ (zeroth derivative) or $\nu = 2$ (with the additional $1/2$ factor from $\int_0^1 (1-v) dv$), and the maximum of $\pi^{\nu}, 3(3\pi)^{\nu}$.

Similarly, the finite-truncation error in the Bernstein rows (from truncating the Taylor series after 40 terms) is bounded by $10^{-6}$. And the quadrature error is $41 \cdot 4 \cdot 5 \cdot 10^{-31} = 8.2 \cdot 10^{-29}$.

### 7.4 Parameter boxes

The parameter $k$ ranges over $[2.36, \infty)$. Divide into five intervals:
$$[2.36, 2.65], \quad [2.65, 3.2], \quad [3.2, 4], \quad [4, 6], \quad [6, \infty).$$
On each, we replace the curve $q(k)$ by a rectangle (containing box) and use the minimum over box vertices. This is valid because the expression below is *affine in each coordinate* of the box.

**Lemma (Containing-box bound).** *For consecutive endpoints $K, J$ of the partition, real row $l = (a_{1}, \dots, a_{5}, b_{1}, \dots, b_{5})$, and $K \le k \le J$ (no upper bound if $J = \infty$), one has*
$$l q(k) \ge [l]_{K, J} := \min_{x \in \{1/K, 1/J\}} \left( a_{1} + x b_{1} + \sum_{j = 2}^{5} \min_{H \in \{K, J\}} \left\{ (a_{j} + x b_{j}) e^{-H(n_{j} - 1)} \right\} \right).$$
*At $J = \infty$ use $1/J = 0$ and $e^{-J(n_{j} - 1)} = 0$ for $j > 1$.*

**Sketch.** The physical vector $q(k)$ lies in the box $x = 1/k \in [1/J, 1/K]$, $w_{j} = e^{-k(n_{j} - 1)} \in [e^{-J(n_{j} - 1)}, e^{-K(n_{j} - 1)}]$. On this larger box $l q$ is affine in each coordinate while others are fixed, so its minimum occurs at a vertex. First minimizing each $w_{j}$ and then the single shared $x$ gives the displayed formula. This is a *containing-box bound*, not sampling of the curve $q(k)$. $\square$

### 7.5 Target support on every half-gap

The right side of the Bernstein polynomial must dominate the deleted quotient of $T_{k}$, not just of $H_{i}$. The target's own deleted quotient has its own Taylor coefficients $a_{r}^{*}$, and these are bounded uniformly over $k$.

**Lemma (Exponential supports).** *For consecutive parameter endpoints $K, J$ and $K \le k \le J$, the target's deleted quotient satisfies, on every prescribed half-gap,*
$$\mathcal D_{m, \nu} T_{k}(m + y x) \ge \sum_{r = 0}^{40} a_{r}^{*} x^{r} - 10^{-6} \qquad (0 \le x \le 1),$$
*with coefficients:*

- *At the origin ($m = 0$, $\nu = 0$): $a_{r}^{*} = (e^{K} / K) \cdot (-Ky)^{r} / r!$.*
- *At node 1 (and other isolated nodes with $m = 1$): $a_{r}^{*} = K \cdot (-Ky)^{r} / (r + 2)!$.*
- *At node 3 (handling signed half-gaps): $a_{0}^{*} = (J/2) e^{-13 J/6}$ when $J < \infty$; more complex expression at other $y$.*

**Sketch.** At the origin, $T_{k}(s) = e^{-k(s-1)} / k$, and $\partial T_{k} / \partial k = T_{k}(s)(1 - s - 1/k) > 0$ for $k \ge 2.36 > 2$ and $s \le 1/2$. So replacing $k$ by $K$ lowers the target. The Taylor polynomial through degree 40 is exactly $(e^{K}/K)(-Ky)^{r}/r!$, giving case 1.

At $m = 1$, the deleted quotient is $\mathcal D_{1, 2} T_{k}(1 + u) = (e^{-ku} - 1 + ku)/(k u^{2})$, with value $k/2$ at $u = 0$. The $k$-derivative is $(1 - (1 + ku) e^{-ku})/(k^{2} u^{2}) \ge 0$ by $e^{x} \ge 1 + x$, so replacing $k$ by $K$ again lowers the target. The Taylor coefficients are $K(-Ky)^{r}/(r + 2)!$, giving case 2.

The exact tail after degree 41 is bounded by $< 10^{-6}$ using scalar estimates and $K \le 6$ in every parameter box (including $J = \infty$). $\square$

### 7.6 The bounded-range certificate

Combine the Bernstein rows with the target support:
$$\Phi_{i}(x; k) = \sum_{r = 0}^{40} \left( l_{i, r} q(k) + \delta_{i} a_{r}^{*} \right) x^{r}.$$
Its Bernstein coefficients are
$$L_{i, d} q(k) + \delta_{i} \sum_{r = 0}^{d} \omega_{d, r} a_{r}^{*}.$$

**Certification condition.** For every parameter box, every half-gap center $m \in \mathcal N_{0} \cap [0, 88]$, both signed half-gaps at each positive center, and every $d \in \{0, \dots, 40\}$, one has
$$L_{i, d} q(k) + \delta_{i} \sum_{r = 0}^{d} \omega_{d, r} a_{r}^{*} > 0.001.$$
Here $q(k)$ is evaluated at the box lower bound $[l]_{K, J}$ (Section 7.4).

This is a finite rational check: $37310$ separate inequalities, one per combination of (parameter box, half-gap, signed direction, $d$). The verification script `verification/numeric_balls.py` runs them.

**Conclusion of bounded-range signs.** When the certification condition holds, the Bernstein polynomial is $> 0.001$ on $[0, 1]$. Apply the three small-error bounds:
$$\mathcal D_{m, \nu} F_{i, \mathrm{slack}}(m + yx) > .001 - .00011 - .000001 - .000001 = .000888 > 0.$$
So $F_{i, \mathrm{slack}} > 0$ on $[0, 89.5]$ off nodes, and $= 0$ at nodes. Done for the bounded range.

### 7.7 The tail: quadratic barrier for the sine product

Beyond $s_{0} = 89.5$, the number of nodes is infinite, so finite-row Bernstein certification cannot work. Instead, the paper uses the following *barrier*.

**Lemma (Quadratic barrier).** For every real $s$ and every nearest point $m \in \mathcal Z$ of the full periodic zero set (so $\lvert s - m \rvert = \min_{n \in \mathcal Z} \lvert s - n \rvert$),
$$P(s) \ge 0.68 (s - m)^{2}.$$

**Sketch.** $P$ is periodic with period $12$, has double zeros at $\mathcal Z$, and no other real zeros. So it suffices to check the inequality on half-gaps. At each half-gap midpoint the exact values are computable:

| gap mod 12 | $P(s)/(s-m)^{2}$ at its midpoint |
|---|---|
| $[0, 1], [3, 4]$ | $12 - 8\sqrt 2 > .6862 > .68$ |
| $[1, 3]$ | $1 > .68$ |
| $[4, 7], [9, 12]$ | $\frac{8}{9}(3 + \sqrt 2 + \sqrt 3 + \sqrt 6) > .68$ |
| $[7, 9]$ | $9 > .68$ |

Logarithmic concavity of $s \mapsto P(s)^{1/2}/(s - m)$ interpolates the endpoint values at nodes. $\square$

### 7.8 The tail signs

For $s \ge s_{0}$, split each $H_{i}$:
$$R_{i, 0}(s) = C_{i} + \sum_{n \in \mathcal N_{0},\; n \le 88} \left( \frac{c_{in}}{(s - n)^{2}} + \frac{d_{in}}{s - n} \right),$$
$$p_{i}(s) = P(s) R_{i, 0}(s) + p_{i, \text{high}}(s),$$
where $p_{i, \text{high}}$ collects nodes $n > 88$ only.

**Step 1: rational part is positive.** Combine the intermediate values from Section 7.6 (with $\min(0, \cdot)$ floors on the lower evaluations), together with the reference correction. One gets
$$\sigma_{i} R_{i, 0}(s) > .002 \qquad (s \ge s_{0}).$$
The denominators $s - n$ and $(s - n)^{2}$ for $n \le 88$ and $s \ge 89.5$ are all $\ge 1.5$, so all reciprocals are $\le 1$ in modulus.

**Step 2: remainder second derivative is small.** Define
$$Z_{i}(s) = \sigma_{i} \left( p_{i, \text{high}}(s) + K_{p_{3 - i}}(s) \right) + \delta_{i} T_{k}(s).$$
A four-term Taylor bound using the derivative estimates of Section 4.6 gives
$$\sup_{s \ge s_{0}} \lvert Z_{i}''(s) \rvert < .0002.$$
The four terms are: direct high contribution ($36\pi^{2} \cdot 2.4 \cdot 10^{-7}$), low-frequency transformed ($3(3\pi)^{2} \cdot 36 \cdot 5 \cdot e^{-.70 s_{0}}$), refined transformed ($1.25(1.32\pi)^{2}(2.54 \cdot 5 / s_{0} + .04) e^{-.187 s_{0}}$), and the target ($2.36 e^{-2.36(s_{0} - 1)}$).

**Step 3: Zed has exact double zeros at high nodes.** By the exact contacts, at every $m \in \mathcal N$, $m \ge 91$, the decomposition
$$F_{i, \mathrm{slack}}(s) = \sigma_{i} P(s) R_{i, 0}(s) + Z_{i}(s)$$
has $Z_{i}(m) = Z_{i}'(m) = 0$ (the $P \cdot R_{i, 0}$ term vanishes to order 2 at $m$, so it must absorb all of $F_{i, \mathrm{slack}}$'s jet there, but $F_{i, \mathrm{slack}}$ also vanishes to order 2 at $m$, so $Z_{i}$ does too).

**Step 4: Taylor remainder bound.** For $s \ge s_{0}$, let $m$ be the nearest point of $\mathcal Z$ with $m \ge 91$ (always possible in the tail by a nearest-neighbor pairing). Taylor's integral formula gives
$$\lvert Z_{i}(s) \rvert \le .0001 (s - m)^{2}.$$

**Step 5: sign.** Combine:
$$\sigma_{i} P(s) R_{i, 0}(s) \ge .68 \cdot .002\, (s - m)^{2} = .00136 (s - m)^{2}.$$
Since $.00136 > .0001$,
$$F_{i, \mathrm{slack}}(s) \ge .00136(s - m)^{2} - .0001(s - m)^{2} = .00126(s - m)^{2} > 0.$$
At a node, the exact zero gives $F_{i, \mathrm{slack}}(m) = 0$ (by construction). Done.

### 7.9 Completing the proof of Theorem 4.2

Combining Sections 7.6 and 7.8:
- On $[0, 89.5]$ off nodes: $F_{i, \mathrm{slack}} > 0$ (Bernstein certificate).
- At nodes: $F_{i, \mathrm{slack}} = 0$ (exact contact).
- On $[89.5, \infty)$: $F_{i, \mathrm{slack}} \ge 0$ (quadratic barrier).

So $H_{1} \le T_{k}$ and $H_{2} \ge 0$ at every $s \ge 0$ and every $k \ge 2.36$. The four contacts of Proposition 6.1 complete the proof of Theorem 4.2 (the normalized Gaussian pair).

---

## 8. From Gaussian pairs to energy comparisons

The normalized pair gave us one specific function $f(x) = e^{-\pi h b \lvert x \rvert^{2}} H_{1}(b \lvert x \rvert^{2})$ with $\widehat f \ge 0$ and $f \le$ (the Gaussian of scale $k$-dependent). We now show how *any* function $f$ with $f \le G_{\alpha}$ (some target Gaussian) and $\widehat f \ge 0$ translates into an energy lower bound.

### 8.1 The density linear program

**Proposition (density LP).** *Let $f : \mathbb R^{2} \to \mathbb R$ be an integrable even Schwartz function with $\widehat f \ge 0$ pointwise. For every locally finite configuration $C$ of centered disk density one,*
$$\liminf_{R \to \infty} \frac{1}{N_{R}} \sum_{\substack{x, y \in C_{R} \\ x \ne y}} f(x - y) \ge \widehat f(0) - f(0).$$
*If also $f \le g$ pointwise for a function $g : (0, \infty) \to [0, \infty)$, then the same bound holds with $g$ on the left:*

$$E_{g}(C) \ge \widehat f(0) - f(0).$$

**Sketch.** The proof is a Cauchy-Schwarz argument plus a background comparison, due to Cohn and Zhao (2014) and Cohn and de Courcy-Ireland (2018).

*Notation.* Let $\mu_{R} = \sum_{x \in C_{R}} \delta_{x}$ be the point measure, and $\nu_{R}$ be the area measure restricted to $B_{(1 + \varepsilon)R}$, of total area $V_{R} = \pi(1 + \varepsilon)^{2} R^{2}$.

*Positivity.* $\widehat f \ge 0$ implies the sesquilinear form $Q_{f}(\mu, \nu) = \iint f(x - y) d\mu(x) d\nu(y)$ is positive semi-definite. So Cauchy-Schwarz applies:
$$Q_{f}(\mu_{R}, \mu_{R}) \cdot Q_{f}(\nu_{R}, \nu_{R}) \ge \lvert Q_{f}(\mu_{R}, \nu_{R}) \rvert^{2}.$$

*Cross-term estimate.* For each $x \in C_{R}$, $\int f(x - y) dv_{R}(y)$ equals $\int f$ plus an error bounded by $\int_{\lvert z \rvert \ge \varepsilon R} \lvert f \rvert$, because the disk of radius $\varepsilon R$ around $x$ is contained in $B_{(1 + \varepsilon)R}$. So
$$Q_{f}(\mu_{R}, \nu_{R}) = N_{R} (\widehat f(0) + o(1))$$
as $R \to \infty$, uniformly in $C$ through the boundary behavior.

*Self-term.* $Q_{f}(\nu_{R}, \nu_{R}) / V_{R} \to \widehat f(0)$, by dominated convergence.

*Diagonal.* $Q_{f}(\mu_{R}, \mu_{R}) = N_{R} f(0) + \sum_{x \ne y \in C_{R}} f(x - y)$.

*Combining.*
$$\frac{Q_{f}(\mu_{R}, \mu_{R})}{N_{R}} \ge \frac{(N_{R} (\widehat f(0) + o(1)))^{2}}{N_{R} \cdot Q_{f}(\nu_{R}, \nu_{R})} = \frac{N_{R}}{V_{R}} \cdot \frac{(\widehat f(0) + o(1))^{2}}{\widehat f(0) + o(1)}.$$
$N_{R}/V_{R} \to (1 + \varepsilon)^{-2}$. Taking $\liminf$ in $R$ and then $\varepsilon \downarrow 0$:
$$\liminf_{R} \frac{Q_{f}(\mu_{R}, \mu_{R})}{N_{R}} \ge \widehat f(0).$$
Subtracting $f(0)$ gives the claimed bound on the off-diagonal sum. $\square$

**The crucial role of "density one".** Without density $\approx 1$, $N_{R}/V_{R}$ can be 0 or $\infty$, and the bound is vacuous. Centered disk density one is what makes the LP close.

### 8.2 Sharp pairs for every Gaussian parameter

The theorem of Section 4.2 gives the normalized Gaussian pair for $k \ge 2.36$. This corresponds to a specific Gaussian target:
$$k = \pi(\alpha / b - h) \qquad \Longrightarrow \qquad \alpha = b(k / \pi + h).$$
So the target $T_{k}$ gives a "$G_{\alpha}$ minorant" for every $\alpha \ge \alpha_{\min} = b(2.36/\pi + 2/5) \approx 1$.

For $0 < \alpha < \alpha_{\min}$, the paper uses **Fourier complementation** (Cohn-Miller duality): the already-constructed large-$\beta$ witness $F_{\beta}$ at $\beta = 1/\alpha > 1$ supplies the small-$\alpha$ minorant $f_{\alpha} = G_{\alpha} - \widehat F_{\beta}$.

**Lemma (Gaussian duality).** *For every $\alpha > 0$, there is a real radial Schwartz function $f_{\alpha}$ satisfying*
$$f_{\alpha}(x) \le e^{-\pi \alpha \lvert x \rvert^{2}} \text{ pointwise}, \qquad \widehat f_{\alpha} \ge 0 \text{ pointwise},$$
$$f_{\alpha}(a) = G_{\alpha}(a) \text{ for every } a \in A \setminus \{0\}, \qquad \widehat f_{\alpha}(w) = 0 \text{ for every } w \in A^{*} \setminus \{0\}.$$

**Proof sketch for $\alpha \ge \alpha_{\min}$.** Let $k = \pi(\alpha / b - h)$ (so $k \ge 2.36$), let $f$ be the function of Theorem 4.2, and set $f_{\alpha} = k e^{-k} f$. This scales $f$ by a positive factor, so $\widehat f_{\alpha} \ge 0$ still. The pointwise bound: with $s = b \lvert x \rvert^{2}$,
$$f_{\alpha}(x) = k e^{-k} e^{-\pi h s} H_{1}(s) \le k e^{-k} e^{-\pi h s} T_{k}(s) = e^{-(k + \pi h) s} = e^{-\pi \alpha \lvert x \rvert^{2}} = G_{\alpha}(x).$$

The direct contacts come from the exact jets in Theorem 4.2, applied at every $n \in \mathcal N$, which contains every nonzero shell of $A$. The dual contacts come from the same contacts applied at every $n \in \mathcal N$, using that $A^{*} = JA$ has the same shells as $A$ (Section 1.1).

**Proof sketch for $0 < \alpha < \alpha_{\min}$ (Fourier complementation).** Set $\beta = 1/\alpha > 1/\alpha_{\min}$. Let $F_{\beta}$ be the large-$\beta$ witness constructed above:
$$F_{\beta} \le G_{\beta}, \quad \widehat F_{\beta} \ge 0, \quad F_{\beta}(a) = G_{\beta}(a) \; (a \ne 0), \quad \widehat F_{\beta}(w) = 0 \; (w \ne 0).$$

Now define
$$f_{\alpha} := G_{\alpha} - \widehat F_{\beta}.$$

This is the Fourier complement: $f_{\alpha}$ is the Fourier transform of $G_{\beta} - F_{\beta}$ (using $\widehat{G_{\alpha}} = G_{1/\alpha} = G_{\beta}$ and $\widehat{\widehat F}_{\beta} = F_{\beta}$). Verify the four required properties.

*$f_{\alpha} \le G_{\alpha}$.* By $\widehat F_{\beta} \ge 0$, pointwise $G_{\alpha} - \widehat F_{\beta} \le G_{\alpha}$. Done.

*$\widehat f_{\alpha} \ge 0$.* By Fourier duality, $\widehat f_{\alpha} = \widehat{G_{\alpha}} - \widehat{\widehat F}_{\beta} = G_{\beta} - F_{\beta}$. By $F_{\beta} \le G_{\beta}$, this is nonnegative. Done.

*Direct contacts.* For $a \in A \setminus \{0\}$:
$$f_{\alpha}(a) = G_{\alpha}(a) - \widehat F_{\beta}(a).$$
We need $\widehat F_{\beta}(a) = 0$. Since $F_{\beta}$ is radial, $\widehat F_{\beta}$ is radial too, so $\widehat F_{\beta}(a)$ depends only on $\lvert a \rvert$. By hypothesis, $\widehat F_{\beta}(w) = 0$ for every $w \in A^{*} \setminus \{0\}$. Since $A^{*} = JA$ and $J$ is a rotation, $A$ and $A^{*}$ have the same set of radial shells. So $\widehat F_{\beta}$ vanishes on all those shells, including every $a \in A \setminus \{0\}$. Hence $f_{\alpha}(a) = G_{\alpha}(a)$. Done.

*Dual contacts.* For $w \in A^{*} \setminus \{0\}$:
$$\widehat f_{\alpha}(w) = G_{\beta}(w) - F_{\beta}(w).$$
By hypothesis, $F_{\beta}(w) = G_{\beta}(w)$ for every $w \in A \setminus \{0\}$. Since $A$ and $A^{*}$ have the same shells and $F_{\beta}$, $G_{\beta}$ are radial, $F_{\beta}(w) = G_{\beta}(w)$ holds for $w \in A^{*} \setminus \{0\}$ as well. So $\widehat f_{\alpha}(w) = 0$. Done.

So $f_{\alpha} = G_{\alpha} - \widehat F_{\beta}$ is the small-$\alpha$ witness. This is the Cohn-Miller duality in dimension two, with the specific technical point that the direct/dual contact transfer between $A$ and $A^{*}$ depends on their having the same *radial* shells.

**Why the same-shells fact is crucial.** The complement trick above transfers contacts from $A^{*}$ (dual, where $F_{\beta}$ has zero Fourier) to $A$ (direct, where $f_{\alpha}$ needs direct contacts). This transfer only works because both lattices are radial (so $\widehat F_{\beta}$ on $A$ equals $\widehat F_{\beta}$ on $A^{*}$ at equal radii) and their shells coincide. In a hypothetical scenario where $A$ and $A^{*}$ had different shells, the small-$\alpha$ case would require a genuinely different construction. $\square$

**Consolidation.** The lemma shows that for every $\alpha > 0$, there is a minorant $f_{\alpha} \le G_{\alpha}$ with $\widehat f_{\alpha} \ge 0$, and with exact contacts. This gives the density-LP inequality.

### 8.3 Poisson summation and the lattice value

The contacts are what make the LP bound *sharp*. Without contacts, the LP bound $\widehat f(0) - f(0)$ would be strictly less than the target lattice sum. With contacts at both $A$ (direct) and $A^{*}$ (dual), the LP bound equals the target.

**Lemma (Poisson and lattice value).** *For any even Schwartz $f : \mathbb R^{2} \to \mathbb R$ with $\widehat f$ absolutely integrable (automatic on Schwartz functions),*
$$\sum_{a \in A} f(a) = \sum_{w \in A^{*}} \widehat f(w).$$

**Consequence for $f_{\alpha}$.** Write the sum splitting off the origin:
$$f_{\alpha}(0) + \sum_{a \ne 0} f_{\alpha}(a) = \widehat f_{\alpha}(0) + \sum_{w \ne 0} \widehat f_{\alpha}(w).$$
The dual contacts give $\sum_{w \ne 0} \widehat f_{\alpha}(w) = 0$. So
$$\sum_{a \ne 0} f_{\alpha}(a) = \widehat f_{\alpha}(0) - f_{\alpha}(0).$$
By the direct contacts $f_{\alpha}(a) = G_{\alpha}(a)$ at every nonzero $a \in A$:
$$\sum_{a \ne 0} G_{\alpha}(a) = \widehat f_{\alpha}(0) - f_{\alpha}(0).$$

So the LP bound is exactly the lattice sum. $\square$

The proof of the paper's Poisson summation lemma: the periodization $p_{f}(x) = \sum_{a \in A} f(x + a)$ is well-defined and smooth by Schwartz decay. Its Fourier coefficients over the torus $\mathbb R^{2}/A$ are exactly $\widehat f(w)$ for $w \in A^{*}$. Evaluating the Fourier series at $x = 0$ gives the identity.

**Only the lattice is periodized.** No periodicity is assumed of the *competitor* $C$; the LP applies to arbitrary locally finite configurations. This is central to the significance of the result. Periodization is used only to establish the target (LHS vs RHS) identity, not to restrict the competitor.

### 8.4 The Gaussian energy comparison

Combining:
1. The LP with $f_{\alpha}$ (Section 8.1),
2. The minorant $f_{\alpha} \le G_{\alpha}$ (Section 8.2),
3. The lattice value identity $\widehat f_{\alpha}(0) - f_{\alpha}(0) = \sum_{a \ne 0} G_{\alpha}(a)$ (Section 8.3),

yields the following.

**Proposition (Gaussian energy comparison).** *For every $\alpha > 0$ and every locally finite configuration $C$ of centered disk density one,*
$$\liminf_{R \to \infty} \frac{1}{N_{R}} \sum_{\substack{x, y \in C_{R} \\ x \ne y}} e^{-\pi \alpha \lvert x - y \rvert^{2}} \ge \sum_{a \in A \setminus \{0\}} e^{-\pi \alpha \lvert a \rvert^{2}}.$$

Equivalently, $E_{G_{\alpha}}(C) \ge E_{G_{\alpha}}(A)$ for every $\alpha > 0$.

This is the **Gaussian case** of the main theorem, proved in full. The extension to all completely monotone kernels takes the rest of Section 9.

### 8.5 The role of each ingredient

The four ingredients of the Gaussian case, and what each contributes:

| Step | Uses | Gives |
|---|---|---|
| LP (8.1) | $\widehat f \ge 0$, density one | Energy lower bound via $\widehat f(0) - f(0)$ |
| Minorant (8.2) | Exact contacts of $f_{\alpha}$ | $\widehat f_{\alpha}(0) - f_{\alpha}(0)$ matches lattice sum |
| Poisson (8.3) | Periodization of lattice | Identity (target-side equality) |
| Contacts (8.2, 6.1) | Interpolation proof | Makes the LP bound *sharp*, not just a lower bound |

None of these pieces is optional. The LP would give a lower bound without contacts; the contacts make the bound tight. The interpolation machinery (Sections 4-7) is what enables the contacts. The Fourier structure (Section 4) is what enables the interpolation. Each is necessary.

---

## 9. From Gaussian to universal: positive mixtures

We now reduce the general completely monotone kernel $g$ to the Gaussian case, by representing $g$ (or rather $g(\varepsilon + \cdot)$) as a positive mixture of Gaussians. This is the classical Hausdorff-Bernstein-Widder theorem.

### 9.1 The positive mixture representation

**Theorem (Hausdorff-Bernstein-Widder).** *Let $F : [0, \infty) \to [0, \infty)$ be smooth and completely monotone: $(-1)^{j} F^{(j)}(t) \ge 0$ for every $j \ge 0$ and every $t \ge 0$. Then there is a finite positive Borel measure $\rho$ on $[0, 1]$ such that*
$$F(t) = \int_{0}^{1} v^{t}\, d\rho(v) \qquad (t \ge 0),$$
*with $\rho(\{0\}) = 0$ (no atom at the excluded endpoint $v=0$).*

The intuition: each $v^{t} = e^{t \log v}$ is itself completely monotone, and the completely monotone cone is exactly the positive cone spanned by exponentials with negative rate (i.e. by $v^{t}$ for $0 < v \le 1$).

**Proof sketch.** Three steps.

*Step 1: discrete moments.* Fix $d > 0$. Repeated integration by parts gives $(-1)^{r} \Delta^{r} F(j d) \ge 0$, where $\Delta$ is the forward difference operator. Define a linear functional $\Lambda$ on polynomials by $\Lambda(x^{j}) = F(j d)$; positivity gives $\Lambda$ a representing measure $\eta_{d}$ on $[0, 1]$ with $\eta_{d}([0, 1]) = F(0)$.

*Step 2: refine to positive real exponents.* For $l \in \mathbb{N}$, push $\eta_{d}$ forward by $x \mapsto x^{l}$ (with $d = 1/l$). The result $\rho_{l}$ has mass $F(0)$ and $\int v^{j/l} d\rho_{l} = F(j/l)$ for $j \ge 0$. Take a weakly convergent subsequence; the limit $\rho = \lim_{n} \rho_{l_{n}}$ exists by compactness.

*Step 3: continuity in the exponent.* For fixed $t > 0$, let $t_{n} = \lceil l_{n} t \rceil / l_{n}$. Then $\lvert v^{t_{n}} - v^{t} \rvert \le (t_{n} - t)/(e t) \to 0$ uniformly in $v \in [0, 1]$. So continuity of $F$ forces $\int v^{t} d\rho = F(t)$.

Finally, the limit is finite on $(0, 1]$ but the mass at $0$ is zero (because $v^{t} \to \mathbf{1}_{(0, 1]}$ as $t \downarrow 0$, and $F(0^+) = F(0)$ gives $\rho((0, 1]) = F(0)$). $\square$

**Conversion to Gaussians.** Each $v^{t}$ for $v \in (0, 1]$ is
$$v^{t} = e^{t \log v} = e^{-\pi \alpha t}, \qquad \alpha = -\frac{\log v}{\pi},$$
which is exactly the Gaussian kernel $G_{\alpha}$. So
$$F(t) = \int_{0}^{1} G_{-\log(v)/\pi}(t)\, d\rho(v)$$
is a positive mixture of Gaussians.

### 9.2 Why we shift: singular kernels

The mixture above requires the representing measure $\rho$ to be finite, and this finite mass equals $F(0)$. For a completely monotone kernel that is singular at $0$, e.g. $g(t) = t^{-s/2}$ with $0 < s < 2$, $F(0) = \infty$ and the representing measure is infinite. So the direct mixture representation breaks down.

The paper circumvents this by *shifting* away from the singularity. Define, for $\varepsilon > 0$,
$$F_{\varepsilon}(t) = g(\varepsilon + t), \qquad t \ge 0.$$
This is smooth and completely monotone on $[0, \infty)$ (being the composition of a completely monotone function with a shift), and $F_{\varepsilon}(0) = g(\varepsilon) < \infty$ for every singular kernel of interest. So $F_{\varepsilon}$ has a finite positive representing measure $\rho_{\varepsilon}$ on $[0, 1]$, with total mass $g(\varepsilon)$.

The shift is removed as $\varepsilon \downarrow 0$, in the limit of the energy argument. This is critical: without the shift, we cannot represent $t^{-s/2}$ as a finite positive mixture of Gaussians. With the shift, we can, and the $\varepsilon \downarrow 0$ limit recovers the singular behavior.

### 9.3 The mixture energy argument

Fix a kernel $g : (0, \infty) \to [0, \infty)$ that is smooth and completely monotone. For $\varepsilon > 0$, let $\rho_{\varepsilon}$ be the representing measure of $F_{\varepsilon}(t) = g(\varepsilon + t)$, mass $g(\varepsilon)$.

For a configuration $C$ of centered disk density one, consider the shifted energy
$$E_{g(\varepsilon + \cdot)}(C) = \liminf_{R} \frac{1}{N_{R}} \sum_{\substack{x, y \in C_{R} \\ x \ne y}} g(\varepsilon + \lvert x - y \rvert^{2}).$$

A subtlety: the shift moves $\lvert x - y \rvert^{2}$ up by $\varepsilon$, which for a pair at distance zero would give $\varepsilon > 0$ instead of $0$. But zero-distance pairs are excluded from the off-diagonal sum; the shift is only felt at pairs with $\lvert x - y \rvert^{2} > 0$. For pairs at distance at least $\delta > 0$, the shifted energy is close to the unshifted one as $\varepsilon \to 0$. But pairs closer than $\delta$ can have large differences; that is why the limit $\varepsilon \downarrow 0$ must be taken *after* taking $R \to \infty$ and the \emph{lower limit} in $\liminf$.

**Key inequality (Fatou's lemma).** Write $E_{g(\varepsilon + \cdot)}(C)$ as an inner product against the positive mixture. Fatou's lemma makes the $\liminf$ of a positive mixture $\ge$ the mixture of $\liminf$s:
$$E_{g(\varepsilon + \cdot)}(C) \ge \int_{0}^{1} E_{G_{\alpha(v)}}(C)\, d\rho_{\varepsilon}(v)$$
where $\alpha(v) = -\log(v)/\pi$. The Gaussian case (Section 8) applies to each $E_{G_{\alpha}}$, giving
$$E_{G_{\alpha(v)}}(C) \ge E_{G_{\alpha(v)}}(A) = \sum_{a \ne 0} G_{\alpha(v)}(\lvert a \rvert^{2}).$$
So
$$E_{g(\varepsilon + \cdot)}(C) \ge \int_{0}^{1} \sum_{a \ne 0} G_{\alpha(v)}(\lvert a \rvert^{2})\, d\rho_{\varepsilon}(v).$$

Tonelli's theorem moves the sum out:
$$= \sum_{a \ne 0} \int_{0}^{1} G_{\alpha(v)}(\lvert a \rvert^{2})\, d\rho_{\varepsilon}(v).$$
Recognizing the inner integral as $F_{\varepsilon}(\lvert a \rvert^{2}) = g(\varepsilon + \lvert a \rvert^{2})$:
$$E_{g(\varepsilon + \cdot)}(C) \ge \sum_{a \ne 0} g(\varepsilon + \lvert a \rvert^{2}).$$

This is the shifted version of the main theorem.

### 9.4 Removing the shift

As $\varepsilon \downarrow 0$, the RHS $\sum_{a \ne 0} g(\varepsilon + \lvert a \rvert^{2})$ increases monotonically to $\sum_{a \ne 0} g(\lvert a \rvert^{2})$. The LHS $E_{g(\varepsilon + \cdot)}(C)$ requires care, but the paper proves that
$$\lim_{\varepsilon \downarrow 0} E_{g(\varepsilon + \cdot)}(C) = E_{g}(C)$$
for every locally finite $C$ of centered disk density one. The argument uses the singularity structure: for smooth completely monotone $g$, the singular behavior of $g$ at zero is dominated by the $\varepsilon \downarrow 0$ limit and can be controlled by the integral $\int_{0}^{1} v^{t} d\rho_{\varepsilon}$ against the pair density.

**Main Theorem (restated).** *Let $g : (0, \infty) \to [0, \infty)$ be smooth and completely monotone. For every locally finite configuration $C$ of centered disk density one,*
$$E_{g}(C) \ge E_{g}(A) = \sum_{a \in A \setminus \{0\}} g(\lvert a \rvert^{2}),$$
*interpreted in the extended nonnegative reals $[0, \infty]$.*

### 9.5 The infinite case

When $g(\lvert a \rvert^{2})$ is infinite for the first shell $a$ (e.g. $g(t) = t^{-1}$), the RHS is infinite. The theorem says: $E_{g}(C) = \infty$ for every density-one $C$. This is vacuous from the inequality direction but non-vacuous as a statement: it says that no configuration has finite energy for a strongly singular kernel, at density one.

The proof of this non-vacuous statement works by taking the limit $\varepsilon \downarrow 0$ explicitly: for singular $g$, the right side becomes infinite as $\varepsilon \downarrow 0$, and the LHS (which is $\ge$ the RHS by the mixture bound) is therefore also infinite.

### 9.6 Conclusion

We have established the full main theorem. To summarize the four-stage structure:

1. **Stage A (Section 4):** construct the entire family of candidate functions $p_{\mathbf{a}}$, $K_{p_{\mathbf{a}}}$, with explicit jet control.
2. **Stage B (Section 5):** reduce to finite data via $M = 168$, $N = 384$ quadrature points, and compute reference lists $X_{i} q(k)$.
3. **Stage C (Sections 6-7):** correct to exact interpolation and prove the two sign conditions.
4. **Transfer (Sections 8-9):** turn the Gaussian pair into an LP lower bound, then mixture into universal bounds.

Each stage individually is elementary or numerically hardened; their combination proves the main theorem.

---

## 10. Renormalized Riesz and jellium energies

The main theorem's Gaussian comparison extends to the logarithmic and Riesz potentials by *renormalization*, a systematic procedure that neutralizes the divergent parts of the raw energy sum. This section develops the renormalized theory and proves Corollary 2.3.

### 10.1 The field energy setup

Fix $0 \le s < 2$ once and for all. The interaction is $g(x) = \lvert x \rvert^{-s}$ for $0 < s < 2$ and $g(x) = -\log \lvert x \rvert$ at $s = 0$. The basic normalizations:
$$\kappa_{0} = 2\pi, \qquad \kappa_{s} = 4\pi \quad (0 < s < 2),$$
$$g_{0}(x) = -\log \lvert x \rvert, \qquad g_{s}(x) = \lvert x \rvert^{-s}, \qquad q_{0} = 2\pi, \qquad q_{s} = \frac{2^{2-s}\pi\Gamma(1-s/2)}{\Gamma(s/2)}.$$

**The weighted extension.** The potential $g_{s}$ in the plane lifts to a *weighted extension problem* in $\mathbb R^{3}$:
$$w_{s}(y) = \lvert y \rvert^{\gamma}, \qquad \gamma = s - 1 \text{ for } 0 < s < 2, \quad \gamma = 0 \text{ (so } w_{0} = 1) \text{ at } s = 0.$$
Points $(a, 0) \in \mathbb R^{2} \times \{0\}$ are identified with their base-plane position $a$; the "charges" live on the plane $y = 0$.

**Admissible fields.** A field is a gradient $E = \nabla h$ in $L^{p}_{\text{loc}}(\mathbb R^{3}; \mathbb R^{3})$ for some $p \in (1, 2)$ (or more generally $p$ in the range fixed at the start; see the paper), with $w_{s} E \in L^{1}_{\text{loc}}$, satisfying the divergence equation
$$-\operatorname{div}(w_{s} E) = \kappa_{s} \left( \sum_{p \in C} \delta_{(p, 0)} - \delta_{\mathbb R^{2}} \right).$$
Here the right side is the charge distribution: point charges at the particles, minus a uniform neutralizing background on the plane $y = 0$.

**Cutoff and truncated field.** For a self-energy cutoff $\eta > 0$, define
$$f_{s, \eta} = (g_{s} - g_{s}(\eta))_{+}, \qquad E_{\eta} = E - \sum_{p \in C} \nabla f_{s, \eta}(X - (p, 0)).$$
The function $f_{s, \eta}$ is the "capped fundamental solution": it equals $g_{s}$ outside the cap and is constantly $g_{s}(\eta)$ inside. Subtracting it makes $E_{\eta}$ finite-energy (removes the singular part of $E$ at each particle).

**Field energy.** For a centered square $K_{R} = [-R/2, R/2]^{2}$:
$$\mathcal W_{s, \eta}(E) = \limsup_{R \to \infty} \left[ \frac{1}{R^{2}} \int_{K_{R} \times \mathbb R^{k}} w_{s} \lvert E_{\eta} \rvert^{2} - \kappa_{s} g_{s}(\eta) \right].$$
A word about order of limits: for each fixed cutoff $\eta$, we take $R \to \infty$ first, then $\eta \downarrow 0$. The constant $\kappa_{s} g_{s}(\eta)$ subtracts the self-energy of one particle (of a specific cutoff).

$$\mathcal W_{s}(E) = \lim_{\eta \downarrow 0} \mathcal W_{s, \eta}(E), \qquad W_{s}(C) = \inf_{E \text{ compatible}} \mathcal W_{s}(E).$$

### 10.2 The periodic field and the Green function

Everything turns on the periodic picture. Fix a torus $T = \mathbb R^{2} / L$ of area $N \in \mathbb{N}$ (where $N = \lvert T \rvert$ is an integer), and $N$ distinct points $a_{1}, \dots, a_{N} \in T$.

**Canonical potential.** There is a mean-zero trace potential $H$ on $T \times \mathbb R^{k}$ satisfying
$$-\operatorname{div}(w_{s} \nabla H) = \kappa_{s} \left( \sum_{i = 1}^{N} \delta_{(a_{i}, 0)} - \delta_{T} \right).$$
Here $\delta_{T}$ denotes the uniform area measure on the base $T$ (the neutralizing background). The identity understood in the weak sense: test against gradients $\nabla \varphi$ with $\varphi$ compactly supported.

**Periodic energy.** Define
$$W_{s, L}(a_{1}, \dots, a_{N}) = \lim_{\eta \downarrow 0} \left[ \frac{1}{N} \int_{T \times \mathbb R^{k}} w_{s} \lvert \nabla H_{\eta} \rvert^{2} - \kappa_{s} g_{s}(\eta) \right],$$
with $H_{\eta}$ the truncated version of $H$ (subtracting $f_{s, \eta}$ at each particle).

**Green function representation.** Let $\mathcal G_{L}$ be the mean-zero unit-source extension Green function:
$$-\operatorname{div}(w_{s} \nabla \mathcal G_{L}) = \delta_{(0, 0)} - N^{-1} \delta_{T}.$$
Integration by parts gives the explicit formula
$$W_{s, L} = \kappa_{s}^{2} \left[ \frac{1}{N} \sum_{i \ne j} \mathcal G_{L}(a_{i} - a_{j}, 0) + \lim_{x \to 0} \left( \mathcal G_{L}(x, 0) - \frac{g_{s}(x)}{\kappa_{s}} \right) \right].$$
The first sum is off-diagonal (all pairs of distinct particles); the second is a limit that removes the diagonal singularity from the Green function's own self-contribution.

This follows the integration-by-parts argument of Petrache and Serfaty (2017, Proposition 1.5). In the log case ($s = 0$), $\mathcal G_{L}$ reduces to the ordinary 2D torus Green function.

**Comparison across period lattices.** A key issue for the corollary is: two periodic configurations $\mathcal C$ and $\mathcal D$ may not share a lattice, so their tori may have different sizes. The next lemma gives a comparison identity that works across arbitrary period lattices.

### 10.3 Heat-kernel comparison

Define $\Psi_{t}(x) = (4\pi t)^{-1} e^{-\lvert x \rvert^{2}/(4t)}$, the two-dimensional heat kernel. For a periodic configuration $\mathcal C = \bigcup_{i} (a_{i} + L)$ with $N = \lvert T \rvert$ classes per unit cell, set
$$P_{t, L}(x) = \sum_{v \in L} \Psi_{t}(x - v), \qquad H_{\mathcal C}(t) = \frac{1}{N} \sum_{i, j = 1}^{N} P_{t, L}(a_{i} - a_{j}) - \Psi_{t}(0).$$

$P_{t, L}$ is the heat-kernel Green function on $T$; $H_{\mathcal C}(t)$ is (up to constants) the "self-energy at heat time $t$."

**Lemma (Heat comparison across period lattices).** *Let $0 \le s < 2$, and let $\mathcal C$ and $\mathcal D$ be simple periodic configurations of density one. Then*
$$V_{s}(\mathcal C) - V_{s}(\mathcal D) = \frac{\kappa_{s} q_{s}}{\Gamma(1 - s/2)} \int_{0}^{\infty} \left( H_{\mathcal C}(t) - H_{\mathcal D}(t) \right) t^{-s/2}\, dt,$$
*where $V_{s}(\mathcal C) = W_{s, L}(a_{1}, \dots, a_{N})$ is the canonical periodic value of $\mathcal C$, and similarly for $\mathcal D$.*

**Proof sketch.** The heat representation of the Green function on any flat torus is
$$\mathcal G_{L}(x) = \frac{1}{\Gamma(\alpha)} \int_{0}^{\infty} \left( P_{t, L}(x) - \frac{1}{N} \right) t^{\alpha - 1}\, dt,$$
where $\alpha = 1 - s/2$ for $0 < s < 2$.

This is a spectral identity on the torus $T$; it applies to any flat torus. The individual $H_{\mathcal C}(t)$ are uniformly bounded: near $t = 0$ by $Ct^{-1} e^{-c/t}$ (from the Gaussian decay of $\Psi_{t}$ at nonzero image displacements), and near $t = \infty$ by $Ce^{-\lambda t}$ where $\lambda$ is the spectral gap of the torus Laplacian.

So the difference $H_{\mathcal C} - H_{\mathcal D}$ is absolutely integrable against $t^{-s/2}$ on $(0, \infty)$, and the identity holds. $\square$

**Why the comparison is useful.** The key point: $H_{\mathcal C}(t)$ is directly comparable to Gaussian energies! Setting $g_{t}(u) = e^{-u/(4t)}$ at heat time $t$, one has
$$H_{\mathcal C}(t) = \frac{1}{4 \pi t} E_{g_{t}}(\mathcal C).$$
So the Gaussian energy comparison (Section 8) applies: $H_{\mathcal C}(t) \ge H_{A}(t)$ for every $t > 0$ and every periodic $\mathcal C$ of density one. Apply the lemma:

$$V_{s}(\mathcal C) \ge V_{s}(A), \qquad \text{for every periodic } \mathcal C \text{ of density one.}$$

This is the periodic version of Corollary 2.3(i).

### 10.4 The infinite-system comparison

For a *nonperiodic* configuration $\mathcal C$ of density one, the periodic heat comparison does not directly apply: there is no period lattice. The paper extends the comparison via a *cubic approximation*: replace $\mathcal C$ by a square-periodic configuration on large cubes, prove the comparison there, and let the cube grow.

**Sketch.** Fix a large cube $\Lambda_{R} = [-R, R]^{2}$. Let $T_{R} = \mathbb R^{2} / (2R \mathbb Z)^{2}$ be the torus of area $(2R)^{2}$. Replace $\mathcal C$ by the periodic configuration $\mathcal C_{R}$ obtained by wrapping points into the torus. This introduces boundary errors of size $O(R)$ (from points near $\partial \Lambda_{R}$), but these vanish relative to $R^{2}$ in the field-energy normalization.

For each $R$, the periodic comparison gives $V_{s}(\mathcal C_{R}) \ge V_{s}(A)$. Take $R \to \infty$: on the left side, $V_{s}(\mathcal C_{R}) \to \mathcal W_{s}(E_{\mathcal C})$ for a canonical field $E_{\mathcal C}$ obtained from the periodic canonical potentials $H_{R}$ by a limit. On the right side, $V_{s}(A/R \mathbb Z) \to W_{s}(A)$ trivially (the lattice repeats).

Combining: $\mathcal W_{s}(E_{\mathcal C}) \ge W_{s}(A)$. Infimizing over $E_{\mathcal C}$ yields
$$W_{s}(\mathcal C) \ge W_{s}(A).$$

**Screening.** The field-energy lower bound ignores the possibility that $\mathcal C$ has "screening" structure: negative charge redistributions on the field side that improve the energy. The paper handles this by *screening*: construct a sequence of fields $E_{j}$ compatible with a density-one configuration $C_{j}$ whose renormalized energies approach $W_{s}(A)$ from above.

The mechanism: take the canonical field $E_{\mathcal C}$ from the periodic limit, and modify it in shells of thickness $\varepsilon R$ around each particle and near the cube boundary, to correct for the "screening defect"
$$V_{s}(\mathcal C) - \mathcal W_{\eta}(E^{\mathrm{can}}_{\mathcal C}) = \frac{\kappa}{\lvert T \rvert} \sum_{p \in T} \sum_{q \ne p} q_{\eta}(p - q) - 2\kappa F_{\eta}, \qquad F_{\eta} = \int_{\mathbb R^{2}} f_{\eta}(x, 0)\, dx.$$

The defect is nonnegative (stationarity), and vanishes as $\eta \downarrow 0$ for minimizers of $\mathcal W$. Screening provides an explicit approximating field whose energy is $V_{s}(\mathcal C)$ up to $o(1)$.

**Corollary (from Section 10.4).** *Under these definitions,*
$$W_{s}(A) = \inf_{\mathcal C \text{ density one}} W_{s}(\mathcal C),$$
*attained at the triangular lattice $A$.*

The full argument for screening (Petrache-Serfaty 2017, Section 6) is technical: it constructs local correction fields in each shell, uses weighted trace and Poincaré inequalities, and proves that the sum of shell contributions is $o(R^{2})$ in the field-energy normalization. The paper writes out the construction in full.

### 10.5 The periodic torus inequality

**Corollary (Torus inequality).** *For every integer $n \ge 1$, and every configuration of $n^{2}$ distinct points in $\mathbb R^{2}/(nA)$,*
$$V_{s}(a_{1}, \dots, a_{n^{2}}) \ge V_{s}(a_{1}^{\triangle}, \dots, a_{n^{2}}^{\triangle}),$$
*where $a_{i}^{\triangle}$ are the $n^{2}$ classes of $A/(nA)$ in the torus.*

**Proof.** Apply the heat comparison lemma (Section 10.3) with $\mathcal D = A$ and $\mathcal C$ the given torus configuration embedded periodically. The heat-kernel identity gives
$$V_{s}(\mathcal C) - V_{s}(A) = \frac{\kappa_{s} q_{s}}{\Gamma(1 - s/2)} \int_{0}^{\infty} (H_{\mathcal C}(t) - H_{A}(t)) t^{-s/2}\, dt.$$
The Gaussian energy comparison (Section 8) gives $H_{\mathcal C}(t) \ge H_{A}(t)$ for every $t > 0$, so the integral is $\ge 0$. $\square$

Note: repeating the argument with $nA$ instead of $A$ (and $H_{A}$ unchanged because $A$ is $nA$-periodic) shows that the inequality holds on the $nA$ torus as well. So the triangular lattice minimizes $W_{s}$ on every larger torus, not just the unit torus.

### 10.6 The jellium thermodynamic limit

The **jellium** is a finite system of $N$ ordered pairs of charges (points plus neutralizing background) on $\sqrt{N} \times \sqrt{N}$ squares, with the *ordered-pair* interaction
$$J_{s, R}(a_{1}, \dots, a_{N}) = \sum_{i \ne j} g_{s}(a_{i} - a_{j}) - 2 \sum_{i = 1}^{N} \int_{K_{R}} g_{s}(a_{i} - x)\, dx + \int_{K_{R}} \int_{K_{R}} g_{s}(x - x')\, dx\, dx'.$$

The three terms are: point-point interaction (excluding self), point-background interaction, and background-background interaction. All at scale $\sqrt{N}$ (so $N$ points in $N$ area, density one).

The jellium thermodynamic minimum is
$$e_{\mathrm{Jel}}^{(s)} = \lim_{N \to \infty} \frac{1}{N} \min_{a_{1}, \dots, a_{N} \in K_{\sqrt{N}}} J_{s, \sqrt{N}}(a_{1}, \dots, a_{N}).$$

**Corollary (Jellium minimum).** *For every $0 \le s < 2$,*
$$e_{\mathrm{Jel}}^{(s)} = \frac{W_{s}(A)}{\kappa_{s}}.$$

The proof has three parts.

**Part 1: Periodic comparison.** For each $N$, minimize $J$ over $N$ points in $K_{\sqrt{N}}$. Confinement: points near the boundary can be moved inward without increasing the energy, because the fundamental solution $g_{s}$ is monotone (attractive / monotone-positive). So the minimum is attained in a compact region and is finite.

**Part 2: Limit exists.** Using the periodic comparison $V_{s} \ge V_{s}(A)$ on the $nA$ torus, one bounds the min-energy from below by $N \cdot W_{s}(A) / \kappa_{s} - o(N)$. Using the periodic comparison plus a "widely spaced" approximation, one bounds it from above by the same thing. So the limit exists and equals $W_{s}(A)/\kappa_{s}$.

**Part 3: Identification.** The classical thermodynamic results of Lewin-Lieb-Seiringer (2019, Section VI, Theorem 2) and Lauritsen (2021, Theorem II.1) supply the necessary existence and equality statements for the thermodynamic limit of the jellium. Combined with the periodic comparison, they identify the limit as $W_{s}(A)/\kappa_{s}$.

This completes Corollary 2.3(ii).

### 10.7 Summary of Section 10

The renormalized/jellium results reduce to the Gaussian energy comparison via three moves:

1. **Heat representation.** Turn the renormalized energy $V_{s}$ into an integral over heat time of Gaussian energies at rate $1/(4\pi t)$, where the Gaussian case (Section 8) is available.
2. **Period lattice transfer.** Cross period lattices via the heat-kernel comparison lemma, using the fact that the heat kernel on any torus has the same spectrum but with different gaps.
3. **Screening + confinement.** Handle infinite configurations via a cubic approximation, and handle the jellium via confinement + periodic comparison + classical thermodynamic identification.

Each of these moves is standard in the crystallization literature; the paper's novelty is that the Gaussian comparison underlying move 1 (Section 8) is now *sharp*, thanks to the construction in Sections 4-7.

### 10.8 The shell-field decomposition

The "screening" move of Section 10.4 is a field-theoretic construction that must be fully specified in any complete proof. Here is the shape of the argument; the paper's Section ~6 gives the full details.

Fix the cutoff $\eta$ and a compatible field $E_\eta$ with $\mathcal W_\eta(E_\eta) \le m + \delta$ (where $m = W_s(A)$). The target is: replace $E_\eta$ by a screened field $\widehat E_\eta$ compatible with a *new* configuration with the same background and with energy no larger (up to $o(R^2)$).

The construction splits the ambient space $\mathbb R^{2+k}$ into five regions.

*Inner region.* $D_0 = K'_R \times [-\ell, \ell]^k$, where $K'_R$ is a concentric square inside $K_R$, its faces at distance $\varepsilon R$ from the outer boundary. Inside $D_0$, the old field $E_\eta$ is retained.

*Cap regions* (only $k = 1$). Above $D_0$: $K_R \times (\ell, R)$; below: $K_R \times (-R, -\ell)$.

*Shell region.* The lateral sides of the shell between $K'_R$ and $K_R$, subdivided into rectangular cells $H_i$ ordered so each has area $|H_i| \in [1/m_i, 1]$ for integer-labeled $m_i$. On $H_i$, integer $n_i = \lfloor m_i |H_i| \rfloor$ particles are placed.

*Outside shell.* Subdivided further into auxiliary "rooted tree" cells, each of side a fixed small width $b$, so that $\eta$-balls around retained centers distribute across these cells in a tree structure.

**The five component fields.**

- $E_\mathrm{cross}$: completes the retained smeared charges across the inner boundary. On each auxiliary tree, the assigned outward flux through the face leading to the root is $-\kappa$ times the tree's charge mass. This makes the cell problems locally compatible.
- $E_\mathrm{vert}$: a harmonic vertical correction $h_2(y) = \frac{\ell^\gamma}{1-\gamma} C^{\pm} |y|^{1-\gamma}$ in the caps, with source $-\ell^\gamma (C^+ + C^-) \delta_{\mathbb R^2}$.
- $E_\mathrm{side}$: on each shell cell $H_i$, a constant-outward-flux solution with source $\kappa(m_i - 1) + \ell^\gamma(C^+ + C^-)$ and datum $-E_\eta \cdot \nu_0 + j$ on the inner face. The $\kappa(m_i - 1)$ corrects for the "missing" background of the shell-cell.
- $E_\mathrm{round}$: a rounding correction with source $\kappa \sum_i \delta_i |H_i|^{-1} \delta_{\mathbb R^2}$ where $\delta_i = n_i - m_i |H_i|$ (the rounding error at each cell). Its $\sum \delta_i = 0$, so it vanishes as a global source.
- $E_\mathrm{new}$: the contributions of the $N_\mathrm{new}$ new particles placed at the barycenters of the $n_i$-decomposed rectangles, each with the truncated fundamental solution field around it.

Each component has an energy bound of the form

$$\int w |E_\mathrm{cross}|^2 \le C(1+g(\eta))(1+M_{R,\eta'}) \varepsilon^{-1} R,$$
$$\int w |E_\mathrm{vert}|^2 \le C e_{\varepsilon,R}\, \varepsilon^3 R^2,$$
$$\int w |E_\mathrm{side}|^2 \le C_{s,\eta'} \varepsilon (M_{R,\eta} + 1 + M_{R,\eta'}) R^2,$$
$$\int w |E_\mathrm{round}|^2 \le o(R^2),$$
$$\int w |E_\mathrm{new}|^2 \le C(1 + g(\eta)) N_\mathrm{new}, \qquad N_\mathrm{new} = O(\varepsilon R^2) + o(R^2).$$

Here $M_{R,\eta} = R^{-2} Q_\eta(E,R)$ is the truncated field energy per unit area. The constants depend on the fixed auxiliary cutoff $\eta'$ but not on $\eta$, $\varepsilon$, $R$, or the realization. Summing these inequalities gives the combined bound

$$R^{-2} \int w |\widehat E_\eta|^2 \le M_{R,\eta} + C_{s,\eta'} \varepsilon(1 + g(\eta))(1 + M_{R,\eta} + M_{R,\eta'}) + o_R(1).$$

The right side can be made as close to $M_{R,\eta}$ as desired (by choosing $\varepsilon$ small), and then $M_{R,\eta}$ is itself close to $m + \delta$ by assumption. Letting $\eta, \delta$ tend to zero gives $\mathcal W_s(\widehat E) = m$.

### 10.9 The two-cutoff defect

The elimination of "coincident particles" (two points at the same location) is a crucial step for making the minimizer a genuine point set. The paper does it via a two-cutoff defect.

Let $\alpha < \eta$ be two cutoffs. Since $E_\eta$ and $E_\alpha$ are both compatible fields (with the same source), they differ by a "cutoff difference":
$$E_\eta - E_\alpha = \sum_p \nabla(f_\eta - f_\alpha)(X - (p, 0)).$$

Stationarity of the natural field-energy functional under translations annihilates the divergence terms in a local energy identity. After separating self, distinct-labeled-pair, and background terms, one obtains

$$D_{\alpha,\eta} := a_\alpha - a_\eta + 2\kappa (F_\eta - F_\alpha) \ge 0,$$

where

- $a_\eta$ is the "mean truncated energy" at cutoff $\eta$, morally $\mathbb E_P[|E_\eta|^2 d\mathbb R^2]$;
- $F_\eta = \int_{\mathbb R^2} f_\eta(x, 0)\, dx$ is the self-energy of one charged particle at cutoff $\eta$;
- $D_{\alpha,\eta}$ counts the energy of the ordered pairs of truncated kernels $q_\eta(p - q)$ for $p \neq q$.

**Collision detection.** Suppose a coincident pair exists at some location, in some configuration contributing to the ergodic average with positive probability. Let $c_\mathrm{coll} = \mathbb E_P \sum_{p} N_p(N_p - 1)$ count the coincidence intensity (number of colliding pairs per unit area, in expectation). Then the coincident-pair part of $D_{\alpha,\eta}$ is $\kappa[g_s(\alpha) - g_s(\eta)] c_\mathrm{coll}$. Since $g_s(\alpha) \to \infty$ as $\alpha \downarrow 0$ for $0 < s < 2$, and the right side of $D_{\alpha,\eta} \ge 0$ stays bounded as $\alpha \downarrow 0$ (with only $\eta$ fixed), one gets $c_\mathrm{coll} = 0$.

**Simplicity.** So the minimizer is *simple*: no two particles coincide. This is exactly the property one wants to reason about point sets, not "clouds of charges at shared locations."

**Identifiability.** Once simplicity is established, the field minimizer $E_\eta$ and its "collision-free" limit $\lim_{\eta \downarrow 0} E_\eta$ can be treated as an actual gradient field compatible with a genuine point configuration. The minimizer is attained and its energy equals $m = W_s(A)$. The selection of a specific realization uses the ergodic theorem: the defect density $D_{\alpha,\eta}$ is a translation-invariant, integer-valued, mean-$D_{\alpha,\eta}$ spatial density; the ergodic theorem gives almost-sure convergence to its mean, and a countable selection of "vanishing-defect" cutoffs gives a single realization with the desired properties.

---

## 11. The arithmetic verification

The paper's numerical inequalities are not folklore; they are certified. Every strict constant comparison has either a rational proof or a machine-checked enclosure. This section describes both.

### 11.1 Rational enclosures of elementary constants

Every terminating decimal is exact rational. Set $\epsilon = 10^{-80}$. Let
$$p_{*} = 3.14159\ldots (80\ \text{digits}), \qquad b_{*} = 0.8660254\ldots (80\ \text{digits}).$$

**Lemma (Decimal enclosures).** *For $b = \sqrt{3}/2$ and $\pi$:*
$$p_{*} < \pi < p_{*} + \epsilon, \qquad b_{*} < b < b_{*} + \epsilon.$$

**Proof sketch for $\pi$.** Machin's identity: $\pi/4 = 4\arctan(1/5) - \arctan(1/239)$. For $0 < x < 1$, the alternating arctangent series has the finite-polynomial bound
$$A_{N}(x) < \arctan x < A_{N}(x) + x^{2N+1}/(2N+1),$$
where $A_{N}(x) = \sum_{j=0}^{N} (-1)^{j} x^{2j+1} / (2j+1)$.

Taking $N = 65$ gives two explicit rational inequalities that enclose $\pi$ within $p_{*} + \epsilon$ from above and $p_{*}$ from below. They are verified by clearing positive denominators in rational arithmetic.

**Proof sketch for $b$.** $b = \sqrt{3}/2 \iff b^{2} = 3/4$. Squaring the rational endpoints gives the pair of rational comparisons $b_{*}^{2} < 3/4 < (b_{*} + \epsilon)^{2}$. $\square$

For the Gaussian duality threshold $k \ge 2.36$, a shorter rational chain suffices:
$$\pi(2/\sqrt{3} - 2/5) > \frac{157}{50} \cdot \frac{1025838}{1732500} = \frac{1025838}{433250} > \frac{59}{25} = 2.36.$$

**Exponential enclosure.** For rational $a \ge 0$ and integer $q \ge 0$, let $P_{q}(a) = \sum_{j=0}^{q} a^{j}/j!$. If $a > 0$ and $q + 2 > a$, then
$$P_{q}(a) < e^{a} < P_{q}(a) + \frac{a^{q+1}}{(q+1)!} \cdot \frac{1}{1 - a/(q+2)}.$$
This is used with $q = 2000$ for scalar comparisons and with $a = 1$, $q = 6$ for $e < 2719/1000$.

### 11.2 The scalar error budgets

Eleven strict inequalities are the "arithmetic backbone" of the sign and correction bounds. Each has a rational proof via Section 11.1.

| Bound | Value |
|---|---|
| $U_{168}(5, .02) + 10^{-170}$ | $< 8.3 \cdot 10^{-14}$ |
| $U_{91}(5, .02) + 10^{-89}$ | $< 2.4 \cdot 10^{-7}$ |
| Target tail $c_{\mathrm{jet}} e^{-2.36 \cdot 168}/(1 - e^{-2.36})$ | $< 10^{-170}$ |
| Target tail $c_{\mathrm{jet}} e^{-2.36 \cdot 90}/(1 - e^{-2.36})$ | $< 10^{-89}$ |
| Half-gap Taylor tails $\mathcal R(m, Y, \nu)$ at four tuples | $< 10^{-8}$ each |
| Target-support Taylor tails (two cases) | $< 10^{-6}$ each |
| Second-derivative $\sup_{s \ge s_{0}} \lvert Z_{i}''(s) \rvert$ | $< .0002$ |
| Product endpoint $12 - 8\sqrt{2}$ | $> .68$ |
| Bernstein margin $.001 - .00011 - .000001 - .000001$ | $= .000888 > 0$ |
| Low-node rational margin $.0021 - 2 \cdot 10^{-8}$ | $> .002$ |
| Tail-barrier product $.68 \cdot .002$ | $> .0001$ |

Each of these is proved by a chain of rational inequalities, with $\pi$ and $\sqrt{3}$ replaced by the rational enclosures $p_{*}$, $b_{*}$ and, where needed, by $157/50 < \pi < 22/7$ and $\sqrt{3} < 1733/1000$.

### 11.3 The floor-rounded alternative

The paper also supplies a completely automatic rational procedure for verifying every finite matrix comparison, at the cost of introducing a rounding floor $\epsilon = 10^{-80}$ at every step of the arithmetic.

**Prescription (Appendix B, in brief).**

- `R(x)` = $\epsilon \cdot \lfloor x/\epsilon \rfloor$, computed separately on real and imaginary parts.
- `e_star(w)` computes $e^{w}$ via 12 squarings of the Taylor polynomial at $w/4096$, with `R(x)` applied at every arithmetic step.
- `W` (2311 × 171): mass matrix, exact rational + `R` floors. Computed once.
- `Ssharp`, `Ysharp`: 168 × 171 and 168 × 10, from `W` and jet projection, with `R` floors.
- `E1, E2`: 3 × 10, from the explicit constants of Section 5.5.
- `U, vp, vm`: power-sum doubling (9 stages), with `R` floors and strict bounds checked against the tables of Section 5.6.
- `X1, X2`: assembly into 171 × 10 reference lists.

The row-comparison procedure for signs (Section 7.6) becomes:

- For each parameter box, each half-gap center $m \in \mathcal N_{0} \cap [0, 88]$, both signed half-gaps at $m > 0$, both lists $i = 1, 2$, and $d \in \{0, \dots, 40\}$:
- Compute `L_{i, d}` as a rational row-vector `lsum = sum_{r=0}^d w_dr l_r` with $w_{dr} = \binom{d}{r}/\binom{40}{r}$, then evaluate `lsum` at the containing-box lower bound `[l]_{K, J}`.
- Compute `tsum` (target support) and add.
- Test `value > .001`.

Similarly, the tail-region check (Section 7.8) becomes: for each parameter box, each $i$, sum the sign-flipped low-node contributions and the reciprocal-scaled $d_{n}, c_{n}$ rows, and test `value > .002`.

**Total comparisons.** 37310 half-gap checks + 10 tail checks. All strict and rational.

**Roundoff budget.** From Appendix B:

- Exponential error: $< 10^{-68}$ (via a 12-stage squaring budget with amplification factor $3^{12} e^{7}$).
- Single-moment error: $< 10^{-48}$ after $41$-term composition.
- Matrix power errors $e_{r}$: $< 10^{-36}$ after 9 stages.
- Partial-sum error $v_{9}$: $< 10^{-34}$.
- Column error (X1, X2 assembled): $< 10^{-29}$.
- Bernstein row entry error: $< 10^{-16}$.
- Lower-evaluation error: $< 10^{-6}$ (94 terms × $2 \cdot 10^{-28}$).

Combined budget: the floor-rounded verification, if it passes, transfers to the looser exact-array bounds with margin. For example, $.00216 - 10^{-6} > .0021$ and $.00116 - 10^{-6} > .001$.

### 11.4 The Arb verifier

Two Python programs verify the finite matrix data at high precision:

- `verification/arithmetic_bounds.py`: uses 512-bit Arb balls and exact rational arithmetic, checks scalar budgets and error-propagation majorants.
- `verification/numeric_balls.py`: uses 256-bit precision, encloses quadrature arrays, checks the power table and the sign tables.

Both check every strict comparison in every table. If all pass, the finite construction is certified.

**What the verifier does not do.** It does not implement the floor-rounded procedure of Section 11.3 end-to-end; the floor-rounding is a separate (redundant) rational proof. The user can run either verifier, and independently run the floor-rounded procedure via a short implementation described in Appendix B.

### 11.5 What certification actually proves

The certifier checks finitely many strict comparisons. Combined with Sections 5.2 (quadrature error $< 10^{-31}$), 4.6 (tail bound $U_{r}$), and 6.3 (operator bounds), this recovers the exact infinite-frequency statement of Proposition 6.1 with the $2 \cdot 10^{-8}$ correction bound.

The four error sources in moving from finite to exact:

| Source | Bound |
|---|---|
| Quadrature → exact spectral | $10^{-31}$ |
| Power-sum truncation $D^{512}$ | $< 10^{-24}$ |
| Omitted finite nodes (target $T_{k}$ loss) | $< 8 \cdot 10^{-12}$ |
| Tail (nodes $n > 168$) | $< 8.3 \cdot 10^{-14}$ |

The largest is $8 \cdot 10^{-12}$; multiplied through the Schur complement ($\lVert \mathcal Q^{-1} \rVert \approx 1$), the finite solution extends to the exact infinite solution with correction $\lVert \ell - \ell^{\mathrm{ref}} \rVert < 2 \cdot 10^{-8}$.

---

## 12. Epilogue

### 12.1 What was proved

The theorem is now fully stated and proved:

**Theorem.** *Let $g : (0, \infty) \to [0, \infty)$ be smooth and completely monotone. For every locally finite configuration $C$ of centered disk density one in $\mathbb R^{2}$,*
$$E_{g}(C) \ge E_{g}(A) = \sum_{a \in A \setminus \{0\}} g(\lvert a \rvert^{2}),$$
*both sides interpreted in $[0, \infty]$. The triangular lattice $A$ is universally optimal.*

**Corollary.** *For $0 \le s < 2$, the triangular lattice minimizes the logarithmic ($s = 0$) and Riesz ($0 < s < 2$) energies with the standard renormalizations, on every torus $\mathbb R^{2}/(nA)$ and in the infinite system with unit volume density. The ordered-pair thermodynamic jellium minimum equals $W_{s}(A)/\kappa_{s}$.*

### 12.2 The proof structure, in one place

Five stages:

1. **Interpolation setup (Section 4).** Construct cardinal functions $p_{\mathbf{a}}$ (sine product times rational poles) with an explicit Fourier transform. Prove the transformed jet operator $S$ has bounded norm with exponentially decaying distant outputs.
2. **Finite truncation (Section 5).** Truncate at $168$ finite nodes and $N = 384$ quadrature points per frequency interval. Compute reference lists $X_{i} q(k)$ via 9-stage power-sum doubling.
3. **Exact correction (Section 6).** Show $X_{i} q(k)$ is within $2 \cdot 10^{-8}$ of the exact infinite interpolation solution, via a Schur complement contraction.
4. **Sign certification (Section 7).** On $[0, 89.5]$, use Bernstein-polynomial certificates on half-gaps (37310 rational inequalities). On $[89.5, \infty)$, use a quadratic barrier for the sine product plus a small second-derivative remainder.
5. **Energy transfer (Sections 8-10).** Gaussian comparison → density LP → positive mixtures → universal kernels → renormalized/jellium via heat-kernel representation.

Each stage is complete in itself. Their combination proves the theorem.

### 12.3 The Cohn-Kumar quantifier subtlety

The Cohn-Kumar conjecture asked for a single auxiliary function $f$, allowed to depend on $g$. This paper constructs *a* sharp auxiliary function $f_{\alpha}$ for each Gaussian scale $\alpha$, then integrates over the positive-mixture representation to obtain the energy bound for general $g$. It does not construct a single $f$ for a mixed kernel like $t^{-1/2}$.

This is a genuine weakening of the Cohn-Kumar per-potential quantifier. It is fully sufficient for the *energy* conclusion, which is the headline. It is not sufficient for stronger conclusions (e.g. positive-definiteness of some associated form).

The paper flags this honestly. The energy theorem is the positive result.

### 12.4 Relationship with the atomic companion

The atomic companion (OpenAI, *An atomic certificate for triangular-lattice universal optimality*, 2026) proves the same energy conclusion by an independent route: finite atomic columns at prescribed node sets, Gaussian damping $h = 17/50$, an auxiliary construction that omits node 24. Its energy conclusion *coincides* with the present paper's, but its auxiliary function pair does not: the node sets differ (the atomic version omits 24; the present version includes it), the damping constants differ ($17/50$ vs $2/5$), and the jet impositions differ. The theorem *statements* are the same; the *proofs* are not inter-derivable at the auxiliary level.

The Sharp-planar companion (OpenAI, *A sharp Fourier certificate for planar circle packing*, 2026) develops the general cardinal/Fourier framework for planar packing that this paper adapts to Gaussian targets. Sections 4-6 of the present paper are in dialogue with Sections 2-4 of the sharp-planar paper.

### 12.5 What's open

Three categories of open questions:

1. **Higher dimensions.** The Cohn-Kumar conjecture asks the same for $E_{8}$ (proved by CKMRV 2022) and the Leech lattice (proved). What about other lattices in intermediate dimensions? What about the hexagonal close packing in 3D and beyond? No general theorem yet.

2. **Stronger quantifiers.** A single auxiliary function for a general mixed kernel is not constructed. This is the Cohn-Kumar per-potential quantifier. Whether it is achievable, and whether it would give something strictly stronger than the energy conclusion, is open.

3. **Sharpness of the numerical constants.** The present proof uses $M = 168$, $N = 384$, $s_{0} = 89.5$, correction $2 \cdot 10^{-8}$, and several "generous" rational bounds. These could in principle be tightened substantially. Whether a *much* shorter certificate (e.g. $M = 24$, $N = 48$) exists is open.

### 12.6 One last summary for a high-school reader

The whole paper is the construction of one clever function $f$, or rather a family $f_{\alpha}$, with three properties:

- **Lower bound:** $f_{\alpha}$ sits *below* the Gaussian $G_{\alpha}$ everywhere.
- **Nonnegative Fourier:** $\widehat f_{\alpha}$ is nowhere negative.
- **Sharp contacts:** $f_{\alpha}$ equals $G_{\alpha}$ at every lattice point, $\widehat f_{\alpha}$ vanishes at every dual lattice point.

The three properties together make $f_{\alpha}$ a *certificate* that $A$ minimizes $G_{\alpha}$ energy against all comers. The rest of the paper is:

- Building $f_{\alpha}$ from "sine product times rational poles" (Sections 4-6).
- Certifying that the built function has the three properties (Section 7).
- Turning certifiicates into energy statements, and then generalizing from Gaussians to the whole completely monotone family (Sections 8-10).
- Verifying that the construction and certification are rigorous, without hand-waving (Section 11).

The triangular lattice's universal optimality follows. Nothing is missing.

---

## Appendix A: Where the numbers come from

The specific values $168$, $384$, $89.5$, and the constants in the tables are not chosen by magic; each is a threshold at which a rational tail estimate crosses a strict bound. The pattern is:

- $M = 168$: the point where $U_{168}$ (tail error at node $168$) crosses $8.3 \cdot 10^{-14}$ for the reference list size used. Any larger $M$ works too; $168$ is a convenient round choice.
- $N = 384 = 168 + 216$: the Fejér quadrature order at which $e^{.083 \pi M} < 10^{83}$ (the Cauchy-derivative bound). Roughly half of $M$; not tight.
- $s_{0} = 89.5$: the smallest $s$ where the "low-node" vs "high-node" cutoff is consistent with the reference-list bound $w < 5$. Determined by maximizing the test ratio $\mu \cdot \lambda / \sup \lvert Z_{i}'' \rvert > 1$ on the interval.
- Correction $2 \cdot 10^{-8}$: the sum of the finite-correction bound $1.37 \cdot 10^{-8}$ (Section 6.7) and the tail-correction bound $8.31 \cdot 10^{-14}$, rounded up.

Every choice is a matter of tradeoff between the analytic simplicity of the bound and the numerical tightness of the tables. None is sharp.

## Appendix B: Reference for the floor-rounded algorithm

The floor-rounded procedure (Section 11.3) is spelled out in full in the paper's Appendix B. Its relation to the exact finite construction is: Appendix B is a redundant rational path that does not go through Arb; a reader who runs it and the Arb verifier gets an independent double-check.

The floor-rounded procedure uses 37310 half-gap comparisons + 10 tail comparisons. Its roundoff budgets are as in Section 11.3. The floor-rounded algorithm is fully written out in mathematical pseudo-code; a machine implementation in any language (Python, Sage, Julia, Magma) should reproduce the comparison outputs within the rigorous bounds.

## Appendix C: Statement of every hypothesis used

For ease of reference, every numerical hypothesis invoked anywhere in the paper:

1. $b_{*} < \sqrt{3}/2 < b_{*} + 10^{-80}$, $p_{*} < \pi < p_{*} + 10^{-80}$ (Lemma 11.1).
2. $157/50 < \pi < 22/7$ (for shorter numerical checks).
3. $\sqrt{3} < 1733/1000$.
4. $e < 2719/1000$.
5. Table values: power bounds (Section 5.6), column bounds, sign bounds, half-gap bounds.
6. The $E_{1}, E_{2}$ specific values from Section 5.5 (extras).
7. $k \ge 2.36$ throughout the construction (from the Gaussian duality threshold).
8. Reference list size $\lVert X_{1} q \rVert + \lVert X_{2} q \rVert < 4.8$.
9. Exact list correction bound $\lVert \ell - \ell^{\mathrm{ref}} \rVert < 2 \cdot 10^{-8}$.
10. Extras constant bound $< .02$.
11. $s_{0} = 89.5$.

Everything else is proved from these.