Difference Bases: Verification Notebook

This notebook is a self-contained mathematical and computational companion to the Difference bases section in the paper.

The task template, including the prompt given to the agents, is available here.

The outline is as follows:

1. Problem Definitions
2. S1. Independent recovery of a record in the Leech–Golay family.
3. Additional findings
4. Appendix: the classical one-frequency lower bound

The required packages are:

None beyond the Python standard library.

1. Problem Definitions

The difference-basis problem asks how few integers are needed so that every distance from $1$ through $n$ occurs between two chosen integers. Its normalized limit measures how efficiently these distances can be covered as $n$ grows.

For a finite set $B\subset\mathbb Z$ , write

B-B=\{b-b':b,b'\in B\}.

The set $B$ is a difference basis for $\{1,\ldots,n\}$ if $\{1,\ldots,n\}\subseteq B-B$ . Let $\Delta(n)$ be the smallest possible cardinality of such a set and let

C=\inf_{n\ge1}\frac{\Delta(n)^2}{n}.

Rédei and Rényi (1949) proved that the sequence $\Delta(n)^2/n$ converges to this infimum. A finite certificate therefore consists of a list of distinct marks, its complete positive-difference set, and the first positive integer absent from that set.

A set $D\subset\mathbb Z_v$ is a cyclic difference cover if every residue modulo $v$ occurs as $d-d'$ for some $d,d'\in D$ . It is a planar difference set when every nonzero residue occurs exactly once as an ordered difference. This notebook verifies that property directly for the stored 90-point set.

The shared routines below are used by all verification cells. They authenticate and load the stored construction, form ordinary and cyclic differences, and identify the first missing positive difference.

Show code

Code cell 3 · In [1]

from __future__ import annotations

import ast
from bisect import bisect_left, bisect_right
import cmath
from collections import Counter
from fractions import Fraction
import hashlib
import json
import math
from pathlib import Path
import struct

DATA_PATH = Path("difference_basis_q89.npy")
METADATA_PATH = Path("construction_metadata.json")
PROVENANCE_PATH = Path("alphaevolve_provenance.json")
EXPECTED_NPY_SHA256 = "7b3032fa66098cf6af4c151bd7303366a05fd41cd442848c504322c0731f852a"
EXPECTED_BASIS_SHA256 = "498eaa873e07f473c282b00c8d7df51a706267d9438a99e1cef1c7abed1bc018"
EXPECTED_INNER_SHA256 = "8ccb0f4af4714bca24725a8acfb533c46b75cd752c3830c9c7e8308fc7df0ca4"
EXPECTED_PROVENANCE_SHA256 = "1daeee992a6346c22abed88665dc31b78dc131a85e0566abe1d0ad88d059ff47"


def sha256_bytes(payload: bytes) -> str:
    return hashlib.sha256(payload).hexdigest()


def canonical_list_sha256(values) -> str:
    payload = json.dumps(list(values), separators=(",", ":"), ensure_ascii=True).encode("ascii")
    return sha256_bytes(payload)


def load_npy_v1_int64_vector(path: Path) -> tuple[int, ...]:
    'Load precisely the bundled small no-pickle NumPy v1.0 format.'
    if not path.is_file():
        raise FileNotFoundError(f"required construction file is missing: {path}")
    payload = path.read_bytes()
    if sha256_bytes(payload) != EXPECTED_NPY_SHA256:
        raise ValueError("construction SHA-256 mismatch")
    if payload[:6] != b"\x93NUMPY" or payload[6:8] != bytes((1, 0)):
        raise ValueError("expected a NumPy v1.0 file")
    header_length = struct.unpack("<H", payload[8:10])[0]
    header_end = 10 + header_length
    header = ast.literal_eval(payload[10:header_end].decode("latin1").strip())
    if header.get("descr") != "<i8" or header.get("fortran_order") is not False:
        raise ValueError("expected a little-endian, C-order int64 array")
    shape = header.get("shape")
    if not isinstance(shape, tuple) or len(shape) != 1:
        raise ValueError("expected a one-dimensional array")
    count = int(shape[0])
    body = payload[header_end:]
    if len(body) != 8 * count:
        raise ValueError("array length does not match its header")
    return tuple(struct.unpack(f"<{count}q", body))


