Flat Autoconvolution: Verification Notebook

This notebook is a self-contained mathematical and computational companion to the Flat autoconvolution 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. Binary step functions preserve the unrestricted supremum.
3. Additional findings

The required packages are:

NumPy
python-flint
SymPy

1. Problem Definitions

The flat-autoconvolution problem asks how closely the self-convolution of a nonnegative function can resemble a flat-topped function: constant on its support and zero outside.

For a nonzero nonnegative function $f\in L^1(\mathbb R)\cap L^2(\mathbb R)$ , define

(f*f)(t)=\int_{\mathbb R}f(x)f(t-x)\,dx, \qquad Q(f)=\frac{\lVert f*f\rVert_2^2}{\lVert f*f\rVert_1\lVert f*f\rVert_\infty},

and let

C_{6.3}=\sup_{0\ne f\ge0}Q(f).

Hölder’s inequality gives $Q(f)\le1$ . For an arbitrary nonnegative output, equality holds precisely when that output is constant on its support. Whether an autoconvolution can approach this ideal closely enough that $C_{6.3}=1$ remains open. The functional was posed by Martin and O’Bryant (2009) and subsequently studied by Matolcsi and Vinuesa (2010). Without the nonnegativity constraint, De Dios Pont and Madrid (2021) proved that the sharp constant is $1$ .

Lemma 1.1 (exact score for equal-grid step functions).

Let $w_0,\ldots,w_{n-1}\ge0$ , let $h>0$ , and set

f_w(x)=\sum_{j=0}^{n-1}w_j\mathbf 1_{[jh,(j+1)h)}(x), \qquad c_k=\sum_{i+j=k}w_iw_j.

Then

Q(f_w)= \frac{2\sum_kc_k^2+\sum_kc_kc_{k+1}} {3\left(\sum_kc_k\right)\max_kc_k} \tag{1}

and $\sum_kc_k=(\sum_jw_j)^2$ .

Proof.

The autoconvolution is piecewise linear. Its interior knot values are $hc_k$ , and its two endpoint values are zero. A linear segment of width $h$ with endpoint values $ha$ and $hb$ contributes

\frac{h^3}{3}(a^2+ab+b^2)

to the squared $L^2$ -norm. Summing over the segments, including the two outer segments, gives

\lVert f_w*f_w\rVert_2^2 =\frac{h^3}{3}\left(2\sum_kc_k^2+\sum_kc_kc_{k+1}\right).

Nonnegativity and Tonelli’s theorem give

\lVert f_w*f_w\rVert_1 =h^2\sum_kc_k =h^2\left(\sum_jw_j\right)^2,

while piecewise linearity gives $\lVert f_w*f_w\rVert_\infty=h\max_kc_k$ . Substitution proves (1). Thus (1) exactly scores the continuous function; it is not quadrature or point sampling. $\square$

2. S1. Binary step functions preserve the unrestricted supremum.

A finite equal-grid step function has the form

f(x)=\sum_{j=0}^{m-1}v_j\mathbf1_{[x_0+jh,x_0+(j+1)h)}(x), \qquad v_j\ge0.

It is binary when every $v_j\in\{0,1\}$ ; equivalently, it is the indicator of a finite union of cells from one grid. Define

\begin{aligned} C&=\sup\{Q(f):0\ne f\ge0,\ f\in L^1(\mathbb R)\cap L^2(\mathbb R)\},\\ C_{\rm step}&=\sup\{Q(f):f\text{ is a nonzero finite equal-grid step function}\},\\ C_{01}&=\sup\{Q(f):f\text{ is a nonzero binary step function}\}. \end{aligned}

Theorem 2.1 (binary-step reduction).

C_{01}=C_{\rm step}=C.

More concretely, every weighted step function can be approximated in score by binary functions on finer grids. The equality $C_{01}=C_{\rm step}=C$ concerns unbounded resolution: it neither asserts equality at a fixed grid size nor gives an efficient rate for a finite search.

Lemma 2.2 (continuity of $Q$ ).

If $u_n,u\ge0$ , $u_n\to u$ in both $L^1$ and $L^2$ , and $u\ne0$ , then $Q(u_n)\to Q(u)$ .

Proof.

Put $d_n=u_n-u$ . Young’s inequality and Cauchy–Schwarz give

\begin{aligned} \lVert u_n*u_n-u*u\rVert_1 &\le(\lVert u_n\rVert_1+\lVert u\rVert_1)\lVert d_n\rVert_1,\\ \lVert u_n*u_n-u*u\rVert_2 &\le(\lVert u_n\rVert_2+\lVert u\rVert_2)\lVert d_n\rVert_1,\\ \lVert u_n*u_n-u*u\rVert_\infty &\le(\lVert u_n\rVert_2+\lVert u\rVert_2)\lVert d_n\rVert_2. \end{aligned}

All three differences tend to zero. The three norms defining $Q$ therefore converge, and the denominator at $u\ne0$ is positive. $\square$

Lemma 2.3 (simultaneous step approximation).

Every nonnegative $f\in L^1(\mathbb R)\cap L^2(\mathbb R)$ admits nonnegative finite equal-grid step functions $g_n$ with

\lVert g_n-f\rVert_1+\lVert g_n-f\rVert_2\longrightarrow0.

Proof.

First truncate $f$ in space and height. The bounded compactly supported functions $\min(f,M)\mathbf1_{[-R,R]}$ converge to $f$ in both norms as $R,M\to\infty$ . For a fixed truncation, replace the function on every grid cell by its cell average. Cell averaging is a contraction in both $L^1$ and $L^2$ by Jensen’s inequality. Continuous compactly supported functions are dense in both spaces, and their cell averages converge uniformly as the mesh tends to zero. A diagonal choice of truncation and mesh proves the lemma. $\square$

Lemma 2.4 (Hoeffding estimate).

If $X_1,\ldots,X_N$ are independent with $a_j\le X_j\le b_j$ , then

\Pr\!\left(\left|\sum_jX_j-\mathbb E\sum_jX_j\right|\ge s\right) \le2\exp\!\left(-\frac{2s^2}{\sum_j(b_j-a_j)^2}\right).

Proof.

Convexity of the exponential on $[a_j,b_j]$ gives $\mathbb E e^{\lambda(X_j-\mathbb EX_j)}\le e^{\lambda^2(b_j-a_j)^2/8}$ . Independence and Markov’s inequality give an upper-tail bound

\exp\!\left(-\lambda s+\frac{\lambda^2}{8}\sum_j(b_j-a_j)^2\right).

Optimizing in $\lambda$ , then applying the same argument to $-X_j$ , proves the two-sided estimate. $\square$

Proof of Theorem 2.1.

The score is invariant under multiplication by a positive constant, translation, and dilation. We may therefore begin with

g(x)=\sum_{i=0}^{m-1}v_i\mathbf1_{[i,i+1)}(x), \qquad0\le v_i\le1.

Divide each unit interval into $T$ microcells. For $0\le a<mT$ , put $p_a=v_{\lfloor a/T\rfloor}$ , choose independent Bernoulli variables $\xi_a$ with $\Pr(\xi_a=1)=p_a$ , and define the binary refinement

F_T(x)=\sum_{a=0}^{mT-1}\xi_a\mathbf1_{[a/T,(a+1)/T)}(x).

At output index $t$ , set

r_t=\sum_{a+b=t}\xi_a\xi_b, \qquad c_t=\sum_{a+b=t}p_ap_b.

The tent-function formula gives

(F_T*F_T)((t+1)/T)=r_t/T, \qquad (g*g)((t+1)/T)=c_t/T. \tag{2}

Off the diagonal, $\mathbb E(\xi_a\xi_b)=p_ap_b$ . At an even $t$ , the one possible diagonal term has expectation $p_a$ , rather than $p_a^2$ . Hence, uniformly in $t$ ,

0\le\mathbb Er_t-c_t\le\frac14.