def positive_differences(values) -> set[int]:
    ordered = sorted(values)
    return {ordered[j] - ordered[i] for j in range(len(ordered)) for i in range(j)}


def first_missing_positive(differences: set[int]) -> int:
    candidate = 1
    while candidate in differences:
        candidate += 1
    return candidate


def cyclic_difference_counts(values, modulus: int) -> list[int]:
    counts = [0] * modulus
    for x in values:
        for y in values:
            counts[(x - y) % modulus] += 1
    return counts


if not METADATA_PATH.is_file():
    raise FileNotFoundError(f"required metadata file is missing: {METADATA_PATH}")
if not PROVENANCE_PATH.is_file():
    raise FileNotFoundError(f"required public provenance file is missing: {PROVENANCE_PATH}")
metadata = json.loads(METADATA_PATH.read_text(encoding="utf-8"))
provenance_payload = PROVENANCE_PATH.read_bytes()
if sha256_bytes(provenance_payload) != EXPECTED_PROVENANCE_SHA256:
    raise ValueError("public provenance SHA-256 mismatch")
public_provenance = json.loads(provenance_payload)
basis = load_npy_v1_int64_vector(DATA_PATH)

assert len(basis) == 360
assert tuple(sorted(set(basis))) == basis
assert canonical_list_sha256(basis) == EXPECTED_BASIS_SHA256
assert metadata["npy_sha256"] == EXPECTED_NPY_SHA256
assert metadata["basis_canonical_json_sha256"] == EXPECTED_BASIS_SHA256

print("Loaded and authenticated the 360-mark certificate.")
print(f"data SHA-256: {EXPECTED_NPY_SHA256}")

Saved output 1

Loaded and authenticated the 360-mark certificate.
data SHA-256: 7b3032fa66098cf6af4c151bd7303366a05fd41cd442848c504322c0731f852a

2. S1. Independent recovery of a record in the Leech–Golay family.

This section certifies the central finding in the main write-up.

Let

q=89,\qquad v=q^2+q+1=8011,\qquad A=\{0,1,4,6\}.

The stored basis has the form

B=\{va+d:a\in A,\ d\in D\},

where $D\subset\{0,\ldots,v-1\}$ has 90 elements.

Lemma 2.1 (product lift). Suppose $D\subset\mathbb Z_v$ is a cyclic difference cover and $A\subset\mathbb Z$ is a difference basis for $\{1,\ldots,s\}$ . For any representatives $\rho(D)\subset\{0,\ldots,v-1\}$ , the product

P=\{va+\rho(d):a\in A,\ d\in D\}

is a difference basis for $\{1,\ldots,sv\}$ .

Proof. Fix $1\le t\le sv$ and write $t=hv+r$ with $0\le r<v$ . If $r=0$ , choose an outer difference $h$ and use the same inner representative twice. If $r>0$ , then $0\le h\le s-1$ . A modular witness for $r$ has ordinary representative difference either $r$ or $r-v$ ; write it as $r-\varepsilon v$ with $\varepsilon\in\{0,1\}$ . The required outer difference is $h+\varepsilon\in\{0,\ldots,s\}$ , where difference zero is obtained by repeating one mark. The resulting product difference is

(r-\varepsilon v)+v(h+\varepsilon)=t.

This includes the endpoint $sv$ . $\square$

Lemma 2.2 (the boundary exposed by a cyclic gap). Let $D$ be planar modulo $v$ , and choose representatives so that the cut lies in a cyclic gap of length $G$ . Then the ordinary positive differences of the representatives contain $1,\ldots,G-1$ and omit $G$ .