For fixed $t$ , group the ordered off-diagonal terms into unordered pairs:

r_t=2\sum_{a<b,\ a+b=t}\xi_a\xi_b +\mathbf1_{\{t\ \operatorname{even}\}}\xi_{t/2}.

These pair variables are independent because the pairs $\{a,t-a\}$ are disjoint. Their squared range lengths sum to at most $2mT+1$ . Lemma 2.4 and a union bound over fewer than $2mT$ output indices give

\Pr\!\left( \max_t\left|\frac{r_t-c_t}{T}\right|> \eta+\frac1{4T}\right) \le4mT\exp\!\left(-\frac{2\eta^2T^2}{2mT+1}\right).

Choose $\eta=T^{-1/4}$ . The right side tends to zero, so for every sufficiently large $T$ at least one deterministic binary refinement satisfies

\max_t|r_t/T-c_t/T|\le T^{-1/4}+1/(4T).

Both convolutions in (2) are linear between consecutive grid points. The last inequality therefore implies

\lVert F_T*F_T-g*g\rVert_\infty\longrightarrow0.

Their supports lie in the fixed interval $[0,2m]$ , so the convergence also holds in $L^1$ and $L^2$ . Therefore $Q(F_T)\to Q(g)$ , and $C_{01}=C_{\rm step}$ . Lemma 2.3 and Lemma 2.2 then give $C_{\rm step}=C$ . $\square$

For a finite set $A\subset\mathbb Z$ , let $f_A(x)=\sum_{a\in A}\mathbf1_{[a,a+1)}(x)$ and let $r_A(k)=\#\{(a,b)\in A^2:a+b=k\}$ . Lemma 1.1 specializes to

Q(f_A)= \frac{2\sum_kr_A(k)^2+\sum_kr_A(k)r_A(k+1)} {3|A|^2\max_kr_A(k)}.

Thus the unrestricted continuous supremum can be studied through finite-set sum-representation functions.

Station Reference. The binary-step reduction originates in Archive #76. The reviewer-response addendum in message archive_76-3 supplies the uniform treatment of diagonal terms, and message archive_76-4 clarifies the unbounded-resolution scope. The finite-set scoring identity was recorded in Archive #27.

Related Work. Cilleruelo, Ruzsa, and Vinuesa (2010) used probabilistic discretization to transfer weighted profiles to sets for a mass-to-maximum functional, and Boyer and Li (2025) derived the finite equal-grid scoring identity. The uniform-profile refinement proving $C_{01}=C_{\rm step}=C$ appears novel.

3. Additional findings

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

3.1 Verified $100{,}000$ -step construction

Here we record the strongest construction found in the Station run: a nonnegative $100{,}000$ -step function satisfying

Q(f)=0.953189630304235\ldots.

We authenticate the complete coefficient vector, convert its binary64 entries to integers over a common dyadic denominator, and evaluate its continuous score by exact integer convolution. This is the explicit witness underlying the Station’s reported bound $C_{6.3}>0.953189$ .

Theorem 3.1 (a $100{,}000$ -step construction).

There exists a nonnegative $100{,}000$ -step function $f_w$ such that

Q(f_w)=0.953189630304235\ldots>0.953189.

Consequently, $C_{6.3}>0.953189$ .

Verification. The file autocorr_6-3_weights.npy contains the coefficients of this function. Every entry is an exact dyadic rational with denominator dividing $2^{79}$ . After multiplying the weights by $2^{79}$ , exact integer polynomial multiplication gives every convolution coefficient. The common scale cancels from (1), so an integer cross-multiplication proves the strict comparison with $0.953189$ .

The next cell authenticates the coefficient array, performs the exact integer convolution, and evaluates the rational score in (1).

Show code

Code cell 5 · In [1]

from pathlib import Path
from fractions import Fraction
from decimal import Decimal, localcontext
import hashlib
import numpy as np
from flint import fmpz_poly

DATA_FILE = Path("autocorr_6-3_weights.npy")
EXPECTED_SHA256 = "0d1f1f2a65b8db8cf7d97a4bca7d49d85a74be4c2f7639d7267cd29c655bdf83"
COMMON_DENOMINATOR = 1 << 79

assert DATA_FILE.is_file(), f"Required construction file is missing: {DATA_FILE}"
assert hashlib.sha256(DATA_FILE.read_bytes()).hexdigest() == EXPECTED_SHA256

weights = np.load(DATA_FILE, allow_pickle=False)
assert weights.shape == (100_000,)
assert weights.dtype == np.float64
assert np.all(np.isfinite(weights))
assert np.all(weights >= 0)
assert np.max(weights) > 0

integer_weights = []
for value in weights:
    numerator, denominator = float(value).as_integer_ratio()
    assert COMMON_DENOMINATOR % denominator == 0
    integer_weights.append(numerator * (COMMON_DENOMINATOR // denominator))

convolution = list(fmpz_poly(integer_weights) ** 2)
weight_sum = sum(integer_weights)
assert all(coefficient >= 0 for coefficient in convolution)
assert sum(convolution) == weight_sum**2

score_numerator = (
    2 * sum(coefficient**2 for coefficient in convolution)
    + sum(left * right for left, right in zip(convolution, convolution[1:]))
)
score_denominator = 3 * weight_sum**2 * max(convolution)
exact_score = Fraction(int(score_numerator), int(score_denominator))
assert exact_score > Fraction(953189, 1_000_000)

print(f"data SHA-256: {EXPECTED_SHA256}")
print(f"number of weights: {weights.size:,}")
print(f"nonzero weights: {np.count_nonzero(weights):,}")
with localcontext() as context:
    context.prec = 24
    exact_decimal = Decimal(exact_score.numerator) / Decimal(exact_score.denominator)
print(f"exact score: {exact_decimal}")
print("certified comparison: Q(f_w) > 0.953189")

Saved output 1

data SHA-256: 0d1f1f2a65b8db8cf7d97a4bca7d49d85a74be4c2f7639d7267cd29c655bdf83
number of weights: 100,000
nonzero weights: 32,683
exact score: 0.953189630304235556430284
certified comparison: Q(f_w) > 0.953189

Station Reference. The vector is row 0 of Research Center Evaluation #1942, equivalently Seed Bank key c000 for Evaluation #1942 attempt 1. Evaluation #1942 records the floating-point value $0.9531896303042352$ ; the exact verification above gives $0.953189630304235556\ldots$ . The exact equal-grid scoring identity was recorded in Archive #27.

Related Work. Georgiev et al. (2025) report $C_{6.3}\ge0.961$ ; their released $50{,}000$ -step construction has score $0.961021\ldots$ . Yüksekgönül et al. (2026) subsequently obtained $0.961206$ , and Ye et al. (2026) obtained $0.962694$ . The vector is an independently generated reproducible construction below these published records.

3.2 The sharp two-interval ceiling

We solve the first disconnected geometric class exactly. If $S$ is a union of at most two bounded intervals, we prove

Q(\mathbf 1_S)\le\frac{20}{27},

with equality only, up to translation, reflection, and dilation, when the shorter interval, intervening gap, and longer interval have proportions $1:1:2$ . The proof divides the two-interval parameter space into seven regions and certifies the required polynomial inequality on each one. Thus no indicator of a union of at most two intervals can have score above $20/27$ .

Theorem 3.2 (indicators of at most two intervals).

If $S\subset\mathbb R$ is a non-null union of at most two bounded intervals, then

Q(\mathbf1_S)\le\frac{20}{27}.

The constant is sharp. Apart from translation, reflection, and positive dilation, equality holds only for

S=[0,1]\cup[2,4],

so the short interval, the intervening gap, and the long interval have lengths in the ratio $1:1:2$ .

Proof.

The score is invariant under translation, reflection, and dilation. A single interval has a triangular autoconvolution and direct integration gives $Q=2/3<20/27$ . Otherwise, after merging overlapping interiors, relabel the two positive-length intervals and normalize their total length to one, so

S=[0,t]\cup[t+g,1+g], \qquad0<t\le\frac12,\quad g\ge0.

Put $b=1-t$ . Define the triangle and short–long overlap profiles

T_L(x)=\max\{0,\min(x,2L-x,L)\},

R_{t,b}(x)=\max\{0,\min(t,x)-\max(0,x-b)\}.

Then direct interval overlap gives

(\mathbf1_S*\mathbf1_S)(x) =T_t(x)+2R_{t,b}(x-t-g)+T_b(x-2t-2g). \tag{3}

Since $|S|=1$ , put $N(t,g)=\lVert\mathbf1_S*\mathbf1_S\rVert_2^2$ and $P(t,g)=\lVert\mathbf1_S*\mathbf1_S\rVert_\infty$ , so $Q=N/P$ . Splitting the piecewise-linear profile (3) at its breakpoints gives the following exhaustive table; boundaries follow by continuity.

Region	Parameter range	$\frac32N(t,g)$	$P(t,g)$
1	$0<t<1/3,\ 0<g<t$	$-g^3+6g^2t+6gt^2-6gt+1$	$1-g$
2	$0<t<1/3,\ t<g<1-2t$	$3g^2t+9gt^2-6gt-t^3+1$	$1-t$
3	$0<t<1/3,\ 1-2t<g<1-t$	$-g^3-3g^2t+3g^2-3gt^2+6gt-3g-9t^3+12t^2-6t+2$	$1-t$
4	$0<t<1/2,\ g>1-t$	$-8t^3+9t^2-3t+1$	$\max(1-t,2t)$
5	$1/3<t<1/2,\ 0<g<1-2t$	same as Region 1	$1-g$
6	$1/3<t<1/2,\ 1-2t<g<t$	$-2g^3+3g^2-6gt^2+6gt-3g-8t^3+12t^2-6t+2$	$2t$
7	$1/3<t<1/2,\ t<g<1-t$	same as Region 3	$2t$

For a linear segment of width $L$ and endpoint values $y_0,y_1$ , the contribution to $N$ is

\frac L3(y_0^2+y_0y_1+y_1^2).

It remains to show $D(t,g):=20P(t,g)-27N(t,g)\ge0$ . Each of Regions 1, 2, 3, 5, 6, and 7 is a triangle. In barycentric coordinates $u,v,w\ge0$ , $u+v+w=1$ , the second exact cell homogenizes $D$ to degree four in Regions 1–2 and degree three in the other four regions. Every coefficient is nonnegative, and exactly one pure-vertex coefficient vanishes. All other coefficients are positive. Hence $D\ge0$ , with equality only at $(t,g)=(1/3,1/3)$ .

Region 4 is even stricter. Its $N$ is independent of $g$ . If $0<t\le1/3$ , then $P=1-t$ and

19P-27N=(3t-1)(48t^2-38t-1)\ge0.

Both factors are nonpositive. If $1/3\le t\le1/2$ , then $P=2t$ and

19P-27N=2(3t-1)(24t^2-19t+9)\ge0;

the quadratic is positive because its discriminant is negative. Thus Region 4 has the stricter ceiling $19/27$ .

Finally, at $(t,g)=(1/3,1/3)$ , direct substitution gives $N=40/81$ , $P=2/3$ , and $Q=20/27$ . The equality conditions in the coefficient certificate give uniqueness modulo the stated symmetries. $\square$

Within each listed region, the breakpoint order and peak branch are constant; the boundary lines are precisely the breakpoint-collision and peak-switch lines.

The next cell derives the seven regional formulas from the interval-overlap profile and checks them at exact rational points.

Show code

Code cell 8 · In [2]

from fractions import Fraction as F
import sympy as sp


def triangle(length, x):
    return max(F(0), min(x, 2 * length - x, length))


def rectangle_overlap(short, long, x):
    return max(F(0), min(short, x) - max(F(0), x - long))


def profile_value(t0, g0, x):
    b0 = 1 - t0
    return (
        triangle(t0, x)
        + 2 * rectangle_overlap(t0, b0, x - t0 - g0)
        + triangle(b0, x - 2 * t0 - 2 * g0)
    )


def exact_N_P(t0, g0):
    b0 = 1 - t0
    knots = sorted({
        F(0), t0, 2*t0,
        t0+g0, 2*t0+g0, 1+g0, 1+t0+g0,
        2*t0+2*g0, 1+t0+2*g0, 2+2*g0,
    })
    values = [profile_value(t0, g0, x) for x in knots]
    N0 = sum(
        (x1-x0) * (y0*y0 + y0*y1 + y1*y1) / 3
        for x0, x1, y0, y1 in zip(knots[:-1], knots[1:], values[:-1], values[1:])
    )
    return N0, max(values), knots, values


t, g = sp.symbols("t g", real=True)
N_expr = {
    1: 2*(-g**3 + 6*g**2*t + 6*g*t**2 - 6*g*t + 1)/3,
    2: -2*(-3*g**2*t - 9*g*t**2 + 6*g*t + t**3 - 1)/3,
    3: -2*(g**3 + 3*g**2*t - 3*g**2 + 3*g*t**2 - 6*g*t + 3*g
           + 9*t**3 - 12*t**2 + 6*t - 2)/3,
    4: -2*(8*t**3 - 9*t**2 + 3*t - 1)/3,
    5: 2*(-g**3 + 6*g**2*t + 6*g*t**2 - 6*g*t + 1)/3,
    6: -2*(2*g**3 - 3*g**2 + 6*g*t**2 - 6*g*t + 3*g
           + 8*t**3 - 12*t**2 + 6*t - 2)/3,
    7: -2*(g**3 + 3*g**2*t - 3*g**2 + 3*g*t**2 - 6*g*t + 3*g
           + 9*t**3 - 12*t**2 + 6*t - 2)/3,
}
P_expr = {1: 1-g, 2: 1-t, 3: 1-t, 5: 1-g, 6: 2*t, 7: 2*t}
vertices = {
    1: [(F(0),F(0)), (F(1,3),F(0)), (F(1,3),F(1,3))],
    2: [(F(0),F(0)), (F(0),F(1)), (F(1,3),F(1,3))],
    3: [(F(0),F(1)), (F(1,3),F(1,3)), (F(1,3),F(2,3))],
    5: [(F(1,3),F(0)), (F(1,3),F(1,3)), (F(1,2),F(0))],
    6: [(F(1,3),F(1,3)), (F(1,2),F(0)), (F(1,2),F(1,2))],
    7: [(F(1,3),F(1,3)), (F(1,3),F(2,3)), (F(1,2),F(1,2))],
}
barycentric_samples = [(F(1,3),F(1,3),F(1,3)),
                       (F(1,6),F(1,3),F(1,2)),
                       (F(1,2),F(1,3),F(1,6))]


def sr(q):
    return sp.Rational(q.numerator, q.denominator)


symbolic_knots = [
    sp.Integer(0), t, 2*t,
    t+g, 2*t+g, 1+g, 1+t+g,
    2*t+2*g, 1+t+2*g, 2+2*g,
]


def symbolic_triangle(length, x, sample):
    # Select the correct linear branch of T_length at one cell interior.
    length_value = sp.Rational(sp.sympify(length).subs(sample))
    x_value = sp.Rational(sp.sympify(x).subs(sample))
    if x_value < 0 or x_value > 2*length_value:
        return sp.Integer(0)
    return x if x_value <= length_value else 2*length-x


def symbolic_rectangle_overlap(x, sample):
    # Select the correct branch of R_{t,1-t}; here t <= 1-t.
    x_value = sp.Rational(x.subs(sample))
    t_value = sp.Rational(t.subs(sample))
    b_value = 1-t_value
    if x_value < 0 or x_value > 1:
        return sp.Integer(0)
    if x_value <= t_value:
        return x
    if x_value <= b_value:
        return t
    return 1-x


def derive_region_formula(sample):
    ordered = sorted(
        symbolic_knots,
        key=lambda x: sp.Rational(sp.sympify(x).subs(sample)),
    )
    values = []
    for x in ordered:
        values.append(sp.expand(
            symbolic_triangle(t, x, sample)
            + 2*symbolic_rectangle_overlap(x-t-g, sample)
            + symbolic_triangle(1-t, x-2*t-2*g, sample)
        ))
    derived_n = sp.expand(sum(
        (x1-x0)*(y0*y0+y0*y1+y1*y1)/3
        for x0, x1, y0, y1
        in zip(ordered[:-1], ordered[1:], values[:-1], values[1:])
    ))
    peak_index = max(
        range(len(values)),
        key=lambda i: sp.Rational(values[i].subs(sample)),
    )
    return derived_n, sp.expand(values[peak_index])


region_samples = {
    1: {t: sp.Rational(1,4), g: sp.Rational(1,8)},
    2: {t: sp.Rational(1,5), g: sp.Rational(1,3)},
    3: {t: sp.Rational(1,5), g: sp.Rational(7,10)},
    5: {t: sp.Rational(2,5), g: sp.Rational(1,10)},
    6: {t: sp.Rational(2,5), g: sp.Rational(3,10)},
    7: {t: sp.Rational(2,5), g: sp.Rational(9,20)},
}
for region, sample in region_samples.items():
    derived_n, derived_p = derive_region_formula(sample)
    assert sp.expand(derived_n-N_expr[region]) == 0
    assert sp.expand(derived_p-P_expr[region]) == 0

# Region 4 has one formula for N and two peak branches.
for sample, expected_p in [
    ({t: sp.Rational(1,4), g: sp.Integer(2)}, 1-t),
    ({t: sp.Rational(2,5), g: sp.Integer(2)}, 2*t),
]:
    derived_n, derived_p = derive_region_formula(sample)
    assert sp.expand(derived_n-N_expr[4]) == 0
    assert sp.expand(derived_p-expected_p) == 0


checks = 0
for region, verts in vertices.items():
    for a, b, c in barycentric_samples:
        tt = a*verts[0][0] + b*verts[1][0] + c*verts[2][0]
        gg = a*verts[0][1] + b*verts[1][1] + c*verts[2][1]
        exact_n, exact_p, _, _ = exact_N_P(tt, gg)
        assert sp.expand(N_expr[region].subs({t: sr(tt), g: sr(gg)})) == sr(exact_n)
        assert sp.expand(P_expr[region].subs({t: sr(tt), g: sr(gg)})) == sr(exact_p)
        checks += 1

# Three checks in the separated region bring the total to 21.
for tt, gg in [(F(1,4),F(2)), (F(2,5),F(2)), (F(1,10),F(1))]:
    exact_n, exact_p, _, _ = exact_N_P(tt, gg)
    expected_p = max(1-tt, 2*tt)
    assert sp.expand(N_expr[4].subs(t, sr(tt))) == sr(exact_n)
    assert exact_p == expected_p
    checks += 1

assert checks == 21

opt_N, opt_P, opt_knots, opt_values = exact_N_P(F(1,3), F(1,3))
assert opt_N == F(40,81)
assert opt_P == F(2,3)
assert opt_N / opt_P == F(20,27)
print("all seven region formulas derived symbolically from the overlap profile")
print(f"additional exact rational-profile checks passed: {checks}")
print(f"equality point: N={opt_N}, P={opt_P}, Q={opt_N/opt_P}")
print("equality profile:", list(zip(opt_knots, opt_values)))

Saved output 1

all seven region formulas derived symbolically from the overlap profile
additional exact rational-profile checks passed: 21
equality point: N=40/81, P=2/3, Q=20/27
equality profile: [(Fraction(0, 1), Fraction(0, 1)), (Fraction(1, 3), Fraction(1, 3)), (Fraction(2, 3), Fraction(0, 1)), (Fraction(1, 1), Fraction(2, 3)), (Fraction(4, 3), Fraction(2, 3)), (Fraction(5, 3), Fraction(1, 3)), (Fraction(2, 1), Fraction(2, 3)), (Fraction(8, 3), Fraction(0, 1))]

The next cell verifies the nonnegative coefficient certificates for the six triangular regions and the factorizations for the separated region.

Show code

Code cell 10 · In [3]

u, v, w = sp.symbols("u v w", nonnegative=True)
Suvw = u + v + w

V = {
    region: [(sr(x), sr(y)) for x, y in verts]
    for region, verts in vertices.items()
}

certificates = {}
for region in (1, 2, 3, 5, 6, 7):
    tt = u*V[region][0][0] + v*V[region][1][0] + w*V[region][2][0]
    gg = u*V[region][0][1] + v*V[region][1][1] + w*V[region][2][1]
    d = sp.Poly(sp.expand((20*P_expr[region] - 27*N_expr[region]).subs({t: tt, g: gg})),
                u, v, w)

    # Homogenize to degree 3 on u+v+w=1. One further degree elevation
    # makes Regions 1 and 2 coefficient-positive.
    homogeneous = 0
    for powers, coefficient in d.terms():
        degree = sum(powers)
        homogeneous += (coefficient * u**powers[0] * v**powers[1] * w**powers[2]
                        * Suvw**(3-degree))
    target_degree = 4 if region in (1, 2) else 3
    if target_degree == 4:
        homogeneous *= Suvw

    poly = sp.Poly(sp.expand(homogeneous), u, v, w)
    coeffs = {
        (i, j, target_degree-i-j):
            poly.coeff_monomial(u**i * v**j * w**(target_degree-i-j))
        for i in range(target_degree + 1)
        for j in range(target_degree - i + 1)
    }
    equality_vertex = {
        1: (0,0,4), 2: (0,0,4),
        3: (0,3,0), 5: (0,3,0),
        6: (3,0,0), 7: (3,0,0),
    }[region]
    assert all(c >= 0 for c in coeffs.values())
    assert [powers for powers, c in coeffs.items() if c == 0] == [equality_vertex]
    assert all(c > 0 for powers, c in coeffs.items() if powers != equality_vertex)
    certificates[region] = coeffs
    print(f"Region {region}: degree {target_degree}, "
          f"{len(coeffs)-1} positive coefficients, zero only at {equality_vertex}")

N4 = N_expr[4]
assert sp.expand(19*(1-t) - 27*N4 - (3*t-1)*(48*t**2-38*t-1)) == 0
assert sp.expand(19*(2*t) - 27*N4 - 2*(3*t-1)*(24*t**2-19*t+9)) == 0
print("Region 4 factor identities verified exactly.")

Saved output 1

Region 1: degree 4, 14 positive coefficients, zero only at (0, 0, 4)
Region 2: degree 4, 14 positive coefficients, zero only at (0, 0, 4)
Region 3: degree 3, 9 positive coefficients, zero only at (0, 3, 0)
Region 5: degree 3, 9 positive coefficients, zero only at (0, 3, 0)
Region 6: degree 3, 9 positive coefficients, zero only at (3, 0, 0)
Region 7: degree 3, 9 positive coefficients, zero only at (3, 0, 0)
Region 4 factor identities verified exactly.

Station Reference. The theorem is stated in Archive #97. Its final seven-region calculation comes from Research Center Evaluation #1583, following preliminary symbolic work in Evaluations #1551, #1562, and #1575. The two cells above recheck the printed formulas at 21 rational interior points and replace the archived stationary-point search by nonnegative-coefficient certificates.

Related Work. The exact classification of two-interval indicators and the unique $1:1:2$ extremizer appears novel. It gives a finite-component counterpart to the global construction problem studied by Matolcsi and Vinuesa (2010).

1. Problem Definitions

2. S1. Binary step functions preserve the unrestricted supremum.

3. Additional findings

3.1 Verified 100,000100{,}000100,000-step construction

3.2 The sharp two-interval ceiling

3.1 Verified $100{,}000$ -step construction