Proof. Rotate the representatives so that the selected empty gap is the wraparound gap. Their ordinary span is $v-G$ . For $1\le r<G$ , the unique modular witness for $r$ cannot wrap: a wrapped witness would have ordinary separation $v-r>v-G$ , larger than the span. Hence $r$ is an ordinary positive difference. The two endpoints of the empty gap give the unique modular witness for residue $G$ , and that witness wraps across the cut. An ordinary difference $G$ would be a second ordered modular witness, contradicting planarity. $\square$

Proposition 2.3 (top-layer boundary). If $A\subset[0,s]$ has maximum difference $s$ , then for $1\le r<v$ ,

sv+r\in P-P \quad\Longleftrightarrow\quad r\in\rho(D)-\rho(D)

as an ordinary positive difference.

Proof. An inner difference lies strictly between $-v$ and $v$ . In a representation of $sv+r$ , any outer difference below $s$ would leave an inner term at least $v+r$ , which is impossible; an outer difference above $s$ is unavailable. Thus the outer difference must be $s$ and the inner difference must be $r$ . $\square$

Theorem 2.4 (the $q=89$ record). There exists a $360$ -element difference basis for $\{1,\ldots,49109\}$ . Consequently,

\Delta(49109)\le360, \qquad C\le\frac{360^2}{49109}=2.639027469506608\ldots.

The exhibited basis misses $49110$ .

Proof. The product $B$ constructed above is the required witness. The exact verification below establishes that $D$ is a $(8011,90,1)$ planar difference set with wraparound gap $1044$ , and that $B$ has $360$ distinct marks, covers every integer through $49109$ , and misses $49110$ . $\square$

The next cell authenticates the stored array, extracts its zero-th product layer, counts all $8011$ modular differences, reconstructs all four product layers, and enumerates every ordinary positive difference.

Show code

Code cell 5 · In [2]

q = 89
v = q * q + q + 1
outer_seed = (0, 1, 4, 6)

assert v == 8011
assert metadata["q"] == q and metadata["v"] == v
assert tuple(metadata["outer_seed"]) == outer_seed

# The zero-th product layer is exactly the inner set.
inner = tuple(x for x in basis if 0 <= x < v)
assert len(inner) == q + 1 == 90
assert canonical_list_sha256(inner) == EXPECTED_INNER_SHA256
assert all(0 <= x < v for x in inner)

# Verify the Singer/planar property without trusting a generator.
modular_counts = cyclic_difference_counts(inner, v)
assert modular_counts[0] == len(inner)
assert all(count == 1 for count in modular_counts[1:])

# Verify the outer ruler and reconstruct all product layers exactly.
outer_differences = positive_differences(outer_seed)
assert set(range(1, 7)).issubset(outer_differences)
reconstructed = tuple(sorted(v * a + d for a in outer_seed for d in inner))
assert reconstructed == basis
assert len(reconstructed) == len(outer_seed) * len(inner) == 360

# Gap and boundary checks.
cyclic_gaps = [inner[i + 1] - inner[i] for i in range(len(inner) - 1)]
cyclic_gaps.append(inner[0] + v - inner[-1])
largest_gap = max(cyclic_gaps)
inner_positive = positive_differences(inner)
inner_first_missing = first_missing_positive(inner_positive)
assert largest_gap == 1044
assert inner[0] + v - inner[-1] == 1044
assert inner_first_missing == 1044
assert set(range(1, 1044)).issubset(inner_positive)

# Complete ordinary-difference verification of the 360 marks.
basis_positive = positive_differences(basis)
first_missing = first_missing_positive(basis_positive)
assert set(range(1, 49110)).issubset(basis_positive)
assert first_missing == 49110
assert 49110 not in basis_positive

# Check the exact top-layer equivalence for every residue, not only the prefix.
assert all(((6 * v + r) in basis_positive) == (r in inner_positive) for r in range(1, v))

score = Fraction(len(basis) ** 2, first_missing - 1)
assert score == Fraction(129600, 49109)
assert basis[-1] - basis[0] == 55033

# Reconstruct the public notebook's product independently from its bundled snapshot.
assert public_provenance["commit"] == "8f447457957deac61e28bf1676746f0753b3b2f8"
assert public_provenance["notebook_path"] == "experiments/difference_bases/difference_bases.ipynb"
assert public_provenance["cell_index_zero_based"] == 1
assert public_provenance["cell_source_sha256"] == "09d0ae46f30326be4abc52dfa9b0f596d70cc6aadc1e19b3946dce6a19a13b8e"
public_outer_seed = tuple(public_provenance["outer_seed"])
public_inner = tuple(public_provenance["inner_set"])
public_basis = tuple(sorted(v * a + d for a in public_outer_seed for d in public_inner))
assert public_outer_seed == outer_seed
assert canonical_list_sha256(public_inner) == EXPECTED_INNER_SHA256
assert canonical_list_sha256(public_basis) == public_provenance["basis_canonical_json_sha256"]
assert public_basis == basis

print("Finite construction certificate: PASS")
print(f"q={q}, v={v}, |D|={len(inner)}, |B|={len(basis)}")
print(f"largest cyclic gap={largest_gap}, ordinary inner prefix={inner_first_missing - 1}")
print(f"covered interval=1..{first_missing - 1}, first missing={first_missing}")
print(f"span={basis[-1] - basis[0]}, exact score={score}, decimal={float(score):.15f}")
print(f"canonical public-list digest={canonical_list_sha256(basis)}")
print("public 360-entry list comparison: PASS")

Saved output 1

Finite construction certificate: PASS
q=89, v=8011, |D|=90, |B|=360
largest cyclic gap=1044, ordinary inner prefix=1043
covered interval=1..49109, first missing=49110
span=55033, exact score=129600/49109, decimal=2.639027469506608
canonical public-list digest=498eaa873e07f473c282b00c8d7df51a706267d9438a99e1cef1c7abed1bc018
public 360-entry list comparison: PASS

Verification provenance. The final digest uses the canonical byte string obtained by JSON-encoding the integer list with no spaces. The bundled provenance snapshot records the public seed and all 90 public inner marks from cell 1 of the commit-pinned AlphaEvolve notebook, together with that cell’s source digest. The cell above reconstructs the public 360-entry product from those marks and compares it entry-for-entry with the local certificate, rather than merely comparing scores.

Station Reference. The accepted 360-mark list and official score come from Eval #88. The product-lift and boundary endpoint proofs are in Archive #1, especially addendum archive_1-3. Evals #120, #132, and #180 independently checked the representative list, boundary prefix, and product mechanism.

Related Work. The Station result is an independent recovery of a known record, not a new result relative to the literature. Rédei and Rényi (1949) established the limiting constant, and Singer (1938) constructed planar cyclic difference sets. Leech (1956) and Golay (1972) combined the four-point basis $\{0,1,4,6\}$ with Singer difference sets to obtain earlier members of this construction family; see also Banakh and Gavrylkiv (2019). Bernshteyn and Tait (2019) record Golay’s preceding bound $C\le 2.6458\ldots$ . The $q=89$ member gives $C\le 2.6390274695\ldots$ , a new numerical upper bound within this classical family. The verification above establishes the Station construction’s entry-for-entry agreement with the construction reported by Georgiev et al. (2025).

3. Additional findings

These are additional results that may be interesting but are not included in the spotlight.

3.1 The Singer boundary advantage is asymptotically transient

For products built from planar Singer difference sets, we prove that the boundary extension is exactly one less than the largest cyclic gap. This identity explains the $q=89$ construction quantitatively: its score falls below $8/3$ precisely because its maximal gap exceeds $6q$ . We then use the flat Fourier spectrum of a Singer set to show that every such gap is $o(v)$ , uniformly over all multipliers. Consequently, the boundary contribution vanishes asymptotically and the four-layer Singer product converges to $8/3$ ; its finite record improvement cannot produce a smaller asymptotic constant.

For a cyclic set $E\subset\mathbb Z_v$ , define a cyclic gap of step $G$ to be the forward distance between two consecutive elements of $E$ ; it contains $G-1$ empty residues. Let

G_q=\max_{u\in\mathbb Z_v^\times}\operatorname{maxgap}(uD), \qquad R_q=G_q-1,

where $D$ is a planar difference set of order $q$ and multiplication is modulo $v=q^2+q+1$ .

Lemma 3.1 (gap identity). The optimized boundary run past $6v$ is exactly $R_q=G_q-1$ .

Proof. Fix a multiplier and place the cut inside a cyclic gap of step $G$ . Any witness arc crossing that cut must span the entire empty gap, so the unique witnesses for residues $1,\ldots,G-1$ do not cross it. The two endpoints of the gap give the unique witness for residue $G$ , and that witness does cross the cut. Thus this cut exposes exactly the ordinary prefix $1,\ldots,G-1$ . Maximizing first over cuts and then over unit multipliers proves the identity. $\square$

Theorem 3.2 (flat-spectrum gap bound). Uniformly over all unit multipliers,

G_q\le 1+\frac{v\sqrt q}{q+1}H_{(v-1)/2} =O(q^{3/2}\log q)=o(v),

where $H_m=\sum_{j=1}^m1/j$ .

Proof. For $f:\mathbb Z_v\to\mathbb C$ , use the Fourier transform

\widehat f(k)=\sum_{x\in\mathbb Z_v}f(x)e^{-2\pi i kx/v}.

The autocorrelation of the indicator $1_D$ equals $q+1$ at zero and $1$ at every nonzero residue. Taking its Fourier transform therefore gives the exact flat spectrum

|\widehat{1_D}(k)|^2=q\qquad(k\ne0).

Multiplication by a unit only permutes the nonzero frequencies, so the same identity holds for every $uD$ .

A gap of step $G$ contains an empty interval $J$ of length $L=G-1$ . Fourier inversion and the triangle inequality give

\frac{(q+1)L}{v} \le\frac{\sqrt q}{v}\sum_{k=1}^{v-1}|\widehat{1_J}(k)|.

The integer $v$ is odd. For $1\le k\le(v-1)/2$ , the geometric-sum formula and concavity of sine on $[0,\pi/2]$ give

|\widehat{1_J}(k)| =\frac{|\sin(\pi kL/v)|}{|\sin(\pi k/v)|} \le\frac{v}{2k}.

Pairing $k$ with $v-k$ yields

\sum_{k=1}^{v-1}|\widehat{1_J}(k)| \le vH_{(v-1)/2}.

Substitution gives $L\le v\sqrt q\,H_{(v-1)/2}/(q+1)$ and hence the displayed bound for $G$ . Since $v\asymp q^2$ and $H_m=O(\log m)$ , this is $O(q^{3/2}\log q)=o(v)$ . $\square$

Corollary 3.3 (vanishing finite dip). For the four-layer product with $A=\{0,1,4,6\}$ ,

|B|=4(q+1),\qquad n=6v+R_q,

and therefore

\frac{|B|^2}{n}\longrightarrow\frac83.

More precisely,

\frac83-\frac{16(q+1)^2}{6v+R_q} =\frac{\frac83R_q-16q}{6v+R_q}.

Proof. Lemma 3.1 gives the exact endpoint $6v+R_q$ , while Theorem 3.2 gives $R_q=o(v)$ . Substitution proves both the identity and the limit. $\square$

The sign of the finite dip is determined by whether $R_q>6q$ .

The next cell exhausts all unit multipliers at $q=89$ , checks the flat spectrum numerically, and verifies the exact finite-dip identity.

Show code

Code cell 9 · In [3]

current_gap_q89 = max(cyclic_gaps)
G_q89 = 0
gap_argmax_units = []
for unit in range(1, v):
    if math.gcd(unit, v) != 1:
        continue
    dilated = sorted((unit * mark) % v for mark in inner)
    unit_gaps = [dilated[i + 1] - dilated[i] for i in range(len(dilated) - 1)]
    unit_gaps.append(dilated[0] + v - dilated[-1])
    unit_gap = max(unit_gaps)
    if unit_gap > G_q89:
        G_q89 = unit_gap
        gap_argmax_units = [unit]
    elif unit_gap == G_q89:
        gap_argmax_units.append(unit)
R_q89 = G_q89 - 1
assert current_gap_q89 == G_q89 == 1044 and R_q89 == 1043
assert 6 * v + R_q89 == 49109

# Numerical evaluation of every nonzero Fourier magnitude; exact planarity was already checked.
max_flat_spectrum_error = 0.0
for frequency in range(1, v):
    transform = sum(
        cmath.exp(-2j * math.pi * frequency * mark / v) for mark in inner
    )
    max_flat_spectrum_error = max(
        max_flat_spectrum_error,
        abs(abs(transform) ** 2 - q),
    )
assert max_flat_spectrum_error < 1e-8

harmonic = math.fsum(1 / j for j in range(1, (v - 1) // 2 + 1))
explicit_gap_upper = 1 + v * math.sqrt(q) * harmonic / (q + 1)
assert G_q89 <= explicit_gap_upper

dip = Fraction(8, 3) - score
dip_identity = (Fraction(8, 3) * R_q89 - 16 * q) / (6 * v + R_q89)
assert dip == dip_identity
assert R_q89 > 6 * q

print("Singer boundary-transience finite checks: PASS")
print(f"G={G_q89}, R=G-1={R_q89}, exact dip={dip}={float(dip):.12f}")
print(f"exhaustive unit scan: {len(gap_argmax_units)} maximizing multipliers")
print(f"maximum flat-spectrum numerical error={max_flat_spectrum_error:.3e}")
print(f"explicit harmonic gap upper bound={explicit_gap_upper:.3f}")

Saved output 1

Singer boundary-transience finite checks: PASS
G=1044, R=G-1=1043, exact dip=4072/147327=0.027639197160
exhaustive unit scan: 6 maximizing multipliers
maximum flat-spectrum numerical error=8.210e-10
explicit harmonic gap upper bound=7451.614

For $q=89$ , the exact finite values are $G=1044$ , $R=1043$ , and $R>6q$ , explaining why this particular member lies below $8/3$ . The harmonic bound suffices to show uniformly that the boundary term is a vanishing fraction of $v$ .

Station Reference. The exact gap law comes from Archives #5 and #6. The score identity, flat-spectrum theorem, Dirichlet-kernel bound, and scope corrections come from Archive #7 and its addendum. Evals #132 and #180 give exhaustive finite checks of the gap identity, Eval #215 checks the exact dip identity, and Eval #228 checks the flat spectrum numerically. Question Room discussion #8 records the later asymptotic scope review.

Related Work. Lucia (2017), answering a question of Banakh, recorded the $O(v^{3/4}\log v)$ bound and noted that the logarithm can be removed. The boundary-gap relation belongs to the classical Leech–Golay analysis; the theorem above gives a self-contained derivation and applies it to the exact score.

4. Appendix: the classical one-frequency lower bound

This appendix records the classical lower-bound argument used for comparison.

Theorem 4.1 (one-frequency lower bound). Let $u_*$ be the unique solution of $\tan u=u$ in $(\pi,3\pi/2)$ . Then

C\ge K:=2\left(1-\frac{\sin u_*}{u_*}\right) =2.434467256422443\ldots.

Proof. Fix $n$ , and let $B$ be any difference basis for $\{1,\ldots,n\}$ with $k=|B|$ . Define

r(t)=\#\{(b,b')\in B^2:b-b'=t\}, \qquad F(\theta)=\sum_{b\in B}e^{ib\theta}.

The covering condition gives $r(t)\ge1$ for $1\le t\le n$ . Expanding the square and using $r(-t)=r(t)$ gives the exact finite identity

k^2-|F(\theta)|^2 =2\sum_{t\ge1}r(t)(1-\cos(t\theta)).

Every summand on the right is nonnegative. Since $|F(\theta)|^2\ge0$ ,

k^2\ge2\sum_{t=1}^{n}(1-\cos(t\theta)).

Set $\theta=u/n$ , where $u>0$ is fixed, and apply the inequality to a minimum basis of size $\Delta(n)$ . Dividing by $n$ and taking a limit inferior, the cosine sum is a Riemann sum, so

\liminf_{n\to\infty}\frac{\Delta(n)^2}{n} \ge2\left(1-\int_0^1\cos(ux)\,dx\right) =2\left(1-\frac{\sin u}{u}\right).

This holds for every $u>0$ . On $(\pi,3\pi/2)$ , the derivative of $\sin u/u$ vanishes exactly when $u\cos u-\sin u=0$ , equivalently $\tan u=u$ . The function $h(u)=u\cos u-\sin u$ satisfies $h'(u)=-u\sin u>0$ there, while $h(\pi)<0<h(3\pi/2)$ ; hence there is one critical point $u_*$ . The sinc function is positive on $(0,\pi)$ . On $(\pi,3\pi/2)$ its derivative has the sign of $h$ , so its value at the unique minimum $u_*$ is strictly below its endpoint value $-2/(3\pi)$ . For $u\ge3\pi/2$ the sinc function is at least $-1/u\ge-2/(3\pi)$ . Thus $u_*$ is its global minimum and gives the stated constant. Rédei–Rényi’s convergence theorem identifies the limit inferior with $C$ . $\square$

The code below locates $u_*$ by floating-point bisection of the monotone function $h$ and checks the constant.

Show code

Code cell 12 · In [4]

def stationary_numerator(u: float) -> float:
    return u * math.cos(u) - math.sin(u)


lo = math.pi
hi = 1.5 * math.pi
assert stationary_numerator(lo) < 0 < stationary_numerator(hi)

for _ in range(200):
    mid = (lo + hi) / 2
    if stationary_numerator(mid) < 0:
        lo = mid
    else:
        hi = mid

u_star = (lo + hi) / 2
sinc_u_star = math.sin(u_star) / u_star
K = 2 * (1 - sinc_u_star)

assert abs(math.tan(u_star) - u_star) < 1e-12
assert abs(sinc_u_star - math.cos(u_star)) < 1e-14
assert abs(K - 2.4344672564224433) < 1e-14

def finite_spectral_bound(n: int, u: float) -> float:
    return 2 * math.fsum(1 - math.cos(t * u / n) for t in range(1, n + 1)) / n

finite_rows = [(n, finite_spectral_bound(n, u_star)) for n in (10, 100, 1000, 10_000)]

print(f"u_*       = {u_star:.15f}")
print(f"sinc(u_*) = {sinc_u_star:.15f}")
print(f"K         = {K:.15f}")
print("finite Riemann-sum bounds:")
for n, value in finite_rows:
    print(f"  n={n:5d}: {value:.12f}")

Saved output 1

u_*       = 4.493409457909063
sinc(u_*) = -0.217233628211222
K         = 2.434467256422443
finite Riemann-sum bounds:
  n=   10: 2.548855725390
  n=  100: 2.446566488491
  n= 1000: 2.435683759033
  n=10000: 2.434588972475

The argument uses one consequence of spectral nonnegativity at one frequency. Averaging this same pointwise inequality over any probability distribution of frequencies cannot improve the result: an average of $2(1-\sin u/u)$ is at most its largest pointwise value $K$ . This obstruction is specific to averaging the one-frequency inequality; Yang and Liao (2022) obtain the stronger published bound $C>2.4421$ from higher-order Fourier information.

Related Work. The argument reproduces the classical one-frequency bound of Leech (1956).