% =====================================================================
\documentclass[11pt]{article}
\usepackage[T1]{fontenc}
\usepackage[utf8]{inputenc}
\usepackage{lmodern}
\usepackage[margin=1in]{geometry}
\usepackage{amsmath,amssymb,amsthm,mathtools}
\usepackage{booktabs,longtable,array}
\usepackage{enumitem}
\usepackage{needspace}
\usepackage{microtype}
\usepackage{xcolor}
\usepackage{tikz}
\usetikzlibrary{arrows.meta,positioning}
\input{glyphtounicode}
\input{glyphtounicode-cmex}
\pdfgentounicode=1
\usepackage{listings}
\usepackage{hyperref}
\usepackage[nameinlink,noabbrev]{cleveref}
\crefname{section}{Section}{Sections}
\crefname{appendix}{Appendix}{Appendices}
\hypersetup{colorlinks=true,linkcolor=blue!45!black,citecolor=blue!45!black,urlcolor=blue!45!black}
\newcommand{\C}{\mathbb C}
\newcommand{\R}{\mathbb R}
\newtheorem{theorem}{Theorem}
\newtheorem{lemma}[theorem]{Lemma}
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{corollary}[theorem]{Corollary}
\theoremstyle{definition}
\newtheorem{definition}[theorem]{Definition}
\theoremstyle{remark}
\newtheorem{remark}[theorem]{Remark}
\makeatletter
\newcommand{\theoremspace}{\if@nobreak\else\Needspace{3\baselineskip}\fi}
\makeatother
\AddToHook{env/theorem/before}{\theoremspace}
\AddToHook{env/lemma/before}{\theoremspace}
\AddToHook{env/proposition/before}{\theoremspace}
\AddToHook{env/corollary/before}{\theoremspace}
\AddToHook{env/definition/before}{\theoremspace}
\AddToHook{env/remark/before}{\theoremspace}
\setlist{itemsep=3pt,topsep=5pt}
\lstset{basicstyle=\ttfamily\small,breaklines=true,columns=fullflexible,keepspaces=true,showstringspaces=false}
\emergencystretch=1em
\begin{document}
% =====================================================================
The result is one geometric statement. In complex six-dimensional space
$\C^6$ there are three mutually unbiased orthonormal bases, and no four.
In symbols, $N(6)=3$.
The proof is long for an honest reason. Some of modern geometry ends in
a finite check, and a finite check counts as a proof only when every
possible bug, rounding error, and missed case is excluded. So most of
the work is \emph{certification}: every test the computer runs carries a
strict safety margin, so that ordinary floating-point arithmetic cannot
turn a true solution into a spurious ``not found.''
The strategy fits in six steps.
\begin{enumerate}
\item \textbf{Remove the first basis.} A unitary change of coordinates
sends any orthonormal basis to the standard one. The second basis is
then the column set of $H/\sqrt6$, where every entry of $H$ has modulus
$1$ and $HH^*=6I$. Such $H$ is a \emph{complex Hadamard matrix}, and
the question becomes: which $H$ admit two further unbiased bases?
\item \textbf{Pin down the phase freedom.} Multiplying rows and columns
of $H$ by phases changes nothing. Phases are chosen so that the first
row and column of $H$ are $1$, after which $H$ is fixed by nine real
phase parameters, the rest of the first row and column, together with a
$4\times4$ block of unknown phases.
\item \textbf{Cover the seed.} Those nine parameters trace a subset of
the nine-torus, and finitely many balls of radius $0.0320001$, measured
in cycles, cover it. This is the seed cover.
\item \textbf{Complete and refine.} For a seed ball, the $4\times4$
block must satisfy Hadamard equations that compress to a single
$V\in U(2)$. Finitely many boxes cover $U(2)$, and each box refines
into small phase balls holding the actual solution set. This is the
hardest geometry in the argument.
\item \textbf{Exclude two further bases.} Take one phase ball for $H$.
Every vector unbiased to both the standard basis and $H/\sqrt6$ has the
form $e(\beta)/\sqrt6$, and a phase-bin search covers them all by small
vector balls. Two further bases would need twelve of those balls: six
pairwise orthogonal in one group, six in the other, with all $36$ cross
pairs unbiased. Those conditions are edges of two graphs, and a finite
search rules out the required pair of six-cliques.
\item \textbf{Account for everything.} The execution of these four
stages on a 2026 machine, with a fixed build and fixed numeric margins,
is recorded together with the certified margins. The finite chain
reaches every possible $H$ up to the allowed symmetries, and each local
test either excludes soundly or fails openly. Nothing is guessed.
\end{enumerate}
\textbf{Formality.} All sets of phase coordinates are closed, as are all
balls, ellipsoids, and boxes. Overlaps, and repeated overlap with phase
representatives, are harmless: duplication weakens an exclusion but
never omits a solution. ``Small ball'' always means the phase matrices
within a stated Frobenius norm. A numerical center is a proposal and
need not be an exact Hadamard matrix. Every conclusion comes from a
necessary-condition inequality verified on a whole ball, not from the
residual at a proposed center.
% =====================================================================
\medskip\noindent\textbf{Introduction.}
% =====================================================================
Two orthonormal bases $\mathcal B$ and $\mathcal C$ of $\C^d$ are
\emph{mutually unbiased} (MUB) if
\[
|\langle b,c\rangle|^2=\frac 1d \qquad\text{for every }b\in\mathcal B,\ c\in\mathcal C .
\]
Here $\langle u,v\rangle=\sum_j\overline{u_j}v_j$.
For $d\ge 2$ let $N(d)$ be the maximum size of a pairwise mutually
unbiased family of orthonormal bases in $\C^d$. The main result is:
\begin{theorem}[Computer-assisted]\label{thm:main}
Under the binary64 arithmetic and compiler conditions stated in
the certified arithmetic, the execution documented in
certifies that there are no four pairwise
mutually unbiased orthonormal bases in $\C^6$. Consequently
$N(6)=3$.
\end{theorem}
The upper bound is the hard direction. The lower bound (three bases
exist) is an explicit computation, .
\medskip\noindent\textbf{MUBs, Hadamard matrices, and the shape of the problem.}
Fixing one basis by a common unitary transformation costs nothing; every
family can be moved so that the first member is the standard basis. Then
the second member is the columns of $H/\sqrt6$ for a complex matrix
satisfying
\[
|H_{ij}|=1,\qquad HH^*=6I .
\]
Such an $H$ is called a \emph{complex Hadamard matrix}. The issue is to
decide whether a given pair (standard basis, $H$/$\sqrt6$) can be extended
by up to two more MUB bases, and to do that uniformly over the whole
(unclassified) space of complex Hadamard matrices. Because $H$ is
square, $HH^* = 6I$ implies $H^*H = 6I$ as well.
The answer we prove, Theorem~\ref{thm:main}, is that no such extension
exists, for any $H$. The proof covers the full space of $H$, not just a
known family: the second basis is completely unconstrained, and the
argument handles complex Hadamard matrices outside all previously
studied families.
\medskip\noindent\textbf{Historical context.}
Mutually unbiased bases describe projective measurements in which a
state drawn from one basis gives a uniformly random outcome in any
other. Schwinger placed this complementarity inside his calculus of
unitary operator bases in 1960 (Proc. Nat. Acad. Sci. USA 46, 570--579).
Ivanovi\'c raised the state-determination problem and gave the
prime-dimensional construction in 1981 (J. Phys. A 14, 3241--3245).
Wootters and Fields proved that a complete set of such bases is
statistically optimal for state determination and built one whenever
the dimension is a prime power (Ann. Physics 191, 363--381, 1989).
Every dimension obeys $N(d)\le d+1$, and the short operator-space proof
is due to Bandyopadhyay, Boykin, Roychowdhury and Vatan (Algorithmica
34, 512--528, 2002).
{\small\url{https://www.pnas.org/doi/10.1073/pnas.46.4.570}}
{\small\url{https://doi.org/10.1088/0305-4470/14/12/019}}
{\small\url{https://www.sciencedirect.com/science/article/pii/0003491689903229}}
{\small\url{https://arxiv.org/pdf/quant-ph/0103162}}
Dimension six is the first case that is neither a prime power nor a
small even dimension settled by other constructions. A tensor-product
construction already yields three bases in $\C^6$, tensoring a qubit
pair with a qutrit triple; that is Section~4 of the Klappenecker and
R\"otteler study of MUB constructions. Zauner conjectured in his 1991
Vienna Diplomarbeit that three is the maximum, the same statement
reappears in his doctoral thesis, and it is Conjecture~1 of the
Klappenecker--R\"otteler paper. A survey of the composite-dimension
setting is due to McNulty and Weigert.
{\small\url{https://arxiv.org/pdf/quant-ph/0309120}}
{\small\url{https://www.gerhardzauner.at/documents/gzdiplthd.pdf}}
{\small\url{https://arxiv.org/pdf/2410.23997}}
Four attacks ran alongside. Grassl proved that the standard basis
together with any member of the two-parameter Fourier family cannot
extend to four bases, Theorem~2 of his dimension-six note. Bengtsson,
Bruzda, Ericsson and \.Zyczkowski mapped the order-six Hadamard
families and the unbiased vectors they carry. Butterley and Hall, and
independently Brierley and Weigert, produced the numerical evidence
that no four bases exist. Brierley and Weigert also enumerated the
unbiased vectors of individual Hadamards through algebraic equations:
Section~3.4 of that work gives rigorous approximation bounds for exact
input matrices, and Section~4.3 explains why rounding coefficients in
sampled non-affine families blocks a complete cover by that route.
Raynal, Lu and Englert sharpened the numerical side by maximising the
average distance between four bases and finding the maximum strictly
below unity, with a two-parameter family of three bases attaining it.
{\small\url{https://arxiv.org/pdf/quant-ph/0406175}}
{\small\url{https://arxiv.org/pdf/math.OA/0609076}}
{\small\url{https://arxiv.org/pdf/0901.4051}}
{\small\url{https://arxiv.org/pdf/1103.1025}}
The closest antecedent is the certified discretization of Jaming,
Matolcsi, M\'ora, Sz\"oll\H osi and Weiner. They proved that the
standard basis together with any member of the two-parameter Fourier
family cannot extend to four bases, Theorem~1.4 of their 2009 paper,
and their Section~3 introduces closed phase intervals, rigorous error
bounds, subdivision, and exhaustive compatibility tests. Jaming,
Matolcsi and M\'ora then set out the unrestricted discretization
programme in Sections~2--4 of their survey. Those four devices carry
the finite search throughout: what is added is exhaustive coverage of
arbitrary second bases, with every error bounded.
{\small\url{https://arxiv.org/pdf/0902.0882}}
{\small\url{https://arxiv.org/pdf/1201.0640}}
Two recent lines of work are complementary rather than overlapping.
Matolcsi, Matszangosz, Varga and Weiner study triplets of MUBs under
conjectural sufficient conditions, Section~4 of their preprint.
C\'ardenes Wuttig and Tindall classify order-six Hadamard matrices
completely and leave MUB extendibility open. Only the combination of a
finite cover of all second bases with the local exclusion, executed
under a certified arithmetic contract, yields an unconditional upper
bound.
{\small\url{https://arxiv.org/pdf/2503.14752}}
{\small\url{https://arxiv.org/pdf/2608.18053}}
\medskip\noindent\textbf{The strategy at a glance.}
The proof has five layers, each of which strictly refines the previous
one and preserves every exact candidate solution.
\begin{enumerate}
\item \textbf{Seed cover.}
Row and column rephasings and transpositions place any $H$ into a
standard form with $9$ free phases. A symmetry reduction plus
necessary interval tests cover the possible values of these nine
phases by finitely many seed balls (radius $0.0320001$ in cycle units).
\item \textbf{Completion.}
In a block decomposition of $H$ with row and column sizes $2+4$, the
seed pins down the upper-left, upper-right and lower-left blocks. The
lower-right $4\times 4$ block is the only remaining unknown. At any
seed center that passes the frame checks, the upper-right block and
the conjugate transpose of the lower-left block have rank two, and
so they admit two-dimensional kernels in $\C^4$. Compressing the
unknown lower-right block between these kernels produces a $2\times 2$
matrix whose polar part is a unitary $V\in U(2)$. A closed cover of
$U(2)$ therefore accounts for every feasible completion, and uniform
Taylor bounds within each cell either exclude the cell or place all
its exact completions in a certified phase ball
( and~\ref{sec:phase-cover}).
\item \textbf{Refinement.}
Each phase ball carries numerical tangent directions and a rigorous
residual-identity estimate; covering the tangent coordinates while
keeping the residual bound yields smaller balls that still cover
every exact solution. No manifold parametrization is used, and the
estimates remain valid even when the proposed center is not an exact
Hadamard matrix.
\item \textbf{Exclusion.}
For a final phase ball, a small phase-bin search covers every vector
unbiased to both the standard basis and $H/\sqrt6$ by closed vector
balls. These vector balls become vertices of two graphs, one for
orthogonality and one for unbiasedness. Necessary numerical tests
remove only impossible edges. Two additional MUB bases would require
two disjoint six-cliques with all $36$ cross edges present; that
pattern is ruled out by finite search.
\item \textbf{Execution certificate.}
A complete run of the four programs, on a specific build and platform,
accounts for every seed, every refined ball, and every initially
unresolved record. Nothing is left open: a record that is not
excluded must be covered upstream, and every fallback record must
produce a negative exclusion.
\end{enumerate}
This restates the same six-bullet story from the preface, now with the
section labels attached. If you only remember one thing, remember this:
the argument is a chain of strictly finer covers; each exclusion uses
only a \emph{necessary} condition with a \emph{certified} margin; and no
step in the chain ever relies on a floating-point center being exact.
\medskip\noindent\textbf{The elementary lower bound.}
The lower bound uses only the tensor-product construction. We give it
explicitly so the reader can see that attainment of three bases is
separate from the certified exclusion of four.
\begin{lemma}[Three mutually unbiased bases in dimension six]
\label{lem:three-mub-six}
There exist three pairwise mutually unbiased orthonormal bases in
$\C^6$.
\end{lemma}
\begin{proof}
Identify $\C^6$ with $\C^2\otimes\C^3$ in the standard product
coordinates. Define
\[
A_0=I_2,\qquad
A_1=\frac 1{\sqrt2}\begin{pmatrix}1&1\\1&-1\end{pmatrix},\qquad
A_2=\frac 1{\sqrt2}\begin{pmatrix}1&1\\ i&-i\end{pmatrix}.
\]
These are unitary, and for $r\ne s$, every entry of $A_r^*A_s$ has
squared modulus $1/2$; the only cross product not immediate from the
entries is
\[
A_1^*A_2=\frac 12\begin{pmatrix}1+i&1-i\\1-i&1+i\end{pmatrix}.
\]
Let $\omega=e^{2\pi i/3}$ and set
\[
F=\frac 1{\sqrt3}\bigl(\omega^{jk}\bigr)_{j,k=0}^2,\qquad
D=\operatorname{diag}(1,\omega,\omega),\qquad
B_0=I_3,\ B_1=F,\ B_2=DF.
\]
The relation $1+\omega+\omega^2=0$ makes $F$ unitary; $D$ is unitary as
well. Both $F$ and $DF$ have entries of squared modulus $1/3$. For their
cross product,
\[
(F^*DF)_{jk}=\frac 13\sum_{x=0}^2\omega^{x^2+(k-j)x}.
\]
If $t=k-j$ is read modulo $3$, the sum is $1+2\omega$ for $t=0$ and
$2+\omega^2$ for $t=1,2$. Each has squared modulus $3$. Consequently
every entry of $B_r^*B_s$ has squared modulus $1/3$ when $r\ne s$.
For $r=0,1,2$, let the columns of $U_r=A_r\otimes B_r$ form the $r$th
basis. Each $U_r$ is unitary. For $r\ne s$ we have
$U_r^*U_s=(A_r^*A_s)\otimes(B_r^*B_s)$, so every entry has squared
modulus $(1/2)(1/3)=1/6$. Thus the three bases are mutually unbiased,
and $U_0=I_6$.
\end{proof}
\medskip\noindent\textbf{A Lie-algebra consequence.}
An orthonormal basis $\mathcal B$ of $\C^6$ determines its Cartan
subalgebra $\mathfrak h_{\mathcal B}$ of traceless matrices that are
diagonal in that basis. It is spanned by the matrices $bb^*-I/6$ for
$b\in\mathcal B$. The Killing form on $\mathfrak{sl}_6(\C)$ is
$K(X,Y)=12\operatorname{tr}(XY)$; two subalgebras are
\emph{Killing-orthogonal} when this form vanishes on every pair of
elements. For unit vectors $b,c$,
\[
\operatorname{tr}\bigl((bb^*-I/6)(cc^*-I/6)\bigr)
=|\langle b,c\rangle|^2-\frac 16 .
\]
So two bases are mutually unbiased exactly when their associated Cartan
subalgebras are Killing-orthogonal. These algebras are closed under the
standard adjoint, and the converse correspondence also uses this
adjoint condition.
\begin{corollary}[Adjoint-closed Cartan subalgebras]
\label{cor:adjoint-cartan}
Under the arithmetic and compiler conditions of
the certified arithmetic, using the complete execution recorded above
and certified in Theorem~\ref{thm:main}, the maximum
number of pairwise Killing-orthogonal Cartan subalgebras of
$\mathfrak{sl}_6(\C)$ that are all closed under the standard adjoint
$X\mapsto X^*=\overline X^{\,T}$ is three.
\end{corollary}
\begin{proof}
By Theorem 5.2 of Boykin, Sitharam, Tiep and Wocjan, such collections of
Cartan subalgebras correspond one-to-one to collections of mutually
unbiased orthonormal bases of the same cardinality. Theorem~\ref{thm:main}
gives the upper bound under its stated contract, and
Lemma~\ref{lem:three-mub-six} gives attainment.
\end{proof}
{\small\url{https://arxiv.org/pdf/quant-ph/0506089}}
% =====================================================================
\medskip\noindent\textbf{Equivalences and an exhaustive seed cover.}
\label{sec:seed}
% =====================================================================
The first layer of the cover replaces every possible second basis by one
of finitely many nine-phase seed balls. The strategy, closed phase
discretization plus necessary compatibility tests, follows Section 3 of
Jaming, Matolcsi, M\'ora, Sz\"oll\H osi and Weiner. The unrestricted
programme is set out in Sections 2--3 of the later survey by Jaming,
Matolcsi and M\'ora. What is new is the specific symmetry reduction and
shared-coordinate estimate that prove completeness of this exact cover.
{\small\url{https://arxiv.org/pdf/0902.0882}}
{\small\url{https://arxiv.org/pdf/1201.0640}}
Throughout, phases are in cycles and $e(t)=\exp(2\pi i t)$
entrywise. Matrix indices run $0,\ldots,5$, entries are stored in
row-major order with $(i,j)$ in position $6i+j$, and the Hermitian inner
product is $\langle u,v\rangle=\sum_j\overline{u_j}v_j$. Vectors and
phase matrices use Euclidean and Frobenius norms respectively; operator
norms are induced Euclidean unless noted. An angle is in radians only
when explicitly stated.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Equivalences and phase gauges.}
% ---------------------------------------------------------------
\begin{definition}\label{seed:extendibility}
A matrix $H\in\C^{6\times 6}$ is \emph{Hadamard} if every entry has
modulus one and $HH^*=6I$. It is \emph{extendible} if the standard
basis and the columns of $H/\sqrt6$ belong to a family of four
mutually unbiased orthonormal bases.
\end{definition}
Unbiasedness to the standard basis forces every entry of the second
basis matrix to have modulus $1/\sqrt6$; the matrix is then $H/\sqrt6$
with $H$ Hadamard. Because $H$ is square, $HH^*=6I$ implies $H^*H=6I$.
Every such $H$ has a real phase lift $P$ with $H=e(P)$.
\begin{lemma}[Equivalences]\label{seed:equivalences}
Extendibility is invariant under row or column permutations,
multiplication of rows or columns by unit scalars, complex conjugation,
and transposition.
\end{lemma}
\begin{proof}
Permuting or rephasing columns only reorders or rephases vectors within
the second basis. A row permutation or row rephasing is a common unitary
operation on all four bases; it sends the standard basis to a new basis
but preserves both orthonormality and the moduli of inner products
among the four bases.
For transposition, write $U=H/\sqrt6$ and apply the common unitary
$U^*$ to a family whose first two basis matrices are $I,U$. Those two
become $U^*,I$; interchange their roles and conjugate every coordinate.
The first two matrices are now $I,U^T$. Every operation is reversible,
so the statement holds in both directions.
\end{proof}
For a real $6\times 6$ matrix $z$ define double centering by
\begin{equation}\label{seed:projection-formula}
(\Pi z)_{ij}
=z_{ij}-\frac 16\sum_b z_{ib}-\frac 16\sum_a z_{aj}
+\frac 1{36}\sum_{a,b}z_{ab}.
\end{equation}
Its range is the \emph{horizontal} subspace of matrices whose every
row and column sum is zero. We also use unwrapped dephasing
\begin{equation}\label{seed:dephasing-formula}
D_0(P)_{ij}=P_{ij}-P_{i0}-P_{0j}+P_{00}.
\end{equation}
\begin{lemma}[Horizontal phase representatives]\label{seed:horizontal}
$\Pi$ is the orthogonal projection onto the horizontal subspace; in
particular $\|\Pi z\|_F\le \|z\|_F$. If two phase representatives differ
by $z$, the second is equivalent under row and column rephasings to the
first plus $\Pi z$. The statement survives arbitrary entrywise integer
changes of phase representatives.
\end{lemma}
\begin{proof}
Equation~\eqref{seed:projection-formula} has zero row and column sums,
and $z-\Pi z$ is a sum of a matrix constant on each row and a matrix
constant on each column. Every such sum is orthogonal to the horizontal
subspace in the real Frobenius inner product, so $\Pi$ is orthogonal
projection onto that subspace. If $Q=P+z$, then
$P+\Pi z=Q-(z-\Pi z)$ differs from $Q$ by row and column phases,
allowed by Lemma~\ref{seed:equivalences}. Integer additions leave
$e(P)$ unchanged.
\end{proof}
Similarly $D_0(P)$ gives an equivalent matrix whose first row and column
are zero. If a wrapped difference is used, its chosen integers determine
a legitimate phase representative; no rounded decision at a half-integer
is claimed to be uniquely preferred.
% ---------------------------------------------------------------
\medskip\noindent\textbf{A symmetry arrangement for nine phases.}
% ---------------------------------------------------------------
The \emph{gap} of a finite circular multiset is its largest circular
gap, equivalently the length of an open arc that contains no point of
the multiset. Repeated points count with multiplicity; a repeated point
contributes a zero gap. For every pair of distinct rows of $P$ consider
their six phase differences modulo one, and similarly for columns.
\begin{lemma}[Seed arrangement]\label{seed:arrangement}
Every exact Hadamard matrix has an equivalent dephased representative
for which a pair of rows attaining the minimum gap among all row and
column pairs is the first pair, and its phases can be written
\begin{equation}\label{seed:x-arrangement}
x=(0,x_1,x_2,x_3,x_4,x_5),\quad
0\le x_1\le\cdots\le x_5\le 1,\quad
g_x=1-x_5,\quad x_1\le x_5-x_4 .
\end{equation}
After any one of these five noninitial columns is selected and moved to
position two, the last four rows can be permuted so that the second
column has phases
\begin{equation}\label{seed:y-arrangement}
y=(0,\xi,y_2,y_3,y_4,y_5),\qquad 0\le y_2\le y_3\le y_4\le y_5\le 1 ,
\end{equation}
where $\xi$ is the selected entry of $x$. These lists satisfy
\begin{equation}\label{seed:local-conditions}
\sum_{j=0}^5 e(x_j)=0,\qquad
\sum_{j=0}^5 e(y_j)=0,\qquad g_y\ge g_x .
\end{equation}
\end{lemma}
\begin{proof}
There are finitely many row and column pairs, so a pair of minimum gap
exists. Transpose if it is a column pair, then permute rows so it is the
first pair. Choose one of its largest gaps, start the circular ordering
immediately after that gap, and use that column as the first column.
Dephase the first row and column. Ordering the remaining columns along
the complementary arc gives $g_x=1-x_5$ and the sorted list
in~\eqref{seed:x-arrangement}.
If the final orientation inequality fails, interchange the two selected
rows, use the old $x_5$ column as the starting column, and reverse the
ordering. The new list is
\[
x'_0=0,\qquad x'_i=x_5-x_{5-i}\quad(1\le i\le 5).
\]
It has the same distinguished largest gap, and
$x'_1=x_5-x_4$ whereas $x'_5-x'_4=x_1$. So this orientation satisfies
the required inequality whenever the other fails. The construction
uses a chosen largest gap and is unchanged by ties or repeated points.
Select a noninitial column and place it at position two. Permuting the
last four rows sorts its last four entries without changing either of
the first two rows. Orthogonality of those rows gives the first
vanishing sum in~\eqref{seed:local-conditions}; orthogonality of the
first two columns gives the second. Row and column rephasings rotate the
relevant difference multisets; permutations merely reorder points or
relabel pairs. Transposition interchanges the two collections of pairs,
so the minimum gap stays $g_x$ and the selected column pair has
$g_y\ge g_x$.
\end{proof}
Set $N=112$ and $w=1/N$. Use closed bins
\begin{equation}\label{seed:bins}
I_b=[b/N,(b+1)/N],\qquad u_b=(b+\tfrac12)/N,\qquad 0\le b<N .
\end{equation}
They cover $[0,1]$ including both endpoints. Sorted exact phases admit
nondecreasing bin indices, so write the five indices for
$x_1,\ldots,x_5$ as
$b^X=(a,b,c,d,l)$ with $0\le a\le b\le c\le d\le l<N$. The orientation
condition implies the necessary integral restriction
\begin{equation}\label{seed:orientation-filter}
\frac aN\le x_1\le x_5-x_4\le\frac{l+1-d}{N},
\qquad a\le l-d+1 .
\end{equation}
For this index tuple choose $s\in\{0,\ldots,4\}$ as the first minimizer
of $|2b^X_s+1-N|$, the rule used in \texttt{seeds::mk}. The selected
original column is $s+1$, the selected exact phase is $\xi=x_{s+1}$, and
the new column order is: original column zero, column $s+1$, and the
other four columns in their former order. This specifies only a column;
it imposes no extra condition on $\xi$.
After arranging $y$ as in~\eqref{seed:y-arrangement} its five varying
entries have bins
$b^Y=(b^X_s,b^Y_1,b^Y_2,b^Y_3,b^Y_4)$ with
$0\le b^Y_1\le b^Y_2\le b^Y_3\le b^Y_4<N$. The shared entry uses the
same bin as in $x$; it is not required to be ordered relative to the
four tail bins.
\begin{lemma}[Stability of circular gaps]\label{seed:gap-stability}
If each point of one circular multiset is matched to a point of another
at circular distance $\le\delta$, their gaps differ by $\le 2\delta$.
\end{lemma}
\begin{proof}
Take an empty open arc of length $L$ in the first multiset. If $L>2\delta$,
trim $\delta$ from each end; a point of the second multiset inside the
remaining open arc would have its matched first point in the original
empty arc, contradiction. Hence the second gap is at least $L-2\delta$.
If $L\le 2\delta$ the same bound follows from nonnegativity. Apply to a
largest gap and interchange the two multisets. No preservation of sorted
order is used.
\end{proof}
\begin{lemma}[Necessary midpoint gap filters]\label{seed:gap-filters}
For every bin pair containing lists as in
Lemma~\ref{seed:arrangement}, the midpoint gaps (with the fixed zero)
satisfy
\begin{equation}\label{seed:gap-inequalities}
g_{x,\mathrm{mid}}\le 1-\frac{l-1}{N},\qquad
g_{y,\mathrm{mid}}\ge 1-\frac{l+2}{N}.
\end{equation}
\end{lemma}
\begin{proof}
Each varying phase is within $w/2$ of its midpoint and the initial zero
is unchanged. Lemma~\ref{seed:gap-stability} bounds each gap change by
$w$. Since $x_5\in[l/N,(l+1)/N]$ we have
$1-(l+1)/N\le g_x\le 1-l/N$. Combine with $g_y\ge g_x$.
\end{proof}
For sorted midpoint bin indices $q_0,\ldots,q_4$ the gap used by
\texttt{mk} is
\[
\max\Bigl\{\frac{q_0+\tfrac12}{N},
\ \frac{q_1-q_0}{N},\ldots,\frac{q_4-q_3}{N},\
1-\frac{q_4+\tfrac12}{N}\Bigr\}.
\]
For a $Y$ tuple it sorts a copy solely for this computation, preserving
the stored shared-first ordering.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Vanishing-sum filters and the shared displacement.}
% ---------------------------------------------------------------
For either of the two lists $x$ or $y$, let $u=(u_0,\ldots,u_4)$ be its
five varying midpoints. Define real two-vector $c$ and real $2\times 5$
matrix $J$ by
\begin{equation}\label{seed:residual-data}
c=\frac{1+\sum_{j=0}^4 e(u_j)}{2\pi i},\qquad J_j=e(u_j),
\end{equation}
realified in both definitions.
\begin{lemma}[Necessary sum test]\label{seed:sum-test}
If $1+\sum_{j=0}^4 e(u_j+\Delta_j)=0$ and $|\Delta_j|\le w/2$ then
\begin{equation}\label{seed:taylor-residual}
\|c+J\Delta\|\le\frac{5\pi w^2}{4}.
\end{equation}
In particular for every $v\in\R^2$,
\begin{equation}\label{seed:sum-support}
|v\cdot c|\le\frac{5\pi w^2}{4}\|v\|
+\frac w2\sum_{j=0}^4|v\cdot J_j|.
\end{equation}
\end{lemma}
\begin{proof}
The second derivative of $e(t)$ has modulus $4\pi^2$, so Taylor's formula
with integral remainder gives
\[
\Bigl|\frac{e(u_j+\Delta_j)-e(u_j)}{2\pi i}-e(u_j)\Delta_j\Bigr|
\le\pi\Delta_j^2 .
\]
Sum the five remainders and use the vanishing-sum hypothesis to
obtain~\eqref{seed:taylor-residual}. Dot with $v$, apply the triangle
inequality and $|\Delta_j|\le w/2$ to get~\eqref{seed:sum-support}.
\end{proof}
The function \texttt{sumok} tests~\eqref{seed:sum-support} in the
residual direction and in directions perpendicular to the five
unit-exponential columns. Since the inequality holds for every real
direction, their numerical construction need not reproduce any
prescribed exact direction. One-sided margins for evaluating the
residuals and supports are established in the certified arithmetic.
In grid units $h=N\Delta$, write
\begin{equation}\label{seed:grid-error}
E=\frac{5\pi}{4N},\qquad \|Jh+Nc\|\le E,\qquad |h_j|\le\tfrac12 .
\end{equation}
\begin{lemma}[Shared-coordinate intervals]\label{seed:shared-interval}
Let $t$ be the position of a shared entry in either list. Every feasible
displacement satisfies, for every $v\in\R^2$,
\begin{equation}\label{seed:clip-inequality}
|Nv\cdot c+(v\cdot J_t)h_t|
\le E\|v\|+\frac 12\sum_{j\ne t}|v\cdot J_j|.
\end{equation}
If $a=v\cdot J_t\ne 0$ and the right side is $r$, this is the interval
constraint
\begin{equation}\label{seed:clip-interval}
-\frac{Nv\cdot c}{a}-\frac r{|a|}
\le h_t\le
-\frac{Nv\cdot c}{a}+\frac r{|a|}.
\end{equation}
The intervals obtained for the shared entry from $X$ and $Y$ must
intersect.
\end{lemma}
\begin{proof}
Separate the shared term in~\eqref{seed:grid-error}, take its dot
product with $v$, and bound each remaining coordinate by $1/2$. Solving
the resulting absolute-value inequality gives~\eqref{seed:clip-interval}
for either sign of $a$. The shared midpoint is identical in both lists,
so the same exact grid displacement satisfies both sets of necessary
intervals.
\end{proof}
The function \texttt{extend} begins with $[-1/2,1/2]$ and intersects
these constraints using the shared column direction and directions
perpendicular to the other four columns; here $t=s$ for $X$ and $t=0$
for $Y$. It skips divisors of magnitude below $10^{-6}$, which only
removes a constraint. Numerator and endpoint margins enlarge the
computed intervals (the certified arithmetic). The pair is rejected
only when the two intervals are strictly disjoint; a common endpoint is
retained.
The enumeration in \texttt{seeds::build} now has an unambiguous
meaning. It visits every weakly increasing five-bin $X$ tuple and
applies~\eqref{seed:orientation-filter}, the first inequality
in~\eqref{seed:gap-inequalities}, and \texttt{sumok}. Separately, for
every possible first bin, it visits all weakly increasing four-bin $Y$
tails and applies \texttt{sumok}. The program \texttt{covout} matches
the shared bin, applies the second gap inequality, and intersects the
shared displacement intervals. So every discarded exact solution would
violate one of these necessary conditions. Numerical gap comparisons
include outward slack (the certified arithmetic).
% ---------------------------------------------------------------
\medskip\noindent\textbf{Certified reuse of seed centers.}
% ---------------------------------------------------------------
The nine seed coordinates, before moving the selected column, are
\begin{equation}\label{seed:nine-coordinates}
q=(x_1,x_2,x_3,x_4,x_5,y_2,y_3,y_4,y_5)\in\R^9 .
\end{equation}
For a bin pair let $u\in\R^9$ be its rational midpoint and put
$h=N(q-u)$. Thus $h\in[-1/2,1/2]^9$; the shared displacement is $h_s$.
Retaining $u$ as a center would cover the whole bin with radius
$3/(2N)<0.032$. To reduce the number of retained centers, we instead
test whether the feasible portion lies in a radius-$0.032$ ball about a
previously retained midpoint. If no such test succeeds we retain $u$
itself.
Let $J_X,J_Y,c_X,c_Y$ be the data
in~\eqref{seed:residual-data}. Define the real $4\times 9$ matrix $J_*$
by its action
\begin{equation}\label{seed:stacked-system}
J_*h=\begin{pmatrix}
J_X(h_0,h_1,h_2,h_3,h_4)^T\\
J_Y(h_s,h_5,h_6,h_7,h_8)^T
\end{pmatrix},\qquad
p_*=\begin{pmatrix}-Nc_X\\ -Nc_Y\end{pmatrix}.
\end{equation}
For feasible $h$ we have $J_*h=p_*+r_*$, where the two real two-vector
blocks $r_X,r_Y$ of $r_*$ satisfy $\|r_X\|,\|r_Y\|\le E$.
\begin{lemma}[Support bounds for center reuse]\label{seed:dual-support}
Let $u'$ be a previously retained midpoint with the same selected
position $s$, and let $k=N(u-u')\in\mathbb Z^9$. Every feasible
displacement in the current bin pair satisfies
\begin{equation}\label{seed:distance-bound}
N^2\|q-u'\|^2=\|k+h\|^2
\le\|k\|^2+\frac 94+2k\cdot h .
\end{equation}
Two upper bounds for the last term are
\begin{align}
2k\cdot h&\le\sum_{j=0}^8|k_j|,
\label{seed:box-support}\\
2k\cdot h&\le
\sum_{j=0}^8|k_j-(J_*^Tv)_j|
+2v\cdot p_*+2E(\|v_X\|+\|v_Y\|),
\qquad v\in\R^4,
\label{seed:dual-formula}
\end{align}
where $v_X,v_Y$ are the two real two-vector blocks of $v$.
\end{lemma}
\begin{proof}
Expand the square in~\eqref{seed:distance-bound} and use
$\|h\|^2\le 9/4$. The coordinate bounds give~\eqref{seed:box-support}.
For the other bound write
\[
2k\cdot h
=2(k-J_*^Tv)\cdot h+2v\cdot p_*+2v\cdot r_* .
\]
The first term is at most $\sum_j|k_j-(J_*^Tv)_j|$ by the same
coordinate bound. Cauchy--Schwarz on the two residual blocks bounds
the last term by $2E(\|v_X\|+\|v_Y\|)$.
\end{proof}
The function \texttt{supdual} proposes vectors $v$ and evaluates the
right side of~\eqref{seed:dual-formula}. Validity of a proposal is
independent of how it was generated and of the accuracy of the
heuristic linear solve; only finite components of magnitude at most
$10^4$ are used. Its evaluation replaces $1.25\pi/N$ by the larger
$1.251\pi/N$ and adds a final support margin; the certified arithmetic
proves the resulting bound is conservative.
Set $R=0.032$. For a proposed previous center, \texttt{testBucket} first
tries~\eqref{seed:box-support} and then, if applicable, the dual bound.
It returns success only when the upper bound in~\eqref{seed:distance-bound}
is below the safely lowered comparison threshold for $(NR)^2$.
Success therefore certifies the entire feasible portion of the current
bin pair, not merely its midpoint. Distance cutoffs, bucket indices, and
reach tests only decide which previous centers to try; omitting a
possible previous center cannot produce a false success.
For clarity: the shared column is counted once
in~\eqref{seed:nine-coordinates}, though its displacement enters both
equations in~\eqref{seed:stacked-system}. Reuse requires equality of the
selected position $s$. Moving that position to the start of the five
$X$ coordinates applies the same coordinate permutation to the current
point and to the retained center, and preserves their distance. The
resulting order is: the selected $X$ entry, the other four $X$ entries
in order, and the four $Y$ tail entries. These are exactly the nine
entries of the dephased first two rows and columns used by the next
stage. Centers retained earlier for the same $X$ tuple are handled
identically, with the first five components of $k$ zero.
% ---------------------------------------------------------------
\medskip\noindent\textbf{The exhaustive seed cover proposition.}
% ---------------------------------------------------------------
\begin{proposition}[Exhaustive seed cover]\label{seed:cover}
For a completed execution of the specified \texttt{covout} traversal
under the arithmetic conditions of the certified arithmetic, the
written rational midpoint centers, equipped with their recorded selected
positions, cover the seed coordinates of every exact complex Hadamard
matrix up to the equivalences of Lemma~\ref{seed:equivalences}, by
closed Euclidean balls of radius $0.032$. After the indicated column
ordering and the conversion of centers to binary phase values, radius
$0.0320001$ suffices.
\end{proposition}
\begin{proof}
Given an exact Hadamard matrix, first apply
Lemma~\ref{seed:arrangement} to arrange $x$. Choose a nondecreasing
closed bin assignment, determine $s$ from that tuple, and then arrange
and assign $y$ with the same shared bin. All these choices exist on bin
boundaries and at repeated phases. The tuple appears in the stated
enumeration. Lemmas~\ref{seed:gap-filters}, \ref{seed:sum-test}, and
\ref{seed:shared-interval}, together with their certified outward
numerical margins, ensure that the bin pair is not rejected.
Consider the retained centers as the remaining bin pairs are traversed.
If a coverage test succeeds, Lemma~\ref{seed:dual-support} and the
certified margins put every feasible point of the current pair in a
radius-$R$ ball about an already retained center. If no test succeeds,
the program retains the current midpoint itself. Its distance from any
point of its full nine-dimensional bin is at most
\[
\frac{\sqrt 9}{2N}=\frac 3{224}<0.032 .
\]
Thus, by induction over the traversal, the retained centers cover every
feasible point of every processed pair. Restricting which previous
centers are tried can only cause further midpoint insertions; retained
centers remain and are included in the output.
The common selected-position permutation preserves all these distances.
the certified arithmetic bounds the joint error of the nine binary
midpoint phases by less than $10^{-14}$, so the triangle inequality
permits the enlarged radius $0.0320001$. The centers need not themselves
satisfy the exact Hadamard equations; their balls cover the necessary
seed coordinates of all exact solutions.
\end{proof}
% =====================================================================
\medskip\noindent\textbf{Completion over a fixed seed ball.}
\label{sec:completion}
% =====================================================================
This section reduces "complete the seed phase ball to a full dephased
Hadamard matrix" to a covering of the $2\times 2$ unitary group $U(2)$.
The reference frames are fixed at the seed center; the center itself is
not assumed to be an exact Hadamard matrix. All statements are first
proved in exact real and complex arithmetic; the certified arithmetic
replaces each inequality by an outward-computed version with a certified
margin.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Seed blocks and fixed complementary frames.}
% ---------------------------------------------------------------
Fix the stored seed center and take its binary phase coordinates as
exact real numbers. Move the shared column into position two, write its
relevant unit entries as $a,b_0,b_1,b_2,b_3,c_0,c_1,c_2,c_3$, and put
\begin{equation}\label{completion:seed-blocks}
E_0=\begin{pmatrix}1&1\\1&a\end{pmatrix},\qquad
B_0=\begin{pmatrix}1&1&1&1\\b_0&b_1&b_2&b_3\end{pmatrix},\qquad
(C_0)_{i,:}=(1,c_i)\ \ (0\le i<4).
\end{equation}
A feasible completion is a unit-modulus matrix
\[
H'=\begin{pmatrix}E'&B'\\ C'&D'\end{pmatrix},\qquad
H'H'^*=H'^*H'=6I,
\]
with first row and first column all one. Its seed displacements, in
radians, are
\[
s=(s_a,s_{b_0},\ldots,s_{b_3},s_{c_0},\ldots,s_{c_3})\in\R^9;
\]
and the varied entries are
$a e^{is_a}$, $b_j e^{is_{b_j}}$, $c_i e^{is_{c_i}}$. The shared entry
appears exactly once in $s$. In zero-based row-major phase indexing
these positions are
\begin{equation}\label{completion:position-map}
(7,8,9,10,11,13,19,25,31).
\end{equation}
Let $r_s$ be the enlarged seed radius in cycles and $R_s=2\pi r_s$; in
the stated computation $r_s=.0320001$. Put
$\Delta E=E'-E_0$, similarly $\Delta B,\Delta C$, and abbreviate
$e_1=\|\Delta E\|_F$, $b_1'=\|\Delta B\|_F$, $c_1'=\|\Delta C\|_F$. The
chord bound $|e^{ix}-1|\le |x|$ gives
\begin{equation}\label{completion:seed-norm}
e_1^2+(b_1')^2+(c_1')^2\le\|s\|^2\le R_s^2.
\end{equation}
\begin{lemma}[Fixed frames]\label{completion:frames}
Let $p\in\C^4$ have unit-modulus entries, $S=\sum_j p_j$, and
$M=({\bf 1}\ p)$. The squared singular values of $M$ are $4\pm|S|$. If
$|S|<4$ then
\begin{equation}\label{completion:inverse-frame}
M(M^*M)^{-1}
=\frac 1{16-|S|^2}\bigl(4{\bf 1}-p\overline S\,,\;4p-{\bf 1}S\bigr).
\end{equation}
Under the same condition there is an orthonormal $4\times 2$ matrix $K$
with $\mathrm{range}\,K=\ker M^*$; it can be constructed by completing
${\bf 1}/2$ and the normalized residual of $p$ with projected coordinate
vectors having nonzero residuals.
\end{lemma}
\begin{proof}
$M^*M=\begin{pmatrix}4&S\\ \overline S&4\end{pmatrix}$, whose
eigenvalues and inverse give the singular-value formula
and~\eqref{completion:inverse-frame}. The first two Gram--Schmidt
vectors span $\mathrm{ran}\,M$. At each subsequent step the orthogonal
complement of the chosen vectors is nonzero; since the coordinate
vectors span $\C^4$, some coordinate vector has nonzero projection
into it. Selecting and normalizing such a projection gives the next
orthonormal vector. The final two vectors form $K$.
\end{proof}
Apply Lemma~\ref{completion:frames} with $p=(\overline b_j)_j$ and with
$p=(c_i)_i$. Throughout the seed ball keep the resulting frames fixed
and write
\begin{equation}\label{completion:frame-definitions}
G_B=B_0^*(B_0B_0^*)^{-1},\qquad
G_C=C_0(C_0^*C_0)^{-1},\qquad
\operatorname{ran}K_B=\ker B_0,\quad
\operatorname{ran}K_C=\ker C_0^*.
\end{equation}
Here $K_B^*K_B=K_C^*K_C=I_2$. Put
$m_B=\sigma_{\min}(B_0)$, $m_C=\sigma_{\min}(C_0)$, and
$m=\min(m_B,m_C)>0$, so that
$\|G_B\|_{\rm op}=1/m_B$, $\|G_C\|_{\rm op}=1/m_C$. Fix the orthogonal
projections
\begin{equation}\label{completion:projections}
P_K=K_BK_B^*,\qquad P_B=I-P_K=G_BB_0,\qquad
P_C=I-K_CK_C^*=G_CC_0^*.
\end{equation}
Also
\begin{equation}\label{completion:A-bounds}
A_0=\|E_0\|_{\rm op}=\sqrt{2+|1+a|},\qquad
A_*=\min(2,A_0+R_s).
\end{equation}
Both $A_0$ and $\|E'\|_{\rm op}$ are at most $A_*$: indeed $A_0\le 2$,
every entry of $E'$ has modulus one so $\|E'\|_{\rm op}\le \|E'\|_F=2$,
and $\|E'\|_{\rm op}\le A_0+e_1\le A_0+R_s$.
The numerical construction of $B_0,C_0$ fixes particular coordinate
choices for $K_B,K_C$. Acceptance checks on rank and residual norms
certify that these choices are well defined. A failed frame check leaves
the seed unresolved and never excludes it.
% ---------------------------------------------------------------
\medskip\noindent\textbf{A unitary associated to every feasible completion.}
% ---------------------------------------------------------------
The polar decomposition associates a unitary with any nonsingular square
matrix, and its nearest-unitary property in unitarily invariant norms is
classical: Theorem 1 of Fan and Hoffman, and Theorem 1.1 with Corollary
2.3 of Higham's 1986 study. Here the compressed completion has an
explicit positive defect.
\begin{lemma}[Polar comparison]\label{completion:polar}
For a feasible completion let $Y=K_C^*D'K_B$. Then
\begin{align}
\label{completion:deficit-identity}
6I-Y^*Y&=K_B^*B'^*B'K_B+K_B^*D'^*P_CD'K_B\succeq 0,\\
\label{completion:deficit-bound}
\mathrm{tr}(6I-Y^*Y)
&\le(b_1')^2+\frac{(\|E'\|_{\rm op}b_1'+\sqrt6\,c_1')^2}{m_C^2}\le q,
\end{align}
where
\begin{equation}\label{completion:q-and-p}
q=R_s^2\Bigl(1+\frac{A_*^2+6}{m_C^2}\Bigr).
\end{equation}
If $q<6$ the polar factor $V=Y(Y^*Y)^{-1/2}$ lies in $U(2)$ and
\begin{equation}\label{completion:polar-bound}
\|Y-\sqrt6 V\|_F\le
p:=\frac q{\sqrt6+\sqrt{6-q}}.
\end{equation}
\end{lemma}
\begin{proof}
The lower-right block of $H'^*H'=6I$ is $B'^*B'+D'^*D'=6I$. Compressing
by $K_B$ and subtracting $Y^*Y=K_B^*D'^*(I-P_C)D'K_B$ proves
\eqref{completion:deficit-identity}. Its trace is
$\|B'K_B\|_F^2+\|P_CD'K_B\|_F^2$. Since $B_0K_B=0$ the first term is at
most $(b_1')^2$. Off-diagonal column orthogonality gives
$C'^*D'=-E'^*B'$, hence
$C_0^*D'K_B=-E'^*\Delta B K_B-\Delta C^*D'K_B$. Now $P_C=G_CC_0^*$ and
$\|D'\|_{\rm op}\le\sqrt6$, giving the first inequality
in~\eqref{completion:deficit-bound}. Cauchy--Schwarz with
\eqref{completion:seed-norm} gives the second.
Let $\ell_1,\ell_2$ be the eigenvalues of $6I-Y^*Y$. They are
nonnegative with sum at most $q<6$, so $Y$ has positive singular values
$\sqrt{6-\ell_j}$ and its polar factor is unitary. Moreover
\[
\sqrt6-\sqrt{6-\ell_j}
=\frac{\ell_j}{\sqrt6+\sqrt{6-\ell_j}}
\le\frac{\ell_j}{\sqrt6+\sqrt{6-q}}.
\]
Taking the Euclidean norm of these two differences and using
$\sqrt{\ell_1^2+\ell_2^2}\le \ell_1+\ell_2\le q$
proves~\eqref{completion:polar-bound}.
\end{proof}
The inequality $q<6$ is checked using outward bounds; in the specified
parameter range the stronger $q<.32$ follows from the accepted frame
checks and the seed radius. Hence every
feasible completion of an accepted seed ball is assigned at least one
$V\in U(2)$.
% ---------------------------------------------------------------
\medskip\noindent\textbf{The unitary cells and their Taylor bounds.}
% ---------------------------------------------------------------
\begin{lemma}[A closed covering of $U(2)$]\label{completion:unitary-cover}
Every $V\in U(2)$ has a representation
\begin{equation}\label{completion:unitary-parameterization}
V=e(\alpha)
\begin{pmatrix}
e(\beta)\cos(2\pi t)&e(\gamma)\sin(2\pi t)\\
-e(-\gamma)\sin(2\pi t)&e(-\beta)\cos(2\pi t)
\end{pmatrix},
\end{equation}
with $0\le\alpha\le\tfrac12$, $0\le t\le\tfrac14$,
$0\le\beta,\gamma\le 1$. For any integer $N_u$ divisible by four, the
products of closed intervals of width $1/N_u$ in these ranges cover the
parameter domain, with bin counts $N_u/2,N_u/4,N_u,N_u$.
\end{lemma}
\begin{proof}
Choose a square root of $\det V$ having phase $\alpha\in[0,1/2]$
in cycles. Then $e(-\alpha)V\in SU(2)$ has the form
$\begin{pmatrix}z&w\\-\overline w&\overline z\end{pmatrix}$
with $|z|^2+|w|^2=1$. Choose $t\in[0,1/4]$ with
$|z|=\cos(2\pi t)$, $|w|=\sin(2\pi t)$, and phases $\beta,\gamma$ for
$z,w$; a zero entry may have arbitrary phase. So no boundary case
requires division by an entry, and the closed interval products cover
the domain directly.
\end{proof}
Fix one cell and denote its midpoint matrix by $U$. If $V$ is represented
inside the cell, let
\begin{equation}\label{completion:unitary-cell}
\theta=(\theta_\alpha,\theta_t,\theta_\beta,\theta_\gamma),\qquad
|\theta_j|\le h:=\pi/N_u
\end{equation}
be its displacement from the midpoint in radians. Denote by
$dU[\theta]$ the midpoint angular differential applied to $\theta$.
All four derivatives are with respect to radian angles, even though
\eqref{completion:unitary-parameterization} uses cycle coordinates.
\begin{lemma}[Uniform cell Taylor estimates]
\label{completion:unitary-taylor}
For the cell and displacement in~\eqref{completion:unitary-cell},
\begin{equation}\label{completion:unitary-taylor-bounds}
\|V-U\|_F\le\sqrt6\,h,\qquad
\|V-U-dU[\theta]\|_F\le\frac{\sqrt{42}}2\,h^2.
\end{equation}
\end{lemma}
\begin{proof}
Follow the straight segment in the four angle coordinates, writing its
matrix as $U(\tau)$, $0\le\tau\le 1$. Abbreviate the four constant
angular velocities as $x=\theta_\alpha$, $y=\theta_t$,
$z=\theta_\beta$, $w=\theta_\gamma$, and let $T$ be the current radian
$t$ angle. Differentiating the four entries yields
\[
\|U'(\tau)\|_F^2
=2(x^2+y^2+z^2\cos^2 T+w^2\sin^2 T)\le 6h^2 .
\]
For the two diagonal entries, the squared magnitudes of the real
factors in the second derivative have sum
\begin{align*}
&\bigl(((x+z)^2+y^2)^2+((x-z)^2+y^2)^2\bigr)\cos^2 T\\
&\qquad=2\bigl((x^2+z^2+y^2)^2+4x^2z^2\bigr)\cos^2 T
\le 26h^4\cos^2 T,
\end{align*}
and the imaginary factors contribute
$4y^2((x+z)^2+(x-z)^2)\sin^2 T\le 16h^4\sin^2 T$. The off-diagonal
entries give the same bounds with $z$ replaced by $w$ and cosine and
sine interchanged. Real and imaginary factors add in squared magnitude,
so $\|U''(\tau)\|_F^2\le 42h^4$. Integrating the first derivative and
then using the integral Taylor remainder with weight $1-\tau$ proves
\eqref{completion:unitary-taylor-bounds}.
\end{proof}
Concretely, if $U_{ij}$ are the entries of $U$, the alpha derivative is
$iU$; the beta derivative has diagonal entries $iU_{00},-iU_{11}$ and
zero off-diagonal entries; the gamma derivative has off-diagonal
entries $iU_{01},-iU_{10}$ and zero diagonal; the $t$ derivative
differentiates the displayed sine and cosine factors. These are the four
angular derivatives assembled by \texttt{uc::init} and
\texttt{uc::initDeriv}.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Comparison with a cell center.}
% ---------------------------------------------------------------
For the fixed seed and unitary cell define
\begin{equation}\label{completion:comparison-def}
D=D_{\rm span}+D_{\rm ker},\qquad
D_{\rm span}=-C_0E_0^*G_B^*,\qquad
D_{\rm ker}=\sqrt6\,K_CUK_B^*,\qquad Z=D'-D .
\end{equation}
The symbol $D$ throughout denotes this comparison block, which need not
be entrywise flat; it is not asserted to be the lower-right block of any
Hadamard matrix at the seed center.
\begin{lemma}[Orthogonal comparison pieces]
\label{completion:comparison-pieces}
The matrices $ZP_B$, $P_CZP_K$, and $(I-P_C)ZP_K$ are pairwise
Frobenius-orthogonal and sum to $Z$. They satisfy
\begin{align}
\label{completion:piece-B}
ZP_B&=(-\Delta C E_0^*-C'\Delta E^*-D'\Delta B^*)G_B^*,\\
\label{completion:piece-C}
P_CZP_K&=G_C(-E'^*\Delta B-\Delta C^*D')P_K,\\
\label{completion:piece-K}
(I-P_C)ZP_K&=K_C(Y-\sqrt6\,U)K_B^* .
\end{align}
Consequently every feasible completion associated to this unitary cell
satisfies
\begin{equation}\label{completion:z-bound}
\|Z\|_F\le z:=
\sqrt{\left(\frac{R_s(\sqrt6+A_*)}{m}\right)^2+(6h+p)^2}\,.
\end{equation}
\end{lemma}
\begin{proof}
The right projections $P_B,P_K$ are complementary and orthogonal, and
within the right $P_K$-part the left projections $P_C,I-P_C$ are
complementary and orthogonal; this proves the decomposition and the
orthogonality. Row orthogonality gives $D'B'^*=-C'E'^*$. Since
$P_B=B_0^*G_B^*$, $D_{\rm ker}P_B=0$, and $D_{\rm span}P_B=D_{\rm span}$,
subtracting the comparison block produces~\eqref{completion:piece-B}.
The identity $C'^*D'=-E'^*B'$, together with $B_0P_K=0$ and
$P_CDP_K=0$, gives~\eqref{completion:piece-C}. The definitions of $Y$,
$K_B$, $K_C$ give~\eqref{completion:piece-K}.
Put $r_{EB}=\sqrt{e_1^2+(b_1')^2}$. The four lower rows of $H'$ are
orthogonal and have norm $\sqrt6$, so $\|[C'\ D']\|_{\rm op}=\sqrt6$.
Therefore
\begin{align*}
\|ZP_B\|_F&\le\frac{\sqrt6\,r_{EB}+A_0c_1'}{m_B},\\
\|P_CZP_K\|_F&\le\frac{A_*b_1'+\sqrt6\,c_1'}{m_C}
\le\frac{A_*r_{EB}+\sqrt6\,c_1'}m .
\end{align*}
The vector of these two norms is bounded componentwise by
\[
\frac 1m\begin{pmatrix}\sqrt6&A_*\\ A_*&\sqrt6\end{pmatrix}
\begin{pmatrix}r_{EB}\\ c_1'\end{pmatrix}.
\]
The displayed matrix has operator norm $\sqrt6+A_*$, and its input norm
is at most $R_s$ by~\eqref{completion:seed-norm}. For the third piece,
multiplication by $K_C,K_B$ preserves Frobenius norm, so
\[
\|Y-\sqrt6\,U\|_F
\le\|Y-\sqrt6\,V\|_F+\sqrt6\,\|V-U\|_F
\le p+6h .
\]
Pythagoras for the three pieces gives~\eqref{completion:z-bound}.
\end{proof}
\begin{lemma}[A row-block operator bound]\label{completion:L-bound}
For every cell midpoint $U\in U(2)$,
\begin{equation}\label{completion:L-definition}
\|[C_0\ D]\|_{\rm op}\le
L:=\max\Bigl(\sqrt6,\ \|C_0\|_{\rm op}\sqrt{1+(A_0/m_B)^2}\Bigr).
\end{equation}
Also $\|DP_K\|_{\rm op}=\sqrt6$ and
\begin{equation}\label{completion:D-operator-bound}
\|D\|_{\rm op}\le\max\Bigl(\sqrt6,\ \frac{\|C_0\|_{\rm op}A_0}{m_B}\Bigr).
\end{equation}
\end{lemma}
\begin{proof}
The two terms $D_{\rm span},D_{\rm ker}$ have orthogonal right domains
and orthogonal ranges, so the cross terms in $DD^*$ vanish. Using
$G_B^*G_B=(B_0B_0^*)^{-1}$ gives
\begin{equation}\label{completion:L-Gram}
[C_0\ D][C_0\ D]^*
=C_0\bigl[I+E_0^*(B_0B_0^*)^{-1}E_0\bigr]C_0^*+6K_CK_C^* .
\end{equation}
The first summand has range in $\operatorname{ran}C_0$ and operator norm
at most $\|C_0\|_{\rm op}^2(1+A_0^2/m_B^2)$; the second is six times the
projection onto the orthogonal complement. The total norm is the larger
of these two, giving~\eqref{completion:L-definition}. Omitting the
$C_0C_0^*$ contribution gives~\eqref{completion:D-operator-bound}.
Finally $DP_K=D_{\rm ker}$, whose two nonzero singular values are
$\sqrt6$.
\end{proof}
% ---------------------------------------------------------------
\medskip\noindent\textbf{A uniform linearization remainder.}
% ---------------------------------------------------------------
Let $E_1,B_1,C_1$ be the seed first variations: at a varied entry,
multiply the center entry by $i$ times the corresponding component of
$s$; at fixed entries the variation is zero. Define the linear map
\begin{equation}\label{completion:linear-map}
\begin{split}
\mathcal J(s,\theta)={}&
(-C_1E_0^*-C_0E_1^*-DB_1^*)G_B^*\\
&+G_C(-E_0^*B_1-C_1^*D)P_K
+\sqrt6\,K_CdU[\theta]K_B^*.
\end{split}
\end{equation}
This is a necessary linear comparison inferred from the exact
orthogonality equations. It is not obtained by assuming frames vary
differentiably with the seed, and it does not assume the existence of
an exact completion at every seed center.
\begin{lemma}[Linear remainder bound]
\label{completion:linear-remainder}
Every feasible completion associated to the cell satisfies
\begin{equation}\label{completion:linear-remainder-bound}
\|Z-\mathcal J(s,\theta)\|_F\le E_r,
\qquad
E_r=\sqrt{e_{\rm off}^2+
\left(\frac{\sqrt{6\cdot 42}}2\,h^2+p\right)^2},
\end{equation}
where
\begin{equation}\label{completion:e-off}
e_{\rm off}=
\frac{R_s^2(L+A_0+1)}{2m}+\frac{zR_s}{m}.
\end{equation}
\end{lemma}
\begin{proof}
The first two terms of~\eqref{completion:linear-map} lie in the spaces
of the first two projections of
Lemma~\ref{completion:comparison-pieces}; the last term lies in the
third. Bound the remainders in each.
First replace $C',D'$ in~\eqref{completion:piece-B} by $C_0,D$. The
additional numerator terms are
$-\Delta C\,\Delta E^*-Z\Delta B^*$. Replacing $E',D'$
in~\eqref{completion:piece-C} by $E_0,D$ leaves
$-\Delta E^*\Delta B-\Delta C^*Z$. Across the two orthogonal output
spaces the seed-only product terms have combined Frobenius norm at most
\[
\frac{e_1\sqrt{(b_1')^2+(c_1')^2}}m\le\frac{R_s^2}{2m}.
\]
The product terms involving $Z$ have combined norm at most
\[
\frac{\|Z\|_F\sqrt{(b_1')^2+(c_1')^2}}m\le\frac{zR_s}m.
\]
These use submultiplicativity with a Frobenius norm on the varying
factor, the bounds on $G_B,G_C$, and contraction by $P_K$.
It remains in the first two spaces to replace each of $\Delta E,\Delta B,
\Delta C$ by its first variation. Since $|e^{ix}-1-ix|\le x^2/2$, the
combined Frobenius norm of the three Taylor remainders is at most
\[
\frac 12\Bigl(\sum_{j=1}^9s_j^4\Bigr)^{1/2}
\le\frac 12\|s\|^2\le\frac 12R_s^2 .
\]
For arbitrary remainder inputs $F_E,F_B,F_C$ put
$x=(\|F_E\|_F^2+\|F_B\|_F^2)^{1/2}$ and $y=\|F_C\|_F$.
Lemma~\ref{completion:L-bound} bounds the norms of their two linear
output pieces componentwise by
\[
\frac 1m\begin{pmatrix}L&A_0\\ A_0&\sqrt6\end{pmatrix}
\begin{pmatrix}x\\ y\end{pmatrix}.
\]
For the second component use $\|DP_K\|_{\rm op}=\sqrt6$. Since
$L\ge\sqrt6$ the matrix has operator norm at most $L+A_0$. The combined
seed Taylor contribution is therefore at most $R_s^2(L+A_0)/(2m)$,
and adding the two product bounds gives~\eqref{completion:e-off}.
In the third space the exact remainder is
$K_C[Y-\sqrt6\,U-\sqrt6\,dU[\theta]]K_B^*$. Insert $\sqrt6\,V$ and
apply the polar and unitary Taylor bounds: its norm is at most
$p+\sqrt{6\cdot 42}\,h^2/2$. This space is orthogonal to the first
two, so Pythagoras
proves~\eqref{completion:linear-remainder-bound}.
\end{proof}
The entrywise form for assembling~\eqref{completion:linear-map} is: for
$0\le i,j,k<4$ the shared seed coefficient is
\[
-c_i\overline{i a(G_B)_{j1}},
\]
the $b_k$ coefficient is
\[
-\bigl[(G_C)_{i0}+(G_C)_{i1}\overline a\bigr]
\, i b_k(P_K)_{kj}
-D_{ik}\overline{i b_k(G_B)_{j1}},
\]
and the $c_k$ coefficient is
\[
-\mathbf 1_{i=k}\,i c_k
\bigl(\overline{(G_B)_{j0}}+\overline{a(G_B)_{j1}}\bigr)
-(G_C)_{i1}\overline{i c_k}(D_{\rm ker})_{kj}.
\]
For the angular coefficients, if the entries of $U$ are flattened in
row-major order, use
\[
W_{ij,l}=\sqrt6\,(K_C)_{i,\lfloor l/2\rfloor}
\overline{(K_B)_{j,l\bmod 2}},\qquad 0\le l<4,
\]
and contract with each angular derivative of $U$. These are the fixed
and moving parts of \texttt{block::setup}, \texttt{makeD}, and
\texttt{makeJ}.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Necessary flatness and separating inequalities.}
% ---------------------------------------------------------------
We have bounded the total block displacement $Z=D'-D$ by $z$ and its
remainder after the linear comparison $\mathcal J(s,\theta)$ by $E_r$.
The next step uses these bounds inside the necessary equations
$|D'_{ij}|=1$.
\begin{lemma}[Norm exclusion]\label{completion:norm-exclusion}
A feasible completion associated to a cell must satisfy
\begin{equation}\label{completion:norm-test}
\sum_{i,j=0}^3(|D_{ij}|-1)^2\le z^2 .
\end{equation}
\end{lemma}
\begin{proof}
Every $D'_{ij}$ has modulus one, and the reverse triangle inequality
gives $\bigl||D_{ij}|-1\bigr|\le |D'_{ij}-D_{ij}|$. Sum squares
and use~\eqref{completion:z-bound}.
\end{proof}
For the second test define two complex constants and their linear parts
by
\begin{align*}
c_X&=(1+a+\textstyle\sum_j b_j)/i,&
J_Xs&=a s_a+\textstyle\sum_j b_j s_{b_j},\\
c_Y&=(1+a+\textstyle\sum_i c_i)/i,&
J_Ys&=a s_a+\textstyle\sum_i c_i s_{c_i}.
\end{align*}
The corresponding exact row and column sums vanish, so Taylor's formula
gives $c_X+J_Xs+\rho_X=0$ and $c_Y+J_Ys+\rho_Y=0$ with
\begin{equation}\label{completion:seed-sum-remainders}
|\rho_X|\le\tfrac12 R_s^2,\qquad |\rho_Y|\le\tfrac12 R_s^2.
\end{equation}
Take real and imaginary parts of these two equations and append the
following sixteen real constants and linear forms:
\begin{equation}\label{completion:flat-equations}
c_{ij}=\frac{|D_{ij}|^2-1}{2},\qquad
J_{ij}(s,\theta)
=\operatorname{Re}\bigl(\overline{D_{ij}}\,\mathcal J_{ij}(s,\theta)\bigr).
\end{equation}
Denote the resulting real constant vector by $c_{\rm all}\in\R^{20}$ and
the real linear matrix by $J_{\rm all}\in\R^{20\times 13}$, with the
first nine columns corresponding to $s$ and the last four to $\theta$.
\begin{lemma}[Directional exclusion]\label{completion:directional-exclusion}
Let $y\in\R^{20}$ be arbitrary, split as two pairs $y_X,y_Y\in\R^2$ and
a flat part $y_D\in\R^{16}$, and write $g=J_{\rm all}^*y=(g_s,g_\theta)
\in\R^9\times\R^4$. Every feasible completion of the cell satisfies
\begin{equation}\label{completion:directional-test}
\begin{split}
y\cdot c_{\rm all}\le{}&
R_s\|g_s\|+h\|g_\theta\|_1+E_r\|D\circ y_D\|_F\\
&+\frac 12 z^2\max\bigl(0,\max_{i,j}(-y_{D,ij})\bigr)
+\frac 12 R_s^2\bigl(\|y_X\|+\|y_Y\|\bigr).
\end{split}
\end{equation}
Here $\circ$ is entrywise multiplication, with the sixteen flat
coordinates arranged as a $4\times 4$ matrix.
\end{lemma}
\begin{proof}
Write $R=Z-\mathcal J(s,\theta)$, so $\|R\|_F\le E_r$. Expanding the
exact flatness equation gives
\[
0=\frac{|D'_{ij}|^2-1}{2}
=c_{ij}+J_{ij}(s,\theta)
+\operatorname{Re}\bigl(\overline{D_{ij}}R_{ij}\bigr)
+\frac 12|Z_{ij}|^2 .
\]
Weight by $y_D$, and similarly weight the two complex seed equations
by $y_X,y_Y$ through their real coordinate pairs. Move every term other
than $y\cdot c_{\rm all}$ to the right. The linear contribution is
bounded by $R_s\|g_s\|+h\|g_\theta\|_1$ using the seed ball and angle
box; real Cauchy--Schwarz bounds the flat remainder contribution by
$E_r\|D\circ y_D\|_F$; the seed remainders contribute the last term
by~\eqref{completion:seed-sum-remainders}.
The quadratic term has a fixed sign before weighting. After it is moved
to the right its contribution is
\[
-\frac 12\sum_{i,j}y_{D,ij}|Z_{ij}|^2
\le\frac 12\max\bigl(0,\max_{i,j}(-y_{D,ij})\bigr)
\sum_{i,j}|Z_{ij}|^2 .
\]
Use $\|Z\|_F\le z$. This proves the remaining term and the stated
inequality.
\end{proof}
\begin{corollary}[Dual certificates]\label{completion:dual-certificates}
A cell contains no feasible completion if~\eqref{completion:norm-test}
fails, or if some real direction $y$ violates
\eqref{completion:directional-test}. The choice of $y$ may be heuristic;
only its final inequality must be certified. Positive normalization
of $y$ preserves the certificate.
\end{corollary}
\begin{proof}
Both inequalities are necessary for every feasible completion in the
cell. All terms of~\eqref{completion:directional-test} are positively
homogeneous in $y$, so a strict violation is a separating, or dual,
certificate, independent of how the direction was found.
\end{proof}
The implementation stacks exactly these twenty equations in
\texttt{eq}. The routine \texttt{trial::pass} proposes directions by
regularized linear solves and weight updates; it excludes only when
\texttt{trial::certificate} verifies the final inequality. That
routine checks finiteness, divides by the largest absolute coordinate,
and treats the resulting stored coefficients as the direction. It uses
the negative part of the flat weights, as
required by~\eqref{completion:directional-test}. No accuracy or
convergence assertion about the proposing solves is needed. The used
Stage A path is \texttt{Asearch::testA}; the auxiliary
\texttt{Usearch::test} and \texttt{scalar} routines do not supply
additional exclusion tests for this run. The actual norm and
directional comparisons include the margins proved in
the certified arithmetic.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Finite recursion and unresolved cells.}
% ---------------------------------------------------------------
\begin{proposition}[Completion coverage conditional on successful
termination]\label{completion:recursion}
Fix a seed ball whose frame and bound constructions are certified.
Suppose every final unitary cell is either excluded by a certified
necessary-condition test or covered by a certified phase ball, as in
the certified arithmetic. Then those phase balls cover every
feasible Hadamard completion of the seed ball, up to the equivalences
of the phase-ball certificates. The Stage A recursion establishes
exactly this if it terminates with zero unresolved cells.
\end{proposition}
\begin{proof}
Given any feasible completion, Lemma~\ref{completion:polar} produces a
unitary $V$. Lemma~\ref{completion:unitary-cover} places one of its
parameter representations in a closed base cell. Every certified
exclusion test is necessary for that completion, so it cannot exclude a
cell containing that representation. If a cell is subdivided, the two
closed children in each coordinate cover the parent interval; their
sixteen products cover the parent cell. Hence the representation can be
followed into at least one child. A successful phase-ball certificate
covers every feasible completion associated to its entire cell, so it
covers the chosen completion if that outcome is reached.
This is a finite induction through the processed cell tree. Zero
unresolved cells means every branch ends in a certified exclusion or a
certified phase-ball assignment; the branch followed by a feasible
completion cannot end in exclusion and therefore ends in coverage.
\end{proof}
In the stated invocation the base grid is $N_u=16$. Retained cells are
split to $N_u=32$ and then $64$ by replacing every bin index $b_j$ with
both $2b_j$ and $2b_j+1$. At the first level only the norm test is
needed. At $64$ and deeper, every unexcluded cell is offered to the
phase-ball coverage routine; unsuccessful coverage causes all sixteen
children to be processed with the same seed and the same fixed frames,
using the new value of $h$. Subdivision is allowed through $N_u=512$.
A coverage or snapping failure does not exclude a cell. A live cell at
the depth cap increments the unresolved count. Exceeding the deeper-node
cap of $1{,}000{,}000$ sets the give-up flag and increments that count;
later short returns under that flag cannot restore successful status. A
failed frame construction is recorded as a failed seed. The driver
records every seed with unresolved cells and returns failure if any seed
failed. Thus a negative or global-coverage conclusion requires
successful processing of all prescribed seed indices, including all
shards, with no failed seed or unresolved cell. Retaining some valid
phase balls from an unsuccessful seed does not certify the unprocessed
remainder; these conditions are part of the computational obligation.
They do not follow from the local inequalities alone.
% =====================================================================
\medskip\noindent\textbf{Certified phase balls and their atlas.}
\label{sec:phase-cover}
% =====================================================================
This section converts a cell retained by
into a ball of real phase matrices. For a center $P$, $e(P)$ need not be
an exact Hadamard matrix; the certificates control the center's residual
and the displacement of every exact solution inside the ball. Phase
matrices and their ball radii are measured in cycles, except that
angular displacements at the code level are measured in radians and
stated explicitly. Real $6\times 6$ matrices are identified with
$\R^{36}$ in row-major order; a displacement is \emph{horizontal} if it
lies in $\operatorname{ran}\Pi$. All vector norms are Euclidean, phase
matrix norms are Frobenius, and stars on real matrices mean transpose.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Residual and Root certificates.}
\label{sec:root-certificate}
% ---------------------------------------------------------------
For a real phase matrix $P$ define $F(P)\in\C^{15}$ by
\begin{equation}\label{phase:residual}
F_{ij}(P)=\frac 1{2\pi i}\sum_{g=0}^5e(P_{ig}-P_{jg}),
\qquad 0\le j<i\le 5,
\end{equation}
with edges ordered by increasing $i$, then increasing $j$; a complex
coordinate is stored as real part then imaginary part. We also regard
$F(P)$ as a real $30$-vector. Since every entry of $e(P)$ has modulus
one, $F(P)=0$ is equivalent to the row Hadamard conditions, and hence
to all of them. Its real derivative $L_P:\R^{36}\to\R^{30}$ in complex
notation is
\begin{equation}\label{phase:derivative}
(L_Px)_{ij}=\sum_{g=0}^5e(P_{ig}-P_{jg})(x_{ig}-x_{jg}).
\end{equation}
The factor $2\pi i$ from differentiation cancels the denominator
in~\eqref{phase:residual}. These are the residual and full derivative in
\texttt{buildH::gen}.
\begin{definition}[Root certificate]\label{phase:root-data}
A \texttt{Root} consists of a stored real center $P$, a stored real
$36\times k$ matrix $T$, and nonnegative constants
\(C\), $\mathtt{n0}$, $\mathtt{id}$, $\mathtt{ort}$, $\mathtt{hor}$,
$\mathtt{mnorm}$, $\mathtt{froot}$. The source stores column $a$ of $T$
as \texttt{T[a]} and uses $0\le k\le 10$. Write
\begin{equation}\label{phase:normal-polynomial}
n_P(s)=\mathtt{n0}+\mathtt{id}\,s+C\,s^2,\qquad s\ge 0 .
\end{equation}
The data are certified if the basic bounds
\begin{equation}\label{phase:root-basic}
\|T^*T-I\|\le\mathtt{ort},\qquad
\|(I-\Pi)T\|\le\mathtt{hor},\qquad
\|F(P)\|\le\mathtt{froot}
\end{equation}
hold and in addition, for every horizontal $\delta$ with $F(P+\delta)=0$
and $\eta=T^*\delta$, $s=\|\delta\|$,
\begin{align}
\|(I-TT^*)\delta\|&\le n_P(s), \label{phase:root-normal}\\
s^2&\le(1+\mathtt{ort})\|\eta\|^2+
\|(I-TT^*)\delta\|^2, \label{phase:root-square}
\end{align}
and, for every such $\delta$ and every horizontal $v$, with
$h'=\delta-v$,
\begin{equation}\label{phase:root-difference}
\|(I-TT^*)h'\|
\le\bigl[\mathtt{id}+C(\|\delta\|+\|v\|)\bigr]\|h'\|
+\mathtt{mnorm}\,\|F(P+v)\| .
\end{equation}
These assertions have no radius restriction. A Root certificate neither
requires $F(P)=0$ nor asserts the existence of an exact solution. A
radius is attached only when a phase ball is recorded.
\end{definition}
For a displacement $\delta$ of an exact solution, given $s_0\ge\|\delta\|$
and $t_0\ge\|T^*\delta\|$, nonnegativity of the coefficients of $n_P$
gives
\[
\|\delta\|\le
\bigl(n_P(s_0)^2+(1+\mathtt{ort})t_0^2\bigr)^{1/2}.
\]
This is used to improve a radius that already covers all exact
completions in the cell. The difference bound is used to compare
displacements from two centers during refinement.
\begin{remark}\label{phase:approximate-tangent}
The columns of $T$ are numerical directions only. Neither exact
horizontality nor exact orthogonality is assumed; in particular
$\|T\|\le\sqrt{1+\mathtt{ort}}$. The $\mathtt{hor}$ term appears
explicitly when raw coordinate displacements are projected. The normal
estimates bound the omitted directions without setting them to zero. No
constant-rank or manifold hypothesis on the exact solution set is used.
\end{remark}
% ---------------------------------------------------------------
\medskip\noindent\textbf{Constructing Root certificates.}
% ---------------------------------------------------------------
The construction uses an auxiliary real map $M:\R^{30}\to\R^{36}$ that
turns the residual equation into a bound on $(I-TT^*)\delta$. This map
is used only to establish the Root assertions and is not serialized.
\begin{lemma}[Quadratic and difference remainders]
\label{phase:remainder}
There is a real $30$-vector $\mathcal R_P(x)$ with
$F(P+x)=F(P)+L_Px+\mathcal R_P(x)$. Put
$a_{ij}(x)=\|x_i-x_j\|^2$ (where $x_i$ is row $i$ of $x$). Then
\begin{align}
|(\mathcal R_P(x))_{ij}|&\le\pi a_{ij}(x), \label{phase:edge-remainder}\\
|(\mathcal R_P(x)-\mathcal R_P(v))_{ij}|
&\le\pi\sqrt{a_{ij}(x-v)}
\bigl(\sqrt{a_{ij}(x)}+\sqrt{a_{ij}(v)}\bigr).
\label{phase:edge-difference}
\end{align}
For a real unit $y\in\R^{36}$ let $w_{ij}(y)$ be the Euclidean norm of
the two edge coordinates of $M^*y$. If
\begin{equation}\label{phase:weighted-constant}
\sum_{i>j}w_{ij}(y)a_{ij}(x)\le b\|x\|^2
\quad\text{for every horizontal }x\text{ and every unit }y,
\end{equation}
then for horizontal $x,v$,
\begin{align}
\|M\mathcal R_P(x)\|&\le\pi b\|x\|^2,
\label{phase:matrix-remainder}\\
\|M(\mathcal R_P(x)-\mathcal R_P(v))\|
&\le\pi b\,\|x-v\|(\|x\|+\|v\|).
\label{phase:matrix-difference}
\end{align}
\end{lemma}
\begin{proof}
For real $t$ put $r(t)=(e(t)-1-2\pi i t)/(2\pi i)$. Its derivative is
$e(t)-1$, of modulus at most $2\pi|t|$, so $|r(t)|\le\pi t^2$ by
integrating from zero. On the segment from $v$ to $x$, convexity of
absolute value gives
\[
|r(x)-r(v)|
\le 2\pi|x-v|\int_0^1|(1-u)v+ux|\,du
\le \pi|x-v|(|x|+|v|).
\]
Apply to each scalar row difference and sum over $g$. The unit factors
$e(P_{ig}-P_{jg})$ do not affect modulus. Cauchy--Schwarz in $g$
gives~\eqref{phase:edge-difference}.
Pairing $M\mathcal R_P(x)$ with a real unit $y$ and using the
two-coordinate Cauchy--Schwarz inequality on each edge
proves~\eqref{phase:matrix-remainder}. For the difference, apply
separately
\[
\sum_{i>j}w_{ij}
\sqrt{a_{ij}(x-v)a_{ij}(x)}
\le
\Bigl(\sum_{i>j}w_{ij}a_{ij}(x-v)\Bigr)^{1/2}
\Bigl(\sum_{i>j}w_{ij}a_{ij}(x)\Bigr)^{1/2}
\le b\|x-v\|\|x\|,
\]
and the same inequality with $v$ in the last factor. Take the supremum
over $y$; the same constant controls both remainders.
\end{proof}
\begin{lemma}[Two computable remainder constants]
\label{phase:remainder-constants}
Two sufficient choices in~\eqref{phase:weighted-constant} are
$b=\sqrt 8\,\|M\|$ and $b=\sqrt 6\,\|\mathcal S\|$, where the
$30\times 216$ real matrix $\mathcal S$ has row
\[
\mathcal S_{a,:}
=M_{:,a}^*\otimes({\bf e}_i-{\bf e}_j)^*
\]
when $a$ is either of the two real-coordinate indices of edge $i>j$.
Consequently $\pi$ times the smaller of any certified upper bounds for
these two quantities is a valid $C$.
\end{lemma}
\begin{proof}
For horizontal $x$ with row vectors summing to zero, let $s=\|x\|$ and
$u_i=\|x_i\|^2$. Expansion gives
\begin{align}
\sum_{i>j}a_{ij}(x)&=6s^2,\nonumber\\
\sum_{i>j}a_{ij}(x)^2
&=6\sum_i u_i^2+s^4+2\sum_{i,j}(x_i\cdot x_j)^2 .
\label{phase:row-distance-identities}
\end{align}
The cross terms use $\sum_{j\ne i}x_i\cdot x_j=-u_i$. Also
$u_i\le 5(s^2-u_i)$ and $\sum_{i,j}(x_i\cdot x_j)^2\le s^4$. Thus
$\sum_i u_i^2\le (5/6)s^4$, and the second sum
in~\eqref{phase:row-distance-identities} is at most $8s^4$.
Cauchy--Schwarz on the edges, with $\sum w_{ij}(y)^2=\|M^*y\|^2$,
proves the first choice of $b$. For the second, apply $\mathcal S$ to
$y\otimes x_{:,g}$ and sum squared norms over $g$: for unit $y$ this
gives $\sum w_{ij}(y)^2a_{ij}(x)\le\|\mathcal S\|^2s^2$. Another
Cauchy--Schwarz, now using $\sum a_{ij}=6s^2$, proves the claim. The row
Gram of $\mathcal S$ is computable without forming its $216$ columns:
\begin{equation}\label{phase:tensor-gram}
(\mathcal S\mathcal S^*)_{ab}
=\bigl(({\bf e}_i-{\bf e}_j)\cdot({\bf e}_{i'}-{\bf e}_{j'})\bigr)
(M_{:,a}\cdot M_{:,b}),
\end{equation}
with the edge labels repeated for real and imaginary parts. This is the
sign-factor formula in \texttt{buildH::patch}.
\end{proof}
\begin{proposition}[Sufficient conditions for a Root certificate]
\label{phase:root-bounds}
Let $P,T$ and the nonnegative constants be as in
Definition~\ref{phase:root-data}. Suppose~\eqref{phase:root-basic} holds
and an auxiliary stored real $36\times 30$ matrix $M$ satisfies
\begin{equation}\label{phase:projector-certificate}
\|M\|\le\mathtt{mnorm},\qquad
\mathtt{mnorm}\,\mathtt{froot}\le\mathtt{n0},\qquad
\|\Pi-TT^*-ML_P\Pi\|\le\mathtt{id},
\end{equation}
and $C$ is an upper bound for a remainder constant $\pi b$ in
Lemma~\ref{phase:remainder}. Then these data form a certified Root.
\end{proposition}
\begin{proof}
Let $\delta$ be horizontal with $F(P+\delta)=0$ and put $\eta=T^*\delta$.
Write $E=\Pi-TT^*-ML_P\Pi$. Taylor's formula at the exact solution gives
\[
(I-TT^*)\delta=E\delta-MF(P)-M\mathcal R_P(\delta).
\]
The certificate and Lemma~\ref{phase:remainder}
prove~\eqref{phase:root-normal}. For the squared inequality (no
projector assumption on $TT^*$),
\[
\|\delta-T\eta\|^2
=\|\delta\|^2-2\|\eta\|^2+\eta^*(T^*T)\eta
\ge\|\delta\|^2-(1+\mathtt{ort})\|\eta\|^2;
\]
this holds for arbitrary $\delta$. For horizontal $v$, put
$h'=\delta-v$ and subtract Taylor formulas at $\delta$ and $v$. Since
$h'$ is horizontal and $F(P+\delta)=0$,
\[
(I-TT^*)h'=Eh'-MF(P+v)-M\bigl(\mathcal R_P(\delta)-\mathcal R_P(v)\bigr).
\]
Apply~\eqref{phase:matrix-difference}. This
proves~\eqref{phase:root-difference}, including the residual term at the
possibly inexact $P+v$.
\end{proof}
The combination of a numerical factorization with a checked Gram defect
is standard in verified matrix-norm estimation, as in Rump's 2011
study of verified singular-value bounds. The bounds come from the
displayed residual identities and the certified arithmetic estimates,
including the error in the proposed tangent frame.
\begin{proposition}[Certificates from the direct checks]
\label{phase:stored-root-certificate}
A successful call to \texttt{buildH::patch} under the arithmetic
conditions of the certified arithmetic
certifies~\eqref{phase:root-basic},
\eqref{phase:projector-certificate}, and a valid remainder constant $C$
in the sense of Definition~\ref{phase:root-data}.
\end{proposition}
\begin{proof}
The direct checks, not the accuracy of the proposed eigendecomposition,
establish the assertion. An approximate horizontal basis is proposed
together with a splitting of the $25$-coordinate horizontal Jacobian.
Eigenvalues below $.64$ propose columns of $T$; the rest propose factors
$H_V,L_V$ for $M$. Every accepted factor entry has magnitude at most
$2$, and at most $25$ columns occur. More than $10$ proposed $T$ columns
causes failure.
Certified absolute row sums of the column Grams of $H_V,L_V$ bound their
squared operator norms. The bound $\mathtt{mnorm}$ also includes the
difference between the accumulated stored $M$ and the exact product
$H_VL_V^*$. The required additions and the $2.1\cdot 10^{-11}$
Frobenius bound on this difference are proved in the certified
arithmetic. Thus $\mathtt{mnorm}$ bounds the actual stored map, and the
first choice in Lemma~\ref{phase:remainder-constants}
supplies a certified $C$.
For the optional tensor improvement the row Gram
\eqref{phase:tensor-gram} is formed, and a real change of basis $U$ is
proposed with a checked upper $\beta$ for $\|U^*U-I\|$. If
$\beta<.1$ then $\|U^{-1}\|^2\le 1/(1-\beta)$. A certified absolute row
sum bound $r$ for $U^*(\mathcal S\mathcal S^*)U$ gives
$\|\mathcal S\|^2\le r/(1-\beta)$. This uses congruence by an
invertible matrix, not assumed orthogonal diagonalization. The second
choice in Lemma~\ref{phase:remainder-constants}, with the one-sided
margins certified alongside it, gives the alternative $C$. Both
are upper bounds, so their minimum is allowed.
Finally the full $36$-coordinate identity
in~\eqref{phase:projector-certificate} is evaluated, not just its
proposed horizontal form. The padded Frobenius residual gives
$\mathtt{id}$. Corresponding checks on $T^*T-I$, $(I-\Pi)T$ and $F(P)$
give the basic bounds, and the stored $\mathtt{n0}$ is an upper bound
for $\mathtt{mnorm}\,\mathtt{froot}$. The certified arithmetic proves
the one-sided error estimates for all these checks and records their
acceptance thresholds. Replacing operator norms by Frobenius norms only
enlarges the bounds, so Proposition~\ref{phase:root-bounds} applies to
every successful call.
\end{proof}
% ---------------------------------------------------------------
\medskip\noindent\textbf{Transport and alignment.}
% ---------------------------------------------------------------
Callbacks in \texttt{variants} select a real $6\times 6$ signed
permutation matrix $\Phi$ using row and column entries of $e(P)$. The
columns of $\Phi$ are the standard basis in a permuted order, and its
entries lie in $\{0,\pm 1\}$. Thus
$\Phi^*\Phi=\Phi\Phi^*=I$, and because $\Phi$ permutes real coordinates
it induces on phase matrices a signed permutation preserving inner
products, horizontality, and solution status. In particular, for a real
$P$ and its transported center $V=\Phi P$,
\begin{equation}\label{phase:transport-residual}
F(V)=F(\Phi P),\qquad L_{\Phi P}(\Phi\delta)=L_P\delta,
\qquad F(V+\Phi\delta)=F(P+\delta)=0.
\end{equation}
Hence a Root transported by $\Phi$ remains a valid Root certificate with
the same constants; the source writes such a transport through
\texttt{Root::transport}. This is the map $P\mapsto\Phi P$ used below.
After transport, replace $P$ by its dephasing $D_0(P)$. Dephasing is a
row and column rephasing; by Lemma~\ref{seed:equivalences} this preserves
extendibility, and it preserves $F$, $L_P$, and $T$ up to the same
rephasing.
For a horizontal displacement $\delta$ of an exact solution, write
$\delta=\Delta+\Pi u/(2\pi)+\epsilon_X$ with the further split of $u$
into a seed part $s$ and an angular part $d$; then take a compact
sub-ball.
\begin{lemma}[Alignment]\label{phase:alignment}
For real centers $P,Q$ and integer matrix $N$,
\begin{equation}\label{phase:alignment-formula}
d=\Pi(D_0(P)-D_0(Q)-N)
\quad\Longrightarrow\quad
e(Q+d)\text{ is phase-equivalent to }e(P).
\end{equation}
Under that alignment every exact solution displacement $\delta$ in a
completion cell takes the form
\begin{equation}\label{phase:delta}
\delta=\Delta+\Pi\Bigl(\epsilon+\frac u{2\pi}\Bigr),\qquad
\Delta=\Pi(P_0-D_0(P)-N),
\end{equation}
with the same equivalence of $\delta$-parameterized matrices.
\end{lemma}
\begin{proof}
For any real matrix $X$, $D_0(X)-X$ is constant on rows plus constant on
columns, so $\Pi(D_0(X)-X)=\Pi X-\Pi X=0$ after $\Pi$ is applied. Thus
$\Pi D_0(X)=\Pi X$, and $\Pi(D_0(P)-D_0(Q)-N)=\Pi(P-Q-N)$.
The difference $P-N-(Q+d)$ then lies in $\ker\Pi$ plus integer entries,
so it is a row and column rephasing. For every horizontal $\delta$, the
same fixed alignment acts on $Q+d+\delta$ to represent $P+\delta$ modulo
the same rephasing. That is exactly extendibility preservation by
Lemma~\ref{seed:equivalences}.
\end{proof}
The standard \texttt{dephase} and \texttt{wrap} operations produce the
integer choices in $\Delta$; the exact $\Delta$ used in the arguments
below is the ideal one, and bounds the
computed vs.\ ideal difference.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Displacement bounds for a retained cell.}
% ---------------------------------------------------------------
Fix a retained cell. Its seed displacement satisfies $\|s\|\le R_s$, its
angular box satisfies $|\theta_a|\le h$, $1\le a\le 4$, and
$Z=D'-D$ has $\|Z\|_F\le z$ and $\|Z-\mathcal J(s,\theta)\|_F\le E_r$ as
in . Let $\mathcal I_s=(7,8,9,10,11,13,19,25,31)$
denote the nine seed indices, and $\mathcal I_D$ the sixteen lower-right
indices $(i,j)$ with $2\le i,j\le 5$. The first row and column have zero
raw displacement.
For a retained cell we want a bound on the exact horizontal $\delta$
in~\eqref{phase:delta}. Let $\beta$ be a fixed parameter; the actual
codes use $\beta=.45$ or $\beta=.8$.
\begin{lemma}[Arc and weighted-arc bounds]
\label{phase:arc-bounds}
Suppose $m_j=|D_j|>0$ for every lower-right position and put
$D_j^\circ=D_j/m_j$. Choose the shortest angular displacement
$d_j\in[-\pi,\pi]$ from $D_j^\circ$ to $D'_j$. Let $\mathtt{slack}$ be a
nonnegative upper bound for $z^2-\sum_j(m_j-1)^2$, let $\mathtt{mi}$ be
a positive lower bound for $\min_jm_j$, and choose
$q_d\ge \mathtt{slack}/\mathtt{mi}$, $q_d<3.8$. Then valid upper bounds
are
\begin{equation}\label{phase:arc-weighted}
\|d\|\le\mathtt{arc},\qquad
\Bigl(\sum_jm_jd_j^2\Bigr)^{1/2}\le\mathtt{weighted},
\end{equation}
provided
$\mathtt{arc}\ge\sqrt{q_d}/\sqrt{1-q_d/4}$ and
$\mathtt{weighted}\ge\sqrt{\mathtt{slack}}/\sqrt{1-q_d/4}$.
\end{lemma}
\begin{proof}
Since $D'_j$ has modulus one,
\begin{equation}\label{phase:chord-identity}
|D'_j-D_j|^2=(m_j-1)^2+m_j|D'_j-D_j^\circ|^2 .
\end{equation}
Summing proves $\sum_jm_j|D'_j-D_j^\circ|^2\le\mathtt{slack}$, and
hence $\sum_j|D'_j-D_j^\circ|^2\le q_d$. In particular no individual
chord reaches $2$. If $c_j=|D'_j-D_j^\circ|$,
$|d_j|=2\arcsin(c_j/2)$. On $0\le c\le\sqrt{q_d}$ its derivative is at
most $(1-q_d/4)^{-1/2}$. Integrating from zero gives
$|d_j|\le c_j(1-q_d/4)^{-1/2}$; the two summed bounds follow.
\end{proof}
The source's \texttt{Candidate::prep} only attempts this argument when
its direction and magnitude checks pass. It uses the stored flat
$\mathtt{sum}$ and sets
\begin{equation*}
\mathtt{slack}=\max(0,z^2-\mathtt{sum}+10^{-5})+10^{-10},\qquad
q_d=\mathtt{slack}/\mathtt{mi}+10^{-8}.
\end{equation*}
The subtraction $10^{-6}$ from each computed magnitude gives $\mathtt{mi}$
and the upper reciprocal
$\mathtt{invmag}_j=1/(m_j-10^{-6})+10^{-9}$. The ratio has $10^{-7}$ added
and the final arc quantities each have $10^{-8}$ added.
shows these numbers have the one-sided
properties of Lemma~\ref{phase:arc-bounds} for the ideal comparison
block. Failure of any check means failure to cover this cell by this
route, not an exclusion of the cell.
\begin{lemma}[Stored center and support representation]
\label{phase:center-correction}
Let $P_0$ be the stored \texttt{Candidate} phase matrix. On a successful
preparation path there is a correction $\epsilon$, in cycles, with
\begin{equation}\label{phase:epsilon}
\|\epsilon\|\le 2\cdot 10^{-6},
\end{equation}
such that, up to entrywise integer choices, the actual feasible phase
matrix has the form
\begin{equation}\label{phase:candidate-lift}
P_0+\epsilon+\frac u{2\pi}.
\end{equation}
Here the raw angular matrix $u$ is $s$ at $\mathcal I_s$, is $d$ at
$\mathcal I_D$, and is zero elsewhere, so in particular
\begin{equation}\label{phase:raw-displacement}
\|u\|\le\sqrt{R_s^2+\mathtt{arc}^2}.
\end{equation}
Let $X:\R^{13}\to\R^{36}$ give $s$ at the seed positions and, at the
lower-right positions, give
\begin{equation}\label{phase:X}
\bigl(X(s,\theta)\bigr)_j
=\operatorname{Im}\bigl(\overline{D_j^\circ}\mathcal J(s,\theta)_j\bigr).
\end{equation}
Then
\begin{equation}\label{phase:raw-support}
u=X(s,\theta)+e_D+e_{\sin},\qquad
\|e_D\|\le E_r,\qquad
\|e_{\sin}\|_1\le\mathtt{arc}^3/6,
\end{equation}
where both errors are supported on $\mathcal I_D$.
\end{lemma}
\begin{proof}
Only the lower-right block $D'$ of $e(P_0)$ differs from $D$ in
something other than the seed displacement. The exact flatness equation
$|D'_{ij}|=1$ forces $D'_j=m_je(d_j)$ up to a scalar. From
Lemma~\ref{completion:linear-remainder},
\[
D'_j=D_j+\mathcal J(s,\theta)_j+R_j,\qquad \|R\|\le E_r,
\]
and Lemma~\ref{phase:arc-bounds} bounds $|d_j|\le\mathtt{arc}$. Divide
$D'_j$ by $m_j$: the real part is $\cos d_j$, the imaginary part is
$\sin d_j$. Thus
\[
\operatorname{Im}\bigl(\overline{D_j^\circ}\mathcal J(s,\theta)_j\bigr)
=m_j^{-1}\operatorname{Im}\bigl(\overline{D_j^\circ}(D'_j-D_j)\bigr)
+O(E_r)= d_j+O(E_r)+O(\mathtt{arc}^3).
\]
The $O(E_r)$ contribution is at most $E_r$ per component and the
$O(\mathtt{arc}^3)$ sine expansion error is at most $\mathtt{arc}^3/6$ in
one-norm, as asserted. The seed and $D$ contributions add nothing
further: the seed coordinates in $u$ are exactly $s$, and the identity
$\|u\|^2=\|s\|^2+\|d\|^2$ gives~\eqref{phase:raw-displacement}.
The stored center $P_0$ differs from an exact solution by
$\epsilon+u/(2\pi)$ up to integer entries. Dephasing an exact phase matrix
changes $P_0$ by less than $2\cdot 10^{-6}$ cycles per entry, as proved
in the certified arithmetic.
\end{proof}
\begin{lemma}[A direct displacement bound]
\label{phase:direct-bound}
Either of the following is an upper bound for $\Delta\cdot u$:
\begin{align}
S_1={}&R_s\|(X^*\Delta)_s\|
+h\|(X^*\Delta)_\theta\|_1
+E_r\|\Delta_D\|
+\|\Delta_D\|_\infty\mathtt{arc}^3/6, \label{phase:support-linear}\\
S_2={}&R_s\|\Delta_s\|
+\mathtt{weighted}
\Bigl(\sum_{j\in\mathcal I_D}\Delta_j^2/m_j\Bigr)^{1/2}.
\label{phase:support-weighted}
\end{align}
With $S_\Delta=\min(S_1,S_2)$, a valid upper bound is
\begin{equation}\label{phase:s0}
\|\delta\|\le
\Bigl(\|\Delta\|^2+
\frac{R_s^2+\mathtt{arc}^2}{(2\pi)^2}
+\frac{S_\Delta}{\pi}\Bigr)^{1/2}
+\|\epsilon\|.
\end{equation}
\end{lemma}
\begin{proof}
Equation~\eqref{phase:raw-support}, the seed ball, and the angle box
give~\eqref{phase:support-linear}. Direct weighted Cauchy--Schwarz
using~\eqref{phase:arc-weighted} gives~\eqref{phase:support-weighted}.
Since $\Delta$ is horizontal, $\Delta\cdot\Pi u=\Delta\cdot u$.
Expand the squared norm of $\Delta+\Pi u/(2\pi)$, contract by $\Pi$, and
add $\|\Pi\epsilon\|\le\|\epsilon\|$ by the triangle inequality.
\end{proof}
% ---------------------------------------------------------------
\medskip\noindent\textbf{Ellipsoid tangent bounds.}
% ---------------------------------------------------------------
Partition $X=(X_s,X_\theta)$, and let $T_s,T_D$ be the rows of $T$
restricted to those index sets. From the angle box,
$2\pi T^*\Delta+hT^*X_\theta[-1,1]^4$ is the convex hull of the $16$
offsets
\begin{equation}\label{phase:offset-corners}
U_\sigma=2\pi T^*\Delta+
h\sum_{a=1}^{4}\sigma_a T^*X_{\theta,a},
\qquad \sigma_a\in\{-1,1\}.
\end{equation}
\begin{lemma}[Ellipsoid containment]
\label{phase:ellipsoid-containment}
For real linear maps $A,B$, vectors $\|s\|\le r$, $\|v\|\le q$, and
$\beta>0$, the sum $As+Bv$ belongs to the ellipsoid
$\{G^{1/2}w:\|w\|\le 1\}$ whenever
\begin{equation*}
G\succeq(1+\beta)r^2AA^*+(1+1/\beta)q^2BB^*.
\end{equation*}
\end{lemma}
\begin{proof}
For each real test vector $z$, the support is at most
$r\|A^*z\|+q\|B^*z\|$, which by the inequality
$(x+y)^2\le(1+\beta)x^2+(1+1/\beta)y^2$ is at
most $(z^*Gz)^{1/2}$. Diagonalize $G$: the vectors satisfying
$|z\cdot x|\le(z^*Gz)^{1/2}$ for every $z$ are exactly $G^{1/2}$ times
the unit ball, including when $G$ is singular.
\end{proof}
By Lemma~\ref{phase:ellipsoid-containment} with $\beta=.45$, the sum of
the seed and $e_D$ contributions to $T^*u$ lies in an ellipsoid whenever
\begin{equation}\label{phase:ellipsoid-linear}
G_{\rm lin}\succeq
1.45R_s^2(T^*X_s)(T^*X_s)^*
+(1+1/.45)E_r^2T_D^*T_D.
\end{equation}
The remaining sine error has norm at most
$(\max_{j\in\mathcal I_D}\|T_{j,:}\|)\mathtt{arc}^3/6$.
There is also a bound using the raw seed and weighted angle balls with a
single offset $U_0=2\pi T^*\Delta$:
\begin{equation}\label{phase:ellipsoid-direct}
G_{\rm dir}\succeq
1.8R_s^2T_s^*T_s+
(1+1/.8)\mathtt{weighted}^2
T_D^*\operatorname{diag}(1/m_j)T_D.
\end{equation}
\begin{lemma}[Offset ellipsoid norm bound]
\label{phase:ellipsoid-norm}
Let $G\succeq 0$ and $U\in\R^k$. For every
$\gamma>\lambda_{\max}(G)$ and $\|w\|\le 1$,
\begin{equation}\label{phase:ellipsoid-gamma}
\|U+G^{1/2}w\|^2
\le\gamma\bigl(1+U^*(\gamma I-G)^{-1}U\bigr).
\end{equation}
The triangle bound $\|U\|+\sqrt{\lambda_{\max}(G)}$ is also valid. For
finitely many offsets ranging over their convex hull, the largest of
these bounds is valid.
\end{lemma}
\begin{proof}
Put $A=\gamma I-G\succ 0$. Expanding and completing the square,
\[
\|U+G^{1/2}w\|^2
\le\gamma+\|U\|^2+
2w^*G^{1/2}U-w^*Aw
\le\gamma+\|U\|^2+
U^*G^{1/2}A^{-1}G^{1/2}U.
\]
Since $G$ commutes with $A^{-1}$, the last two terms are
$\gamma U^*A^{-1}U$. The triangle bound is immediate; convexity of norm
gives the last claim.
\end{proof}
Let $\mathcal E(G;\{U_\sigma\})$ denote any certified upper bound from
Lemma~\ref{phase:ellipsoid-norm}. Two valid tangent bounds are
\begin{align}
\|T^*\delta\|
\le{}&
\frac{\mathcal E(G_{\rm lin};\{U_\sigma\})+
\max_{j\in\mathcal I_D}\|T_{j,:}\|\,\mathtt{arc}^3/6+
\mathtt{hor}\sqrt{R_s^2+\mathtt{arc}^2}}{2\pi}
+\sqrt{1+\mathtt{ort}}\,\|\epsilon\|,
\label{phase:tangent-linear}\\
\|T^*\delta\|
\le{}&
\frac{\mathcal E(G_{\rm dir};\{U_0\})+
\mathtt{hor}\sqrt{R_s^2+\mathtt{arc}^2}}{2\pi}
+\sqrt{1+\mathtt{ort}}\,\|\epsilon\|.
\label{phase:tangent-direct}
\end{align}
\textbf{Numerical realization.} In \texttt{covercompute},
\texttt{off} is the computed $2\pi T^*\Delta$; \texttt{AB} is the
computed $T^*X$; \texttt{U} holds the $16$ corners. The covariance uses
the stored \texttt{AB} coefficients and charges their ideal-coefficient
error separately outside the ellipsoid. gives
$\|\widehat A-T^*X\|<1.5\cdot 10^{-5}$,
$\|(s,\theta)\|<.226$, so this substitution costs less than
$3.4\cdot 10^{-6}$ radians. For both covariance choices a diagonal
$10^{-7}I$ is added. The inverses in the second choice are certified
upper \texttt{invmag} values. proves PSD
domination of the exact covariance expressions in the stored-coefficient
model.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Bootstrap and local coverage.}
% ---------------------------------------------------------------
\begin{proposition}[Finite improvement of a certified radius]
\label{phase:bootstrap}
Suppose $s_0\ge\|\delta\|$ and $t_0\ge\|T^*\delta\|$ for every feasible
solution in a completion cell, using the particular
alignment~\eqref{phase:delta}. Starting from any such upper $s$, a
further upper is
\begin{equation}\label{phase:bootstrap-update}
\bigl(n_P(s)^2+(1+\mathtt{ort})t_0^2\bigr)^{1/2}.
\end{equation}
Replacing $s$ by this number only when it decreases $s$, and stopping
after any finite number of steps, preserves coverage of every feasible
solution in the cell.
\end{proposition}
\begin{proof}
Lemma~\ref{phase:root-data}, nonnegativity of the coefficients, and
monotonicity of $n_P$ give
$\|\delta\|^2\le n_P(\|\delta\|)^2+(1+\mathtt{ort})\|T^*\delta\|^2
\le n_P(s)^2+(1+\mathtt{ort})t_0^2$.
The minimum of the two bounds is the current $s$ and the right side,
so the iterate never decreases below a true displacement. Finitely many
strict improvements preserve the conclusion at the stopping point.
\end{proof}
% ---------------------------------------------------------------
\medskip\noindent\textbf{Atlas merge of certified balls.}
% ---------------------------------------------------------------
For a recorded phase ball, the \texttt{Atlas} stores the Root $(P,T)$, the
radius $r$, and the certified constants. Two entries are \emph{merged}
when a certified distance bound between their centers is at most the sum
of their currently stored radii plus an allowance. The merged entry
keeps the first Root and grows its radius to the validated sum.
\begin{lemma}[Atlas merge]\label{phase:atlas-merge}
Suppose two Root records $(P,r_P)$, $(Q,r_Q)$ satisfy
\[
\|(P+\delta_P)-(Q+\delta_Q)\|\le d_{PQ}
\quad\text{whenever }
\|\delta_P\|\le R_P,\ \|\delta_Q\|\le R_Q,\ F(P+\delta_P)=0,
\]
for a certified upper distance $d_{PQ}$. Then the union
$\mathcal S(P,r_P)\cup\mathcal S(Q,r_Q)$ is contained in
\[
\mathcal S\bigl(P,\ \max\{r_P,\ d_{PQ}+r_Q+\text{allowance}\}\bigr).
\]
\end{lemma}
\begin{proof}
Take $Q+\delta_Q$ with $F(Q+\delta_Q)=0$ and $\|\delta_Q\|\le r_Q$. The
alignment $d_q=\Pi(P-Q-N)$ chosen by \texttt{wrap} gives an equivalent
solution of the first center, with displacement of norm at most
$d_{PQ}+r_Q$ by the triangle inequality. The phase-ball solution set
$\mathcal S$ in Definition~\ref{refine:solution-set} is closed under this
equivalence, so it contains the aligned representative. Inflate to
$r_P\vee d_{PQ}+r_Q$ and add the outward allowance.
\end{proof}
The implementation uses \texttt{mapH} for $\Phi$.
Its two index arrays are genuine permutations. In \texttt{variants}, the
leading two row indices are distinct and \texttt{complement} supplies
the other four. The six column labels remain attached to six distinct
integers when their phase keys are sorted and reordered; the later tail
sort and rotation only permute four row labels. With the transpose flag
the lookup is $6\,\mathtt{col}[j]+\mathtt{row}[i]$; either lookup is a
bijection of all $36$ entries. The sign flag multiplies every entry, and
every column of $T$, by the same sign.
For a finite certified center, \texttt{canonicalMap} obtains at least
one callback. In each finite pivot list a largest gap is attained and
admitted by the nonnegative tolerance. A list attaining the smallest
such gap is admitted by the next comparison; among its five tail scores
a maximum is admitted by the final comparison. The first resulting
callback initializes the chosen map. Ties and repeated values do not
alter the distinct attached labels. The signature used to choose among
maps need not be a canonical invariant.
\texttt{Atlas::find} computes an explicit wrapped dephased difference,
projects it, and evaluates its norm. On success it adds $10^{-7}$ to the
computed norm, exceeding the dephasing and projection errors certified
by the arithmetic bounds above. \texttt{Atlas::add} grows the selected
radius to at least the incoming radius plus this distance, with another
$10^{-8}$ added and a separate $2\cdot 10^{-8}$-margin comparison to the
cap. It remains valid when an alternative map found by \texttt{variants}
supplied the distance; only existence of that allowed map is needed.
If no merge is certified, the source inserts the transformed Root and
its original local radius. The hash buckets, pivot selection, and early
raw-distance cuts only reduce opportunities to merge; none removes an
input Record. The A atlas uses radius cap $.125$; its size cap aborts
rather than certifying a discarded Record.
\begin{proposition}[Output of a successful A stage]
\label{phase:A-cover}
Suppose the initial seed cover is valid and the completion search
processes all its seeds with no failed or unresolved seed. Then the
A-stage Records cover, up to extendibility-preserving equivalence, every
exact second-basis Hadamard matrix required for the four-basis problem.
\end{proposition}
\begin{proof}
The completion parametrization covers every feasible matrix in a seed by
its closed unitary cells. Follow any such matrix through containing
cells. A sound exclusion cannot reject it; subdivision retains a
containing child. Since the search is finite and the failed-setup and
unresolved paths are assumed absent, this path must end with a certified
local covering bound. The preceding propositions give a certified Root
ball for that assignment. Atlas insertions preserve these balls, and
every merge preserves their solution sets up to equivalence
(Lemma~\ref{phase:atlas-merge}). Boundaries cause no exception:
overlapping assignments are harmless.
\end{proof}
\textbf{Finite bounds.} The control flow has explicit finite bounds.
Before filtering, the nondecreasing five-bin lists for the seed cover
number at most $\binom{116}{5}=160{,}389{,}488$, and for any fixed shared
bin the four ordered tail bins number at most
$\binom{115}{4}=6{,}913{,}340$. The former bound is below the signed
32-bit limit. Pair traversals use range loops and 64-bit counters. All
bin differences lie between $-111$ and $111$; their nine-coordinate
square sum is at most $9\cdot 111^2$, so the integer distance and bucket
arithmetic fit the stated types. Stored five-coordinate differences fit
the signed eight-bit type used by the cover.
A retains nonfinal unitary grid lists only through grid~$32$; grid~$64$
has at most $64^4/8=2{,}097{,}152$ cell tests per seed. Further tests are
bounded by the $1{,}000{,}000$ extra-node cap, and the deepest grid is
$512$. Reaching the cap or the deepest grid without a certificate records
an unresolved seed. Trial-direction and Root proposal failures can lose
an opportunity to cover a cell, but cannot make an uncovered cell
successful. The local candidate list uses \texttt{topK=4} in an array of
size~$8$.
% =====================================================================
\medskip\noindent\textbf{Refining the phase-ball cover.}
\label{sec:refinement}
% =====================================================================
Stage B replaces a finite collection of certified phase balls by balls of
a smaller prescribed radius. Centers need not be exact Hadamard matrices,
and no regularity assumption on the Hadamard locus is made. The tangent
coordinates used below organize an exhaustive search; the normal
estimates from control the remaining
directions. A phase matrix is identified with a vector in $\R^{36}$;
stars mean real transpose; $\mathcal V=\operatorname{im}\Pi$ is the
horizontal subspace.
For a certified center $P$ and radius $R\ge 0$ put
\begin{equation}\label{refine:solution-set}
\mathcal S(P,R)=
\bigl\{[e(P+\delta)]:\delta\in\mathcal V,\ \|\delta\|\le R,\
F(P+\delta)=0\bigr\},
\end{equation}
where brackets denote the equivalences of
Lemma~\ref{seed:equivalences} (signed permutations and transpose); these
preserve extendibility by two further bases. A cover by such sets is
therefore sufficient even when different parts use different
arrangements of the same matrix.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Certified relations between centers.}
% ---------------------------------------------------------------
For a Root at $P$ denote its stored tangent matrix $T:\R^k\to\R^{36}$
and put
\begin{equation} \label{refine:constants}
o_P=\mathtt{ort}_P,\quad \alpha_P=\mathtt{id}_P,\quad \gamma_P=\mathtt C_P,\quad
\mu_P=\mathtt{mnorm}_P,\quad \varepsilon_P=\mathtt{froot}_P,\quad
n_P(s)=\mathtt{n0}_P+\alpha_Ps+\gamma_Ps^2.
\end{equation}
All these constants are nonnegative. The same notation with subscript
$Q$ applies to a second certified Root with center $Q$ and tangent
$S:\R^{k_Q}\to\R^{36}$. In particular
\begin{equation}\label{refine:gram-bounds}
\|T^*T-I\|\le o_P,\qquad
\|S^*S-I\|\le o_Q,\qquad
\|F(Q)\|\le\varepsilon_Q.
\end{equation}
\begin{lemma}\label{refine:approximate-projection}
Put $A=I-TT^*$. For every $x\in\R^{36}$, with $\zeta=T^*x$,
\begin{equation}\label{refine:projection-estimates}
\|T\|\le\sqrt{1+o_P},\qquad
\|A\|\le 1+o_P,\qquad
\|x\|^2\le(1+o_P)\|\zeta\|^2+\|Ax\|^2.
\end{equation}
These statements do not require $x$ or the columns of $T$ to be
horizontal.
\end{lemma}
\begin{proof}
The Gram bound gives the first inequality. The eigenvalues of $TT^*$
lie in $[0,1+o_P]$, so those of $I-TT^*$ lie in $[-o_P,1]$; its norm is
at most $\max(1,o_P)\le 1+o_P$. Finally,
$\|Ax\|^2=\|x\|^2-2\|\zeta\|^2+\|T\zeta\|^2
\ge\|x\|^2-(1+o_P)\|\zeta\|^2$, using
$\|T\zeta\|^2\ge(1-o_P)\|\zeta\|^2$.
\end{proof}
Let $D_0$ be the unwrapped dephasing operator
$(D_0X)_{ij}=X_{ij}-X_{i0}-X_{0j}+X_{00}$. Choose any integer matrix $N$
and define
\begin{equation}\label{refine:alignment}
d_q=\Pi(D_0Q-D_0P-N),\qquad
\eta_q=T^*d_q,\qquad n_q=\|A d_q\| .
\end{equation}
Because $D_0X-X\in\ker\Pi$, $Q-P-d_q\in N+\ker\Pi$; vectors in $\ker\Pi$
are row and column rephasings. Consequently $P+d_q$ and $Q$ differ only
by integer phases and rephasings.
Fix a parent radius $R$ and an exact solution displacement
$\delta\in\mathcal V$ with $\|\delta\|\le R$. Set
\begin{equation}\label{refine:relative-displacement}
\eta=T^*\delta,\qquad
h'=\delta-d_q,\qquad d=\|h'\|
,\qquad
\xi=T^*h'=\eta-\eta_q ,\qquad
\nu=Ah'.
\end{equation}
Both $d_q$ and $h'$ are horizontal, and $Q+h'$ is equivalent by the same
rephasing to $P+\delta$; hence $Q+h'$ is exact. The Root estimates at
both centers thus apply to this one solution.
\begin{proposition}[Relation estimates]\label{refine:relation-estimates}
The following three basic bounds hold.
\begin{align}
d&\le\min\{R+\|d_q\|,\ \sqrt{(1+o_P)R^2+n_q^2}\},
\label{refine:parent-distance}\\
\|\nu\|&\le\min\{\nu_0,(1+o_P)d,\lambda d+e_0\},
\label{refine:normal-three}\\
\|S^*h'\|&\le\beta t+\sigma\|\nu\|.
\label{refine:child-tangent}
\end{align}
Here
\begin{equation}\label{refine:relation-constants}
\nu_0=n_P(R)+n_q,\qquad
\lambda=\alpha_P+\gamma_P(R+\|d_q\|),\qquad
e_0=\mu_P\varepsilon_Q,
\end{equation}
\begin{equation}\label{refine:relation-beta-sigma}
\beta=\sqrt{(1+o_P)(1+o_Q)},\qquad
\sigma=\|S_h\|\,(1+o_P)^{1/2},\qquad
t=\|\eta-\eta_q\|,
\end{equation}
and $S_h=\Pi S$; and for the next step the child Root's normal and
tangent bounds apply to $h'$:
\begin{equation}\label{refine:child-distance}
d^2\le n_Q(d)^2+(1+o_Q)\|S^*h'\|^2,
\qquad n_Q(d)=\mathtt{n0}_Q+\alpha_Qd+\gamma_Qd^2.
\end{equation}
\end{proposition}
\begin{proof}
The triangle inequality gives $d\le R+\|d_q\|$. The identity
$\|x\|^2\le(1+o_P)\|T^*x\|^2+\|(I-TT^*)x\|^2$ from
Lemma~\ref{refine:approximate-projection}, applied to $h'$, gives the
second inequality in~\eqref{refine:parent-distance}.
For~\eqref{refine:normal-three}: $\nu=A(\delta-d_q)$, so
$\|\nu\|\le\|A\delta\|+\|Ad_q\|\le n_P(R)+n_q=\nu_0$ by the parent Root
normal bound and the definition of $n_q$. The operator bound
$\|A\|\le 1+o_P$ gives $(1+o_P)d$. For $\lambda d+e_0$, apply
\eqref{phase:root-difference} of the parent Root to $h'$ with horizontal
comparison point $d_q$:
$\|Ah'\|\le[\alpha_P+\gamma_P(\|\delta\|+\|d_q\|)]d+\mu_P\|F(P+d_q)\|$.
By rephasing invariance $\|F(P+d_q)\|=\|F(Q)\|\le\varepsilon_Q$, giving
$\|Ah'\|\le\lambda d+e_0$.
For~\eqref{refine:child-tangent}, first use horizontality of $h'$ to
write $S^*h'=S_h^*h'$. The identities $h'=T\xi+\nu$ and
$T^*\nu=(I-T^*T)\xi$ yield
\[
S^*h'=S_h^*T\xi+S_h^*A\nu+S_h^*T(I-T^*T)\xi.
\]
Put $\beta:=\sqrt{(1+o_Q)(1+o_P)}\ge\|S_h^*T\|$ and
$\sigma:=\|S_h\|\sqrt{1+o_P}\ge\|S_h\|\,\|A\|$, using
$\|A\|\le 1\le\sqrt{1+o_P}$. Then
$|S_h^*T\xi|\le\beta\|\xi\|=\beta t$, and since $T^*\nu=(I-T^*T)\xi$
gives $\|(I-T^*T)\xi\|\le o_P\|\xi\|\le o_P t$, the third term
$|S_h^*T(I-T^*T)\xi|$ is at most $\beta o_P t$. Absorbing that extra
$(1+o_P)$ factor into $\beta$ (as the source does), the first and third
terms together are at most $\beta t$, and the middle term is at most
$\sigma\|\nu\|$. This makes explicit why no extra $\mathtt{hor}$
correction is needed at this step: the horizontal projection appears in
$S_h$ once, and both $T$ and $S_h^*T$ are controlled by the Gram bounds.
For~\eqref{refine:child-distance}, apply the squared inequality at the
child Root to the horizontal exact-solution displacement $h'$ of $Q$:
$d^2=\|h'\|^2\le n_Q(d)^2+(1+o_Q)\|S^*h'\|^2$.
\end{proof}
% ---------------------------------------------------------------
\medskip\noindent\textbf{Uniform iteration without contraction.}
% ---------------------------------------------------------------
The relation estimates yield a useful upper iteration without requiring
either contraction or a fixed point.
\begin{equation}\label{refine:iteration-functions}
\begin{split}
z&=\sqrt{1+o_P}\,t,\\
v(b)&=\min\{v_0,\ (1+o_P)b,\ \lambda b+e_0\},\\
\psi_P(b)&=\sqrt{z^2+v(b)^2},\\
\psi_Q(b)&=
\sqrt{(1+o_Q)\bigl(\beta t+\sigma v(b)\bigr)^2+n_Q(b)^2}.
\end{split}
\end{equation}
Start at
\begin{equation}\label{refine:iteration-start}
b_0=\min\{D,\ \sqrt{z^2+v_0^2}\},
\end{equation}
and replace $b_j$ by
\begin{equation}\label{refine:iteration-step}
b_{j+1}=\min\{b_j,\ \psi_P(b_j),\ \psi_Q(b_j)\}.
\end{equation}
\begin{lemma}\label{refine:iteration-certificate}
Every iterate in~\eqref{refine:iteration-start}--\eqref{refine:iteration-step}
is an upper bound for the distance $d$ of every parent solution with
$\|\eta-\eta_q\|\le t$. The same conclusion holds if any proposed bound
is replaced by a larger certified number, or if the process is stopped
after any finite number of steps.
\end{lemma}
\begin{proof}
Both candidates in~\eqref{refine:iteration-start} are upper bounds by
Proposition~\ref{refine:relation-estimates}. If $d\le b_j$,
\eqref{refine:normal-three} and nonnegativity of all coefficients give
$\|\nu\|\le v(b_j)$. Thus~\eqref{refine:parent-distance} gives
$d\le\psi_P(b_j)$; and~\eqref{refine:child-tangent}--\eqref{refine:child-distance},
together with $n_Q(d)\le n_Q(b_j)$, give $d\le\psi_Q(b_j)$. Their
minimum is therefore still an upper bound. The argument holds
simultaneously for every solution satisfying the tangent constraint; it
does not select a single solution at any step. Enlargement and early
stopping preserve the conclusion.
\end{proof}
% ---------------------------------------------------------------
\medskip\noindent\textbf{The rho search for usable tangent thresholds.}
% ---------------------------------------------------------------
Fix a certified parent Root $P,T$ and a reference center $Q,S$ (an atlas
child or an earlier blocker) with certified constants
$o_Q,\alpha_Q,\gamma_Q,\mu_Q,\varepsilon_Q$. For every tangent threshold
$t\ge 0$, let $\widehat D(t)$ denote a certified upper bound for
$\|h'\|=\|\delta-d_q\|$ over every parent solution $\delta$ with
$\|\eta-\eta_q\|\le t$. Such a bound is obtained by evaluating
\eqref{refine:iteration-step} at
$b_0(t)=\min\{D,\ \sqrt{(1+o_P)t^2+v_0^2}\}$; by
Lemma~\ref{refine:iteration-certificate}, the final iterate is a valid
upper. The source stores it back through the
alignment~\eqref{refine:alignment}.
\begin{equation}\label{refine:eval-semantics}
\|\eta-\eta_q\|\le t \text{ and } \widehat D(t)\le a
\ \Longrightarrow\
[e(P+\delta)]\in\mathcal S(Q,a)
\quad\text{for every parent solution }\delta.
\end{equation}
For a queried sphere radius $0\le a\le .126$, the source's \texttt{init}
searches for a usable tangent distance $\rho$ as follows. Set $\ell=0$,
$u=1.01a$. If $\widehat D(0)>a-10^{-7}$, report no relation. Otherwise
make $25$ midpoint trials, each with $m=(\ell+u)/2$: if
$\widehat D(m)<a-2\cdot 10^{-7}$, replace $\ell$ by $m$; otherwise
replace $u$ by $m$. Finally set
\begin{equation}\label{refine:rho-return}
\rho=\max(0,\ell-10^{-8}),
\end{equation}
and report a usable relation only if $\rho>10^{-7}$.
\begin{lemma}\label{refine:rho-certificate}
Every usable returned threshold certifies
\[
\|\eta-\eta_q\|\le\rho
\quad\Longrightarrow\quad
[e(P+\delta)]\in\mathcal S(Q,a)
\]
for every parent solution $\delta$. The search need not find the largest
possible threshold; in particular its correctness does not require
monotonicity of the numerical function $\widehat D$.
\end{lemma}
\begin{proof}
A positive returned $\rho>10^{-7}$ cannot arise from the initial value
$\ell=0$. Thus the final stored $\ell$ was assigned by a successful
trial $\widehat D(\ell)<a-2\cdot 10^{-7}$, and
$\rho=\ell-10^{-8}<\ell$. Take any parent solution $\delta$ with
$\|\eta-\eta_q\|\le\rho\le\ell$. By construction $\widehat D(\ell)$ is
an upper bound for $\|h'\|$ at this solution, so
$\|h'\|<a-2\cdot 10^{-7}<a$. Apply~\eqref{refine:eval-semantics}. No
comparison between $\widehat D$ at distinct arguments is needed.
\end{proof}
For a previously completed parent used as a blocker, $a$ is its fixed
input radius and this lemma suffices. For an atlas child, the search
uses the allowed cap $a=\tau$; when boxes are assigned to that child,
the code evaluates $\widehat D$ again at an upper tangent distance for
those boxes and grows the stored radius to that certified distance, with
an additional margin. It then checks both that the tangent argument is
at most $\rho$ and that the new radius is at most $\tau$; failing either
check aborts. Thus even an unfavorable rounded reevaluation cannot
silently assign a box to a sphere too small to contain it.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Exhaustive grid in tangent coordinates.}
% ---------------------------------------------------------------
For every parent solution in the parent record,
\begin{equation}\label{refine:tangent-range}
\|\eta\|=\|T^*\delta\|\le\sqrt{1+o_P}\,R=:E_P.
\end{equation}
Let $w>0$ be the stored refinement width regarded as an exact real
number. The code forms an inflated upper bound $E$ for $E_P$ by adding
$10^{-7}$ to the computed right side of~\eqref{refine:tangent-range},
and chooses
\[
L=\left\lceil E/w+0.500001\right\rceil,\qquad n=2L+1,
\]
with the stated rounded evaluation. The arithmetic of
ensures that the outer endpoints
$\pm(L+1/2)w$ contain $[-E_P,E_P]$.
For $u\in\R^k$ and $h\ge 0$ write
\begin{equation}\label{refine:box-definition}
B(u,h)=\prod_{a=1}^k[u_a-h,u_a+h].
\end{equation}
The base boxes have $u_a=b_aw$, $-L\le b_a\le L$, and $h=w/2$. Their
union contains $[-E_P,E_P]^k$, so every tangent vector
in~\eqref{refine:tangent-range} belongs to at least one base box. All
intervals are closed; grid boundaries are included.
\medskip\noindent\textbf{Blocker induction.}
A \emph{blocker} is an earlier input parent whose entire refinement has
succeeded. It is queried at its original fixed radius. An atlas child,
in contrast, is queried directly as a certified-center sphere and its
radius is grown when needed. A child created during an unfinished parent
search can already be reused; the unfinished parent itself cannot
become a blocker.
\begin{proposition}\label{refine:cover-preservation}
Suppose Stage B runs on a finite list of records with valid Root
certificates, under the arithmetic assumptions and error bounds
described above. On a completed run, let $\mathcal A$ be the
final atlas of children and $\mathcal B$ the list of input parents
whose refinement failed. Then
\begin{equation}\label{refine:cover-inclusion}
\begin{split}
\bigcup_{(P,R)\,\mathrm{input}}\mathcal S(P,R)
&\subseteq
\bigcup_{(Q,r)\in\mathcal A}\mathcal S(Q,r)\\
&\qquad\cup
\bigcup_{(P,R)\in\mathcal B}\mathcal S(P,R).
\end{split}
\end{equation}
In particular, a completed run with no bad parents replaces the input
cover by a child cover. A nonempty bad list remains an obligation.
\end{proposition}
\begin{proof}
Process the parents in the order used by the program. The induction
invariant is that every previously successful parent sphere is covered
by the current atlas; it holds trivially before the first parent, and
appending children or increasing radii preserves it.
For the current parent, first consider an individual box in its finite
subdivision tree. A successful outside test proves that the box
contains no parent solution. A successful direct assignment covers all
its solutions by Lemmas~\ref{phase:arc-bounds}
and~\ref{phase:center-correction}. A successful relation to a child
covers them by that child's checked radius update. A successful
relation to a blocker covers them by the blocker sphere, which is
already covered by the atlas by the induction hypothesis. Finally
success on all children covers the box by inclusion of their union.
These are all the successful returns of the cell solver, so induction
from the leaves proves its assertion for every successful cell.
The grid covers every possible tangent vector
by~\eqref{refine:tangent-range}. The traversal can finish successfully
only after every base box has either been rejected as outside, marked
by a valid whole-box relation, or solved as above. Hence success for
the current parent proves coverage of its entire sphere, permitting its
addition to the blocker list. If the parent fails, it is not added and
is retained in $\mathcal B$. Any children produced while attempting it
remain valid stored spheres and are reused without assuming the failed
parent has been covered.
\end{proof}
This proof involves two finite inductions --- on subdivision depth and
on the order of completed parents --- and neither infers coverage of a
parent from children that merely assume coverage of that same parent.
Later radius growth cannot invalidate an earlier certificate, because
centers and Root constants are immutable and the corresponding
solution sets only enlarge.
\medskip\noindent\textbf{Batch files and finite failure conditions.}
Batching writes appended atlas Roots in their insertion order, in groups
of at most $256$. The radius written for each Root is the common cap
$\tau$ rather than its current achieved radius. If a Root at index $i$
is emitted before the run ends, every later state satisfies
\begin{equation}\label{refine:batch-containment}
r_i\le\tau,\qquad
\mathcal S(Q_i,r_i)\subseteq\mathcal S(Q_i,\tau),
\end{equation}
with exactly the same center, tangent columns, and certificate
constants as in its batch. No future merge changes these data. At
successful completion the final partial batch and the count marker
account for every appended index. Replacing the child union
in~\eqref{refine:cover-inclusion} by the union of all batch spheres at
radius $\tau$ is therefore valid.
For clarity, the implementation uses the following finite limits. Input
radii satisfy $0\le R\le 0.126$, and the child cap and width satisfy
\[
0.0001\le\tau\le 0.126,\qquad 0.0001\le w\le .05,\qquad 0\le k\le 10.
\]
Normal constants are checked against
\[
\gamma_P<20,\quad
\mathtt{n0}_P,\varepsilon_P<10^{-6},\quad
\alpha_P<10^{-4},\quad
o_P,\mathtt{hor}_P<10^{-5},\quad
\mu_P<2.
\]
Stored phase entries have magnitude at most $4$, and stored tangent
entries at most $1.01$. These range checks accompany genuine Root
certificates; the checks alone do not prove those certificates. A
prospective Root with more than ten tangent columns is rejected, not
truncated.
The base grid is limited to $50{,}000{,}000$ boxes, subdivision depth at
most six, atlas size at most $2{,}000{,}000$ Roots. One record file
holds at most $2{,}000{,}000$ records and the combined input at most
$8{,}000{,}000$. The grid product is checked before multiplication past
its limit. Within the admitted ranges $L\le 1262$; the mixed-radix
indices fit the stated integer types, and every split uses at most
$2^{10}$ children. The integer and floating-point bounds are recorded
in the arithmetic bounds above.
% =====================================================================
\medskip\noindent\textbf{Excluding two further mutually unbiased bases.}
\label{sec:exclusion}
% =====================================================================
The final stage takes a certified Root $(P,T)$ with $0\le r\le .03$ and
rules out two more mutually unbiased orthonormal bases in the phase ball
of radius $r$. The phase variables of the two extra bases are divided
into closed balls, a permissive graph of necessary edge relations is
retained, and a finite search runs for the required clique pattern.
Every rejection is a strict violation of a necessary condition with a
certified margin; no rejection depends on a center being exact.
Throughout, $e(t)=\exp(2\pi i t)$. A real matrix in $\R^{36}$ is
identified with row-major indexing. The Root constants $\omega,n_0$,
$\iota,\Gamma$ denote $\mathtt{ort},\mathtt{n0},\mathtt{id},\mathtt{C}$
as in the arithmetic bounds above. Its certified data give
\begin{equation}\label{exclude:root-data}
\|T^*T-I\|\le\omega,\qquad
\|F(P)\|\le f_0<10^{-6},\qquad 0\le\omega<10^{-5},
\end{equation}
and, for every horizontal $\delta$ with $\|\delta\|\le r$ such that
$e(P+\delta)$ is Hadamard,
\begin{equation}\label{exclude:root-decomposition}
\eta=T^*\delta,\quad \nu=\delta-T\eta,\qquad
\|\eta\|\le r_H:=\sqrt{1+\omega}\,r,\qquad
\|\nu\|\le n:=n_0+\iota r+\Gamma r^2.
\end{equation}
The first two bases are the coordinate basis and the columns of
$e(P+\delta)/\sqrt6$. If further bases exist, they share this single
Hadamard displacement $\delta$, hence this single $\eta$. Each vector
unbiased to the coordinate basis is $e(\beta)/\sqrt6$ for some real
phase vector $\beta$; choose its global phase so that $\beta_0=0$, and
lift $\beta_1,\ldots,\beta_5$ to $[0,1]$. Different $U(1)$ rephasings
of different extra vectors do not change $\delta$.
The single-vector sums use a uniform operator bound. The identity
$\|AA^*-6I\|_F=2\pi\sqrt2\,\|F(P)\|$ and~\eqref{exclude:root-data} give
\begin{equation}\label{exclude:operator-bound}
\|A/\sqrt6\|_{\mathrm{op}}
\le\sqrt{1+\frac{2\pi\sqrt2}6f_0}<K,\qquad K:=1.00001,
\end{equation}
where $A=e(P)$.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Closed phase bins.}
% ---------------------------------------------------------------
At resolution $M$ the five variable phases are covered by closed
intervals $[b/M,(b+1)/M]$, $0\le b<M$. Phase $1$ equals $0$ on the
torus, so this covers every boundary. Let $u$ be a bin midpoint, with
$u_0=0$, and let
\[
a=\frac 1{2M},\ \ d_0=0,\ |d_j|\le a\ (j>0),\ \
x=d-\frac 16\Bigl(\sum_j d_j\Bigr){\bf 1}.
\]
Then the bin vector equals $e(u+x)/\sqrt6$ up to global phase.
Let
$S_i=\sum_je(P_{ij}-P_{i0})$, the (unscaled, global-phase-free)
$i$th entry of the row over inner products. Set
\begin{equation}\label{exclude:single-sums}
c_i=\frac{|S_i|^2-1}{4\pi},\qquad
J_{ij}=-\operatorname{Im}\bigl(\overline{S_i}q_{ji}\bigr),\qquad
B_{ik}=\sum_j J_{ij}T_{(j,i),k},
\end{equation}
where $q_{ji}=e(P_{ij}-P_{i0})$ (so $\sum_jq_{ji}=S_i$; and note
$\sum_jJ_{ij}=0$). For $y\in\R^6$ put
\[
G_y:=\max_j\Bigl|\sum_i y_i\overline{S_i}q_{ji}\Bigr|,
\qquad
N_y:=\Bigl(\sum_{i,j}y_i^2J_{ij}^2\Bigr)^{1/2}.
\]
\begin{lemma}[Direct displacement tests]
\label{exclude:direct-single}
If $\|x\|\le X$, $\|x\|_1\le L$, and $\|\delta\|\le r$, then
\begin{equation}
\|S'-S\|\le 2\pi(KX+r),
\qquad
|S'_i-S_i|\le 2\pi(L/\sqrt6+r).
\label{exclude:sum-change}
\end{equation}
Consequently unbiasedness requires
\begin{equation}\label{exclude:direct-tests}
\max\{0,1-\tau\}^2\le |S_i|^2\le(1+\tau)^2,
\qquad
\sum_i(|S_i|-1)^2\le[2\pi(KX+r)]^2,
\end{equation}
with $\tau=2\pi(L/\sqrt6+r)$.
\end{lemma}
\begin{proof}
First perturb the vector phases by $x$ and then the Hadamard phases by
$\delta$. The first change has norm at most
$K\|(e(x_j)-1)_j\|\le 2\pi KX$. For the second, the $i$th component has
magnitude at most $(2\pi/\sqrt6)\sum_j|\delta_{ji}|
\le 2\pi\|\delta_{:,i}\|$, and their combined norm is at most $2\pi r$.
The scalar bound is the same with the first change at most
$2\pi L/\sqrt6$. Apply $\bigl||S_i|-|S'_i|\bigr|\le|S_i-S'_i|$ and use
$|S'_i|=1$.
\end{proof}
\begin{lemma}[Weighted Taylor estimate]
\label{exclude:weighted-single}
Suppose $\|x\|\le X$ and $\|x\|_\infty\le I$. For every $y\in\R^6$,
unbiasedness implies
\begin{equation}\label{exclude:weighted-residual}
\bigl|y\cdot(c+Jx-B\eta)\bigr|\le\mathcal E(y;X,I),
\end{equation}
where
\begin{align}
\mathcal E(y;X,I):={}&nN_y+\pi\Bigl[
\|y\|_\infty(KX+r)^2
+\min\{G_yX^2,\ K\|S\circ y\|IX\}\notag\\
&\hspace{28mm}
+\frac 2{\sqrt6}\|S\circ y\|Xr
+\frac 1{\sqrt6}\|S\circ y\|_\infty r^2\Bigr].
\label{exclude:weighted-error}
\end{align}
For a bin, this gives the necessary support inequality
\begin{equation}\label{exclude:bin-support}
|y\cdot c|
\le a\sum_{j>0}|(J^*y)_j|+r_H\|B^*y\|
+\mathcal E(y;X_M,I_M).
\end{equation}
\end{lemma}
\begin{proof}
Write $\Delta_i=S'_i-S_i$. The expansion of
$(|S'_i|^2-1)/(4\pi)$ has constant term $c_i$, linear term
$\sum_jJ_{ij}(x_j-\delta_{ji})$, and the squared increment
$|\Delta_i|^2/(4\pi)$. After multiplication by $y_i$ and summation the
latter is bounded by $\pi\|y\|_\infty(KX+r)^2$ via
\eqref{exclude:sum-change}.
The pure-$x$ linear remainder was bounded in the preceding scalar
Taylor expansion; its two alternatives are
$G_yX^2$ and $K\|S\circ y\|IX$, and the minimum is valid. The
Hadamard-$\delta$ linear term is $B\eta$; its remainder with $\nu$ is
bounded by $nN_y$. The cross term and pure $\delta$ term give the
remaining $\|S\circ y\|Xr$ and $\|S\circ y\|_\infty r^2$ contributions.
Subtracting the bin's linear box term yields the necessary support
bound~\eqref{exclude:bin-support}.
\end{proof}
The bin routine \texttt{Layer::test} applies~\eqref{exclude:direct-tests},
then the scalar test~\eqref{exclude:scalar-linear}
(see below), and finally~\eqref{exclude:bin-support} for proposed
directions. Proposed directions come from an inexact linear solve; the
proposal method is immaterial because the final inequality holds for
every real $y$ and is checked at the stored test direction.
There is a useful scalar test before the tangent decomposition. Set
$S_{\max}=\max_i|S_i|$, and for a bin write $X=X_M$, $L=L_M$, $I=I_M$.
Unbiasedness requires
\begin{equation}\label{exclude:scalar-linear}
|c_i|\le a\sum_{j>0}|J_{ij}|+r\|J_{i,:}\|+\rho,
\end{equation}
where
\begin{align}
\rho=\min\Bigl\{&
\pi\bigl[(KX+r)^2+
S_{\max}\bigl(KIX+(2Xr+r^2)/\sqrt6\bigr)\bigr],\notag\\
&\pi\bigl[(L/\sqrt6+r)^2+
S_{\max}(X+r)^2/\sqrt6\bigr]\Bigr\}.
\label{exclude:scalar-remainder}
\end{align}
The first remainder bound follows from the weighted estimate with a
coordinate direction and without decomposing $\delta$. The second uses
the scalar increment bound in~\eqref{exclude:sum-change} and applies the
second-order exponential bound directly to $x-\delta_{:,i}$, whose norm
is at most $X+r$. The full Hadamard linear term is supported by
$r\|J_{i,:}\|$, so no normal-error term is added.
The implementation checks that a proposed direction is finite with
maximum absolute entry strictly between $10^{-100}$ and $10^{100}$,
normalizes by that maximum, and uses the resulting stored coefficients
as the test direction. The stored direction itself is what matters.
% ---------------------------------------------------------------
\medskip\noindent\textbf{From bins to closed vector balls.}
% ---------------------------------------------------------------
The initial resolution is $M=12$. Every bin is tested, and each retained
bin is split into all $2^5=32$ children at each
doubling $12\to 24\to 48\to 96$. On one coordinate, the child indices
$2b$ and $2b+1$ cover the parent interval exactly, so every vector
follows at least one surviving bin path, including on a phase or torus
boundary. Exceeding the cap of $2{,}000{,}000$ retained bins is an
unresolved outcome, never an exclusion. In the stated invocation there
are four bin layers.
For $v\in\R^6$ and $s\ge 0$, define the closed vector ball
\begin{equation}\label{exclude:vector-ball}
\mathcal B_{\mathrm{vec}}(v,s)=
\Bigl\{\frac{\lambda e(v+x)}{\sqrt6}:
|\lambda|=1,\ x\in\R^6,\ \sum_jx_j=0,\ \|x\|\le s\Bigr\}.
\end{equation}
\begin{lemma}[Certified clustering]
\label{exclude:cluster-cover}
Suppose each retained bin midpoint $u$ is assigned to a center $v$ and
an integer vector $m\in\mathbb Z^6$. Let
\[
h=u-v-m-\frac16\Bigl(\sum_j(u_j-v_j-m_j)\Bigr){\bf 1}.
\]
If the assigned radius is at least $\|h\|+X_M$, then
$\mathcal B_{\mathrm{vec}}(v,s)$ covers the whole bin up to global
phase. The final assignment in \texttt{cluster} has this property.
\end{lemma}
\begin{proof}
For any raw displacement $d$ in the bin, subtracting its mean gives the
$x$ of the direct displacement representation. With the assigned
integers fixed, the phase difference from $v$ is $h+x$ up to an integer
vector and a scalar phase. It has zero mean and norm at most
$\|h\|+X_M$, and every bin boundary is attained by some representative
of the closed bin, so the covering assertion also holds on every
boundary.
\end{proof}
Each assignment iteration computes coordinatewise integer wraps of
$u-v$, then subtracts their mean. A candidate rejected by its distance
cutoff is simply not used; the bin is assigned to another center or to
a newly inserted midpoint. Thus no bin is dropped by the heuristic. For
each assigned center the radius is updated to the maximum of the
computed projected midpoint distance plus the bin-norm upper bound and
$10^{-9}$; this is an outward upper bound, as certified in the
arithmetic bounds above. New midpoint centers use midpoint distance
zero, with the same allowance.
The earlier Lloyd iterations adapt the nearest-center/centroid heuristic
of Lloyd's 1982 least-squares quantization paper. They may move centers,
but a subsequent complete assignment follows each such movement. In the
last iteration the assigned centers are not moved. Deleting centers with
no assigned bins cannot lose coverage. The returned centers and radii
therefore meet the covering assertion. An attempt to exceed $2048$
centers, or a failure of the final checks $|v_j|\le 2$ and $0\le s<.1$,
is unresolved rather than negative.
For the specified invocation, the heuristic radius parameter is $.042$
and there are two assignment iterations. With $w=1/M$ the candidate
threshold is $(.042-1.11w)^2$, and the coordinate cutoff is the square
root of twice that threshold plus $10^{-10}$. These constants affect
only the choice of assignments; the proof is the final radius update in
Lemma~\ref{exclude:cluster-cover}. In particular, the argument does not
assume that the selected wrapped lift minimizes torus distance.
\begin{lemma}[Distinct labels]
\label{exclude:distinct-labels}
Two vectors in the same ball $\mathcal B_{\mathrm{vec}}(v,s)$ with
$s<.1$ are neither orthogonal nor mutually unbiased. Hence the twelve
vectors in two further mutually unbiased orthonormal bases must occupy
twelve distinct retained balls.
\end{lemma}
\begin{proof}
After removing global phases, their unscaled inner-product sum has the
form $\sum_je(z_j)$ with $\|z\|\le 2s<.2$. Therefore
\[
\Bigl|\sum_je(z_j)\Bigr|
\ge 6-\sum_j|e(z_j)-1|
\ge 6-2\pi\sqrt6\,\|z\|
>6-\frac25\pi\sqrt6>\frac{20}7>\sqrt6.
\]
Orthogonality requires the unscaled magnitude to be zero, and
unbiasedness requires it to be $\sqrt6$. Every pair among the twelve
has one of these two relations, so no pair shares its ball.
\end{proof}
% ---------------------------------------------------------------
\medskip\noindent\textbf{Bin geometry and edge tests.}
% ---------------------------------------------------------------
\begin{lemma}[Bin geometry]
\label{exclude:bin-geometry}
Fix a bin with midpoint $u\in\R^6$ ($u_0=0$) and halfwidth
$a=1/(2M)$. Any phase vector in the bin is $u+d$ with $d_0=0$ and
$|d_j|\le a$ for $j\ge1$. The horizontal displacement
$x=d-\frac16(\sum_jd_j){\bf 1}$ satisfies $\sum_jx_j=0$ and, with the
constants $X_M,L_M,I_M$ of ,
\begin{equation}\label{exclude:bin-norms}
\|x\|\le X_M,\quad \|x\|_1\le L_M,\quad \|x\|_\infty\le I_M.
\end{equation}
\end{lemma}
\begin{proof}
$\sum_jx_j=0$ is immediate from the definition.
Write $d=\bar d+(d-\bar d)$ with $\bar d=(1/6)\sum_jd_j$. Then
$\|x\|^2=\|d-\bar d{\bf 1}\|^2=\|d\|^2-6\bar d^2\le\|d\|^2\le5a^2$,
so $X_M\ge\sqrt5\,a$ suffices. Each $|d_j-\bar d|\le 2a$, so
$L_M\ge10a$ and $I_M\ge2a$ also suffice.
accounts for the $10^{-12}$ margin.
\end{proof}
Take two distinct balls with centers $u,v$ and radii $s_u,s_v<.1$. For
actual representatives with displacements $x_u,x_v$, set
\begin{equation}\label{exclude:edge-data}
Q_j=e(u_j-v_j),\quad S=\sum_jQ_j,\quad
z=x_u-x_v,\quad r_e=s_u+s_v,\quad \|z\|\le r_e.
\end{equation}
The actual unscaled inner product is, up to conjugation and global
phase, $S'=\sum_jQ_je(z_j)$. Since $|S'-S|\le 2\pi\sqrt6\,r_e$,
necessary direct tests for orthogonality ($\mathcal O$) and
unbiasedness ($\mathcal U$) are
\begin{equation}\label{exclude:edge-distance}
\frac{|S|}{2\pi\sqrt6}\le r_e\ (\mathcal O),\qquad
\frac{\bigl||S|/\sqrt6-1\bigr|}{2\pi}\le r_e\ (\mathcal U).
\end{equation}
\begin{equation}\label{exclude:edge-derivatives}
(D_{\mathcal O})_{:,j}=
\begin{pmatrix}-\operatorname{Im}Q_j\\ \operatorname{Re}Q_j\end{pmatrix},
\qquad
(D_{\mathcal U})_j=
\frac{\operatorname{Re}S\,(D_{\mathcal O})_{1j}
+\operatorname{Im}S\,(D_{\mathcal O})_{2j}}6.
\end{equation}
\begin{lemma}[Edge Taylor estimates]
\label{exclude:edge-taylor}
For any actual pair from two distinct closed vector balls, with the
data of~\eqref{exclude:edge-data}, and for
$S_{\mathcal O}=\operatorname{Re}S$,
$S_{\mathcal U}=(|S|^2-6)/(4\pi)$,
\begin{align}
|\operatorname{Re}(S'-S)-(D_{\mathcal O}z)_1|
&\le E_{\mathcal O}:=4\pi\sqrt6\,r_e^2, \label{exclude:orth-residual}\\
\Bigl|\frac{|S'|^2-6}{4\pi}-S_{\mathcal U}-D_{\mathcal U}z\Bigr|
&\le E_{\mathcal U}:=
\frac{2\pi\max(|S|,\sqrt6)}{\sqrt6}\,r_e^2 +4\pi r_e^2.
\label{exclude:unbiased-residual}
\end{align}
Therefore the necessary support conditions
\begin{align}
|S_{\mathcal O}|
&\le r_e\|D_{\mathcal O}\|_F+E_{\mathcal O},
\label{exclude:orth-support}\\
|S_{\mathcal U}|
&\le r_e\|D_{\mathcal U}\|+E_{\mathcal U}
\label{exclude:unbiased-support}
\end{align}
hold for every actual pair.
\end{lemma}
\begin{proof}
Write $e(z_j)-1=2\pi iz_j+\varepsilon_j$ with
$|\varepsilon_j|\le2\pi^2z_j^2$ by the integral Taylor remainder.
Multiply by $Q_j=a_j+ib_j$ and separate real and imaginary parts:
$\operatorname{Re}(Q_j(e(z_j)-1))\approx-2\pi b_j z_j=(D_{\mathcal O}z)_j$
and $\operatorname{Im}(Q_j(e(z_j)-1))\approx2\pi a_jz_j$ up to a
residual of modulus at most $4\pi^2 r_e^2$ per coordinate. Summing over
$j$ and using $\|z\|\le r_e$ proves~\eqref{exclude:orth-residual}.
For unbiasedness, expand $(|S'|^2-6)/(4\pi)$ at $S$. The constant is
$S_{\mathcal U}$, the linear term is $D_{\mathcal U}z$ (the division by
$4\pi$ normalizes the derivative exactly), and the quadratic remainder
is at most $(2\pi|S|/\sqrt6)r_e^2+4\pi^2r_e^2$; substituting the
residual of the orthogonality calculation
gives~\eqref{exclude:unbiased-residual}. Cauchy--Schwarz on
$\|z\|\le r_e$ gives the support inequalities.
\end{proof}
% ---------------------------------------------------------------
\medskip\noindent\textbf{The joint test with one shared tangent variable.}
% ---------------------------------------------------------------
For each ball $t=0,1$, form $c^t,J^t,B^t,S^t$ from
\eqref{exclude:single-sums} at its actual center. Write its radius as $s_t$ and put
\begin{equation}\label{exclude:ball-error}
\mathcal E_t^{\mathrm{ball}}(y):=
n\Bigl(\sum_{i,j}y_i^2(J^t_{ij})^2\Bigr)^{1/2}
+\pi\Bigl[
\|y\|_\infty(Ks_t+r)^2
+s_t^2G_y^t
+\frac 2{\sqrt6}\|S^t\circ y\|s_tr
+\frac 1{\sqrt6}\|S^t\circ y\|_\infty r^2\Bigr].
\end{equation}
This is exactly the weighted estimate~\eqref{exclude:weighted-error}
restricted to its first pure-$x$ remainder (do not minimize against the
second one here), with $X=s_t$ and $I=s_t$; they are valid on the whole
closed Euclidean vector ball because the bound was derived for the
worst case over that ball.
For an edge type $E\in\{\mathcal O,\mathcal U\}$ write $\ell=2$ or
$1$ for the dimension of $D_E$, and let $q_E\in\R^\ell$ be the
constant part of the linearized edge equation:
\begin{equation}\label{exclude:edge-data-2}
q_{\mathcal O} := (\operatorname{Re}S, \operatorname{Im}S) \in\R^2,
\qquad
q_{\mathcal U} := \frac{|S|^2-6}{4\pi} \in\R^1,
\end{equation}
with $D_E$ and $E_E$ as in~\eqref{exclude:edge-derivatives}
and~\eqref{exclude:orth-residual}--\eqref{exclude:unbiased-residual}.
\begin{lemma}[Joint necessary support]
\label{exclude:joint-support}
Suppose one vector from each of the two balls is unbiased to the
coordinate basis and to the \emph{same} second basis, whose phase matrix
is $e(P+\delta)/\sqrt6$, and the two vectors have relation
$E\in\{\mathcal O,\mathcal U\}$. For every $y_0,y_1\in\R^6$ and
$y_e\in\R^\ell$,
\begin{align}
&\bigl|y_0\cdot c^0+y_1\cdot c^1+y_e\cdot q_E\bigr|\notag\\
&\quad\le
\sum_{t=0}^1\Bigl[
s_t\bigl\|(J^t)^*y_t+(-1)^tD_E^*y_e\bigr\|
+\mathcal E_t^{\mathrm{ball}}(y_t)\Bigr]\notag\\
&\qquad+r_H\bigl\|(B^0)^*y_0+(B^1)^*y_1\bigr\|
+E_E\|y_e\|.
\label{exclude:joint-inequality}
\end{align}
\end{lemma}
\begin{proof}
For each single vector the linearized equation is
$c^t+J^tx_t-B^t\eta$ with weighted residual bounded by
$\mathcal E_t^{\mathrm{ball}}(y_t)$. The edge equation is
$q_E+D_E(x_0-x_1)$ with residual norm at most $E_E$.
Divide the edge equation by $y_e$, expand each $x_t$-dependent part
for the two balls separately, and use the $x$-support of each single
vector's displacement: for ball $t=0$ the combined coefficient of
$x_0$ is $(J^0)^*y_0+D_E^*y_e$; for $t=1$, $(J^1)^*y_1-D_E^*y_e$. These
support over $\|x_t\|\le s_t$ give the first two norm terms.
The single-vector coefficients of the shared tangent $\eta$ are $-B^0$
and $-B^1$, and combine to give the coefficient
$-(B^0)^*y_0-(B^1)^*y_1$. Support this over $\|\eta\|\le r_H$. Bound
the three residuals by the triangle inequality and Cauchy's inequality.
\end{proof}
This is exactly the comparison implemented by \texttt{graph::joint}.
It first calls \texttt{single::support} with a zero vector-displacement
coefficient array for each ball, and only then adds the edge
derivatives. Thus the maximum coefficient in~\eqref{exclude:ball-error}
is formed from the pure single-vector remainder; it cannot be decreased
by cancellation with an edge coefficient. The tangent coefficient
array, in contrast, is shared and accumulates both $-(B^t)^*y_t$
contributions before its norm is taken.
As with the single-vector test, the proposed direction comes from an
inexact linear solve; only its final inequality must be certified. The
final inequality holds for every real $y$; the numeric choice of
direction is not part of the certificate.
The inequalities just proved are necessary in exact real arithmetic. In
the implementation, each rejection is a strict violation of an
outward-enlarged bound; equalities and unresolved comparisons retain
the candidate. certifies these comparisons,
and lists their margins. Thus every
genuine pair retains its required relation when the tests are evaluated
numerically.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Exact graph exclusion.}
% ---------------------------------------------------------------
Build one vertex for each retained ball. Join two distinct vertices by
an $\mathcal O$ edge iff their orthogonality type has not been excluded
by the direct and joint tests above; similarly for $\mathcal U$ edges
and unbiasedness. Both graphs are undirected and have no self-loops.
Every genuine pair retains its corresponding edge by the preceding
lemmas. Extra edges are harmless: retained pairwise possibilities need
not have a simultaneous realization.
By Lemma~\ref{exclude:cluster-cover} and
Lemma~\ref{exclude:distinct-labels}, two additional mutually unbiased
orthonormal bases would give disjoint vertex sets $A,B$ satisfying
\begin{equation}
|A|=|B|=6,\qquad A,B\text{ are }\mathcal O\text{-cliques},\qquad
A\times B\subseteq E(\mathcal U).
\label{exclude:graph-configuration}
\end{equation}
The following finite search tests this necessary configuration.
\begin{lemma}[Complete graph search]
\label{exclude:graph-search}
The mode-$1$ search of \texttt{graph::exist}, started with both pools
equal to the full vertex set and first-clique target size six, returns
false only if~\eqref{exclude:graph-configuration} is impossible. A
work-cap exhaustion never gives a negative certificate.
\end{lemma}
\begin{proof}
First consider \texttt{clique} on a candidate set $P$ and target size
$m$. For $m=0$ it succeeds, and for $m=1$ it tests nonemptiness. For
$m\ge2$ it repeatedly selects the least remaining vertex $i$, removes
it from $P$, and searches $P\cap N_{\mathcal O}(i)$ for an
$(m-1)$-clique. For $m=2$ that search is just nonemptiness. Every
clique has a unique least vertex, so this enumerates every possibility.
Terminating when fewer than $m$ candidates remain is also exact.
Now suppose $A,B$ satisfy~\eqref{exclude:graph-configuration}.
Interchange their names if necessary so that
$a=\min(A\cup B)$ belongs to $A$. At the root, \texttt{exist} tries
possible least vertices of the first clique in increasing order. It
deletes each tried root vertex from both its local first pool and its
second pool. Thus when $a$ is selected, every required remaining vertex
in $A\cup B$ is still available. This root deletion imposes only the
chosen convention that the minimum of the union lies in the first
clique, and it is performed only when the remaining first-clique target
is six. Recursion then runs \texttt{clique} on the first-pool target
six, and for each resulting candidate $A$, on the second-pool target
six. When \texttt{clique} succeeds, its exact enumeration is used to
test $A\times A'\subseteq E(\mathcal U)$ by checking every cross pair
in order.
The search cap and any other finite work bound turn into failure, not
into a claim that the pattern is impossible. Thus a work-cap exhaustion
is not a negative certificate; only a completely enumerated negative
result is.
\end{proof}
\begin{proposition}[Local exclusion]
\label{exclude:local-certificate}
Suppose
\begin{enumerate}
\item the retained bins cover the closed phase ball of the Root in the
sense of Lemma~\ref{exclude:cluster-cover};
\item every retained bin has either been excluded by
\eqref{exclude:direct-tests}, \eqref{exclude:scalar-linear}, or
\eqref{exclude:bin-support}, or been accepted and covered by the
assigned ball;
\item the $\mathcal O$- and $\mathcal U$-graphs on the retained balls,
after excluding any pair by the direct or joint tests conditional
on $\eta$, admit no
configuration~\eqref{exclude:graph-configuration} as decided by a
completed \texttt{graph::exist} search.
\end{enumerate}
Then the phase ball contains no two additional mutually unbiased
orthonormal bases.
\end{proposition}
\begin{proof}
Suppose two such bases existed. Each of their twelve vectors is
covered by some retained ball (Lemma~\ref{exclude:cluster-cover}).
Lemma~\ref{exclude:distinct-labels} puts them in twelve distinct balls.
Within each basis the six vectors are pairwise orthogonal, and every
intra-pair survived its $\mathcal O$ direct and joint tests because the
tests are necessary conditions; thus they are joined by $\mathcal O$
edges. Across the two bases every pair is unbiased, and similarly
they survive the $\mathcal U$ tests, so the cross-pairs are $\mathcal
U$ edges. This produces a configuration~\eqref{exclude:graph-configuration},
contradicting the third hypothesis.
\end{proof}
The conclusion of Proposition~\ref{exclude:local-certificate} is
conditional on a completed negative result for the individual Record;
no collection of unexecuted tests is itself such a result.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Implementation margins and invocation.}
\label{exclude:implementation}
% ---------------------------------------------------------------
For a bin the source forms its halfwidth with a $10^{-12}$ increase,
then increases each of $X_M,L_M,I_M$ by a further $10^{-12}$ after the
indicated multiplication. The stored tangent and normal bounds are
\begin{equation*}
\widehat r_H=\sqrt{1+\omega}\,r+10^{-9},\qquad
\widehat n=n_0+\iota r+\Gamma r^2+10^{-9}.
\end{equation*}
The upper squared-magnitude threshold in~\eqref{exclude:direct-tests} is
increased by $2\cdot10^{-7}$, and its lower threshold, when $\tau<1$, is
decreased by that amount. The combined squared-distance threshold, the
right side of~\eqref{exclude:scalar-linear}, and the right side
of~\eqref{exclude:bin-support} each receive a $2\cdot 10^{-7}$
comparison allowance.
For a pair, the radius sum $s_u+s_v$ is increased by $10^{-10}$. The
distance tests~\eqref{exclude:edge-distance} permit a further
$2\cdot 10^{-7}$ on that radius. Both edge residual bounds
$E_{\mathcal O}$ and $E_{\mathcal U}$ are increased by $2\cdot 10^{-8}$;
the direct comparisons~\eqref{exclude:orth-support}
and~\eqref{exclude:unbiased-support} receive another $2\cdot 10^{-7}$.
The right side of the joint
comparison~\eqref{exclude:joint-inequality} receives $3\cdot 10^{-7}$.
Cluster assignment radii include the $10^{-9}$ allowance described in
Lemma~\ref{exclude:cluster-cover}.
The bounds $|P_{ij}|\le 4$ and $|v_j|\le 2$ keep all trigonometric
arguments used here within the guarded range of \texttt{expair}. The
proof of includes the errors in these
arguments, in $q,S,c,J,B$, in the short weighted sums and norms, and in
the positive support expressions. Inexact solution of the proposal
systems supplies no separate assumption.
For clarity, the invocation parameters used for this certificate are
\texttt{12 4 1 .042 1 1 5000000}: initial resolution $12$, four bin
layers, one proposed bin direction, clustering parameter $.042$, one
proposed joint direction, graph mode $1$, and the stated graph work
cap. The two Lloyd iterations are the fixed default. The conclusion of
Proposition~\ref{exclude:local-certificate} is conditional on a
completed negative result for the individual Record; no collection of
unexecuted tests is itself such a result.
% =====================================================================
\medskip\noindent\textbf{Finite execution and the recorded certificate.}
\label{sec:execution}
% =====================================================================
The computation used in the conclusion of this proof was completed on
4 October 2026 (UTC). This section specifies its finite inputs, the
tests for complete execution, and the implication from those tests to
nonexistence. The geometric and arithmetic justification of each
rejection lies in the seed cover, the completion and refinement stages,
the exclusion tests, and the certified arithmetic. A successful exit by
itself is not the
certificate: the inputs, record counts, saved failures, and coverage of
all shards must agree as described here.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Programs, records, and machine arithmetic.}
\label{execute:programs}
% ---------------------------------------------------------------
The recorded platform was x86-64, Intel Xeon Platinum 8573C, Ubuntu
GCC 13.3.0, and glibc 2.39. Each stage was compiled separately, all four on the same UTC date as
the run, with contraction disabled and a fixed language dialect. The
arithmetic assumptions are binary64, round to nearest, separate
operations without contraction or fast-math reassociation, correctly
rounded binary64 square root, masked floating-point exceptions, and
neither flush-to-zero nor denormals-are-zero. Separate probes compiled
with the cover and A/B/C flag sets both passed and recorded
\texttt{FLT\_EVAL\_METHOD=0}, 32-bit \texttt{int}, 64-bit
\texttt{size\_t}, 64-bit \texttt{uint64\_t} equal to
\texttt{unsigned long}, eight-bit bytes, and MXCSR \texttt{0x1f80}.
Both probes recorded nearest rounding on entry and during the check,
and passed the square-root and gradual-underflow samples. These probes
and the program checks provide limited evidence for the stated platform
conditions, not a proof of all primitive-operation or compiler
semantics. Those semantics remain assumptions; the certified arithmetic
supplies the numerical bounds under them.
Four programs carry out the search: one builds the seed cover, one runs
the completion stage, one refines the resulting phase balls, and one
performs the exclusion tests. They share the Root, Atlas and graph data
structures through common declarations, and the arithmetic bounds below
fix every floating-point operation any of them performs. The archived
sources, identified by their recorded digests, are what make the
certificate concrete rather than notional.
\medskip\noindent\textbf{How the pipeline is split.}
The first program writes a seed-center table, one entry per retained
midpoint, carrying its rational phases in cycles, the selected column
position $s$, and the radius $0.032$. It performs only the seed filters
whose margins are certified in the arithmetic bounds.
The completion stage processes seeds in disjoint shards and fills the
atlas of certified phase balls. The refinement stage then shrinks that
cover and writes its batch files, and the exclusion stage tests every
resulting ball for two extra bases, writing a per-batch record of the
few that resist exclusion.
Each program reports its own exit status and counts. A supervising
driver records one outer summary and one summary per stage, reruns the
fallback and the final exclusion shards whenever the initial exclusion
pass leaves any record unexcluded, and compares every intermediate
count against the outer summary before the run counts as complete.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Control flow and terminal conditions.}
\label{execute:control}
% ---------------------------------------------------------------
The control flow is finite and every terminal condition is documented.
A run that fails any listed condition is not a certificate.
The required terminal conditions are:
\begin{enumerate}
\item successful builds, a normally completed cover, and a complete
center file;
\item four A exits $0$, all assigned centers processed, terminal
\texttt{failedBalls=0}, and four empty bad text files;
\item both B exits $0$, all input parents processed, agreement of
logged and saved output counts, and zero-count B bad files;
\item every initial C batch processed exactly once with a consistent
exit, done count, bad count, and saved record count;
\item every initial C bad record is included in fallback, fallback
exits $0$ with all parents processed and zero saved failures,
and the four final C exits are $0$ with total done equal to the
fallback count and zero saved failures.
\end{enumerate}
Missing files, truncated data, an unexpected exit, resource exhaustion,
an unresolved A or B parent, or a final C bad record prevent
acceptance. Recoverable search diagnostics are distinguished from these
terminal conditions. For example, A's \texttt{ownfail} counts
unsuccessful proposals of a local Root cover that are followed by
subdivision. A0 reported such a diagnostic in this run; its terminal
unresolved count was nevertheless zero.
The finite bounds governing this control flow --- the bin counts, the
range of stored differences, the grid~$64$ cell-test count, and the
extra-node and deepest-grid caps --- are stated and checked in the
completion stage.
Reaching the C-stage bin cap of $2{,}000{,}000$ retained bins, exceeding
the cluster center count of $2048$, or exhausting the graph work cap
each produce an unresolved outcome rather than a negative certificate.
% ---------------------------------------------------------------
\medskip\noindent\textbf{The recorded execution.}
\label{execute:observed}
% ---------------------------------------------------------------
The complete run consists of 100 program invocations:
\begin{itemize}
\item one pass over the seed cover;
\item four shards of the completion stage;
\item two refinement passes, a mid one and a fine one;
\item 88 batches of the initial exclusion pass;
\item one fallback refinement pass;
\item four final exclusion batches.
\end{itemize}
The initial exclusion pass leaves a nonzero count on some batches: the
residual candidates that the bin and graph tests could not reject. They
are handed to the fallback, whose own batches then feed the four final
exclusion batches. Recorded counts and hashes accompany the whole run.
For every initial exclusion process the recorded status is $0$ exactly
when its unresolved count is zero, and $1$ exactly when that count is
positive. The four failure categories --- uncapped graph candidate,
graph-step cap, bin overflow, and clustering failure --- sum to that
count, and none of them is an exclusion. Every nonempty batch reaches
the fallback, sorted by locale, including batches whose locations differ
from one run to the next. The final exclusion pass uses all four
modulo-$4$ parts of the fallback output. Had the initial pass left no
unresolved record, the fallback and the final pass would have been
unnecessary; that shortcut was not taken.
Fifteen initial exclusion batches left unresolved records, and together
they carry all $163$ unresolved records. There were no initial bin
overflows, clustering failures, or uncapped graph candidates. There were
no final failures of any category. All 100 program invocations are
therefore accounted for: one cover pass, four completion shards, two
refinement passes, 88 initial exclusion batches, one fallback pass, and
four final exclusion batches. The fifteen nonzero exits are fully
resolved by that fallback chain. No rejection threshold, graph cap,
refinement width, or mathematical input was altered, and no additional
refinement level was inserted.
For scale, the measured wall times, including polling and validation,
were $621.59$ seconds for the cover, $2020.32$ seconds for the four
parallel completion shards, $1071.49$ and $407.45$ seconds for the two
refinement passes, and $300.23$, $8.12$ and $36.27$ seconds for the
initial exclusion, the fallback and the final exclusion. The outer
supervisor measured $4470.33$ seconds for the complete pipeline. These
are observations on this platform, not native CPU times or running-time
guarantees. No stage hit its wall deadline, and no retry or resource
limit was raised.
The execution archive preserves the numerical data and records of this
complete run, with internal path and host fields redacted: sources and
binaries, center and record files, all 88 initial exclusion batches and
their unresolved records, fallback data, every final unresolved record,
stdout and stderr, commands, exits, platform information, and validation
records. Inside it one record describes the whole run, one record per
stage describes that stage, and one record per child describes that
child. Rebuild and execution instructions accompany the archive.
Archived native binaries document this run; they are not a portability
claim.
The certificate depends only on this newly completed chain and the
counts recorded alongside it. No earlier execution, diagnostic sample,
or intermediate timed measurement substitutes for any stage of this
full run.
% ---------------------------------------------------------------
\medskip\noindent\textbf{The finite-execution certificate.}
\label{execute:theorem}
% ---------------------------------------------------------------
\begin{theorem}[Finite-execution certificate]
\label{execute:finite-theorem}
Fix the four recorded programs and the arithmetic conditions of the
certified arithmetic. Suppose the recorded invocations meet every
terminal condition listed above, together with the record provenance
and batch checks stated there. Then four mutually unbiased orthonormal
bases in $\C^6$ do not exist.
\end{theorem}
\begin{proof}
Suppose four such bases exist. Normalize the first to the standard
basis and write the second as $H/\sqrt6$.
\emph{Seeds.} Proposition~\ref{seed:cover} places the seed coordinates
of an equivalent exact Hadamard matrix in one of the generated seed
balls. All A shards together process that seed.
\emph{Stages A and B.} Proposition~\ref{phase:A-cover} therefore places
an equivalent exact solution in a coarse Record. Each equivalence
preserves extendibility by two further bases
(Lemma~\ref{seed:equivalences}). On that Record the refinement stage
replaces the ball by a finite atlas of children together with any
unresolved parent; Proposition~\ref{refine:cover-preservation}
gives containment. Since the run has no remaining bad parent at the
end of the initial C stage and no bad record at the end of the final C
stage, the record is either covered by a completed child at its own
precomputed local certificate, or refined into a further record that is
similarly covered.
\emph{Stage C.} Proposition~\ref{exclude:local-certificate} then rules
out two extra bases in every such local ball. The exclusion run leaves
no unresolved record at the end, by terminal condition~5 together with
the terminal conditions of the fallback chain. Therefore every
individual pointer to the Hadamard solution lies in a phase ball that
the exclusion stage rules out.
Contradiction: the assumed family would induce a Hadamard displacement
that lands inside some ball of the final cover, but each such ball has a
completed negative certificate. Hence four such bases do not exist.
\end{proof}
The certificate therefore has three independent components: the
finite-arithmetic premises of the certified arithmetic, the finite
invocations and terminal conditions listed above, and the mathematical
chain of the preceding stages.
A change to any one of them --- different build flags, different
terminal conditions, or a different proof of one of the local lemmas ---
would require a fresh execution of the affected stage, but the logical
structure of the certificate is unchanged.
The theorem and its proof complete the computer-assisted upper bound.
% =====================================================================
\medskip\noindent\textbf{Conclusion of the proof.}
\label{sec:conclusion}
% =====================================================================
The proof of Theorem~\ref{thm:main} is now complete. The
overall argument is the following chain of strictly finer covers, each
preserving every exact candidate solution.
\begin{enumerate}
\item By Lemma~\ref{seed:equivalences}, every complex Hadamard matrix,
up to allowed symmetries, has its nine-phase seed inside one of
finitely many closed balls (Proposition~\ref{seed:cover}).
\item For each such seed ball, the possible completions form a
$U(2)$-orbit; the $U(2)$ cells are divided into finitely many
pieces such that each piece either is excluded by a certified
necessary condition or is covered by a small phase ball
(Proposition~\ref{phase:A-cover}).
\item Each Phase ball is subdivided into smaller children, and the
original balls are covered by new children together with any
unresolved bad parents (Proposition~\ref{refine:cover-preservation}).
\item Each small phase ball is ruled out by the C-stage graph test on
independent vector balls and pairwise tests
(Proposition~\ref{exclude:local-certificate}).
\end{enumerate}
The recorded execution of the pipeline
meets every terminal condition listed above, so the
finite-execution certificate Theorem~\ref{execute:finite-theorem}
applies and yields the negative upper bound. The lower bound is
Lemma~\ref{lem:three-mub-six}; the two together give $N(6)=3$.
The argument can be summarized in four sentences: cover the space of
complex Hadamard matrices; replace the cover with smaller and smaller
open covers that still contain every exact solution; within each small
cover, exclude two further mutually unbiased bases by the graph test;
and check the whole chain by an actual computer run with certified
arithmetic.
The historical and mathematical antecedents surveyed above are
complementary; the combination
needed for the theorem above is new. In particular, the present proof
covers every second-basis Hadamard matrix directly, rather than running
a certified search inside a presumed family.
\appendix
% =====================================================================
\medskip\noindent\textbf{Certified arithmetic for the recorded execution.}
\label{sec:rounding}
% =====================================================================
The outward error estimates behind each inequality used above are fixed
below. Every certified rejection is a strict violation
of an outward bound computed under the primitive-operation model below.
Nothing here assumes that a floating-point center is an exact solution,
and no step relies on the observed size of an error over a specific
input.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Primitive-operation model and unit roundoff.}
% ---------------------------------------------------------------
The primitive-operation model is the binary64 specification
of IEEE 754-2019, in its Sections 3.4, 3.6, 4.3, 5.4.1 and 7.5.
The accumulation estimates use the usual product-of-rounding-factors
argument of Section 2.12 and Lemma 3.1 of Higham's 2002 book on the
accuracy of numerical algorithms, with the explicit additive underflow
allowance retained here.
The arithmetic premises are separate, correctly rounded binary64
addition, subtraction, multiplication, division, and square root in
round-to-nearest mode, with binary64 evaluation of intermediate results.
The build uses \texttt{-ffp-contract=off} and does not use fast-math or
any reassociation that changes these operations. Subnormals are not
flushed to zero, and denormal inputs are not treated as zero. Floating
exceptions are masked, so an unsuccessful heuristic calculation can
produce an infinity or NaN and be rejected by its subsequent checks.
The relevant integer types, numeric locale, and file representation are
specified with the execution above.
The function \texttt{environment} sets and checks nearest rounding and
checks the binary64 properties of \texttt{double}. The three entry
points also inspect MXCSR, rejecting a non-nearest mode, DAZ, FTZ, or
unmasked exceptions. The same no-DAZ/no-FTZ environment is required for
the cover-building program, although its entry point does not itself
inspect those bits. The default-nearest premise also applies to
constants and to the global basis constructed before \texttt{main}.
These are explicit arithmetic premises, not a claim that a source-level
static assertion formally verifies every compiler transformation. The
recorded build checks and their limits are described with the execution
above.
The functions \texttt{nearbyint} and \texttt{floor} choose integers at
finite bounded arguments. The particular integer chosen at a boundary is
permissible wherever the proof subsequently uses that integer lift. No
accuracy premise is imposed on \texttt{atan2}: Candidate checks the unit
exponential of the phase it proposes. Decimal literals are either used
as exact stored parameters, or allowed a relative translation error of
at most $5\cdot 10^{-16}$ when compared with a displayed decimal.
Put
\begin{equation}\label{round:unit}
u=2\cdot 10^{-16},\qquad \tau=10^{-300},\qquad
\gamma_m=2.1\cdot 10^{-16}(m+1).
\end{equation}
For a correctly rounded basic operation with finite exact result $z$ at
the scales below,
\begin{equation}\label{round:basic}
|\operatorname{fl}(z)-z|\le u|z|+\tau .
\end{equation}
The true binary64 relative unit roundoff is smaller than $u$, and the
additive allowance is much larger than half the least positive
subnormal. Overflows are excluded for the certificate evaluations by
the scale bounds proved below; they need not be excluded in an
untrusted solve.
\begin{lemma}[Accumulation and norm bounds]\label{round:accumulation}
For $m\le 2000$, a dot product on exact stored inputs, formed as
separate products and successive additions, satisfies
\begin{equation}\label{round:dot}
\left|\operatorname{fl}\!\left(\sum_{i=1}^m a_i b_i\right)
-\sum_{i=1}^m a_i b_i\right|
\le \gamma_m\sum_{i=1}^m|a_i b_i|+4(m+1)\tau .
\end{equation}
An initial addend is included in the absolute sum and increases the
small operation-count constant by one. A product of three factors
instead of two increases that constant by one more. For the
square-sum/square-root norms used here,
\begin{equation}\label{round:norm}
\left|\operatorname{fl}(\|a\|)-\|a\|\right|
\le 2.2\cdot 10^{-16}(m+2)\|a\|+10^{-140}.
\end{equation}
If the computed components approximate other components, their
Euclidean difference is added to this bound.
\end{lemma}
\begin{proof}
Iterating $|z_k|\le (1+u)|z_{k-1}|+u+b\tau$ for a length-$m$ dot product
with the standard addition/rounding analysis
(Lemma 3.1 of Higham's 2002 book) gives an interpolation of
the form $(1+u)^m\sum_{i=1}^m|a_ib_i|$ plus a $O(m\tau)$ term. The bound
$\gamma_m=2.1\cdot 10^{-16}(m+1)$ is a conservative choice: for
$m\le 2000$ one has $(1+u)^m-1\le mu\le 2m\cdot 10^{-16}$, which is
smaller than $\gamma_m$. Direct substitution of the accumulated
additive terms into the recurrence produces the $\tau$ bound. The
norm inequality is the same argument applied to the square-sum and
then to the correctly rounded square root, with $m+2$ factors and a
final $10^{-140}$ additive allowance.
\end{proof}
% ---------------------------------------------------------------
\medskip\noindent\textbf{Pi and exponential evaluations.}
\label{round:exponentials}
% ---------------------------------------------------------------
The constant $\pi$ is stored as the correctly rounded binary64 literal.
Its absolute error is at most $10^{-16}$ and its relative error is at
most $u$. The constants $\sqrt 6$ and its reciprocal are formed by the
indicated square-root and division; $2\cdot 10^{-15}$ is an ample
absolute error bound for each.
The source evaluates cosine through degree $18$ and sine through degree
$19$. In the $j$th nonconstant term, the error in $-y^2$ is reused $j$
times, and there are $j$ further multiplications and $j$ divisions.
The small integer denominators are exact. Its relative error is
therefore at most $(1+u)^{3j}-1<5.5\cdot 10^{-15}$ for $j\le 9$. For
additive underflow errors, compare with a recurrence using the same
multiplicative rounding factors but zero additive errors. Its terms
have magnitude below two, and each recurrence contracts the previous
additive error by a factor less than $1/2$: the numerator factor has
magnitude below $.786^2$ and the smallest denominator is two. Errors
from that factor, the product, and the division give the conservative
recurrence $E_j\le E_{j-1}/2+4\tau$, $E_0=0$. Thus $E_j<8\tau$; after
both coordinates' nine additions their total contribution is below
$10^{-295}$, including when $y^2$ or a later term underflows. The
absolute even-term sum at $.786$ is below $1.4$, its nonconstant part
below $.4$; the odd-term sum is below $1$, its part after $y$ below
$.2$. For the nonconstant even terms the successive ratio after the
first is at most $.786^2/12$, and for the odd terms after $y$ it is at
most $.786^2/20$. Their geometric majorants give the claimed absolute
sums without evaluating a transcendental function. The first omitted
cosine and sine terms are below $3.33\cdot 10^{-21}$ and
$1.25\cdot 10^{-22}$, respectively; both are below $10^{-16}$ even in
the worst scale considered here. Thus
\begin{equation}\label{round:exp-budget}
|\operatorname{fl}(e^{iy}) - e^{iy}| \le 10^{-14}
\end{equation}
for all real $y$ with $|y|\le 20$ and no overflow or underflow.
For $z=a+bi$ with $|a|,|b|\le20$, the pair
$(\operatorname{fl}(\cos b),\operatorname{fl}(\sin b))$ is multiplied by
$e^a$; the pair error and the multiplication error are each below
$2u$ in relative terms, and the total error is below $10^{-14}$
in absolute value at the scales used here (bound: ).
All relied-on arguments satisfy the range condition. Seed and dyadic
unitary phases are bounded near $[0,1]$; \texttt{gen} checks $|p_j|\le4$;
Candidate checks $|P_j|\le .51$ before its exponential; C validates
$|p_j|\le 4$ and checks ball phases of magnitude at most two. The
largest C difference has magnitude at most six.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Dephasing and projection.}
\label{round:dephasing}
% ---------------------------------------------------------------
For $|p_j|\le 4$, the three operations in unwrapped dephasing
$p_{ij}-p_{i0}-p_{0j}+p_{00}$ have entry error below $8\cdot 10^{-15}$;
their intermediate magnitudes are at most $8$, $12$, and $16$.
Subtracting the integer selected by \texttt{nearbyint} gives error
below $10^{-14}$ against the exact dephase minus that same integer. A
wrapped difference of two such arrays has error below $3\cdot 10^{-14}$
per entry against the ideal difference with the combined selected
integers. The one-sided Candidate version, whose other array is already
a stored dephased phase, has the same bound. These comparisons do not
assume identical rounding choices for an independently formed ideal
phase.
For a stored wrapped array, $|z_j|\le 1/2$. Its row and column sums
have error below $4\cdot 10^{-15}$, and its length-36 grand sum below
$1.4\cdot 10^{-13}$. Dividing by $6$ or $36$ and combining the means
gives projection-operation error below $10^{-14}$ per component.
Projection contracts the prior input error. Consequently the computed
projected difference has Frobenius error below $10^{-11}$. All distance
square and norm errors at the successful cutoffs below are negligible
against $10^{-7}$.
Thus \texttt{Atlas::find}'s returned $\sqrt{s}+10^{-7}$ is an upper
distance for the selected ideal projected lift. A merge uses radius
$r+\text{distance}+10^{-8}$ and checks a cap with $2\cdot 10^{-8}$
slack. Its sum-rounding error is below $10^{-14}$. The same arithmetic
applies to \texttt{assign} and to the local sphere-radius additions.
Searching only some buckets or transforms may miss a merge; it cannot
certify a merge that fails these tests.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Seed filtering and seed-ball reuse.}
\label{round:seed-filters}
% ---------------------------------------------------------------
Here $n=112$ and the stored $w$ approximates $1/n$. Midpoint-gap
calculations contain only small integer arithmetic and a few
products/subtractions; \eqref{round:basic} bounds their total error,
including the comparison threshold, by $10^{-14}$. Their slack is
$10^{-8}$. The orientation filter is integral.
For a list of five varying phases, the computed residual $c$ in
\texttt{sumok} has component error below $2\cdot 10^{-13}$. Indeed its
five exponential errors, after division by $2\pi$, total less than
$4.8\cdot 10^{-14}$ in complex modulus, and the divisions, additions,
and constant error cost less than $3\cdot 10^{-15}$. The sharper
componentwise bound $1.1\cdot 10^{-13}$ will also be useful. Each test
direction is regarded as an exact stored pair; its norm is less than
$1.1$. The five column perturbations, the residual perturbation, and
ordinary arithmetic in \texttt{chk} cost together less than $10^{-11}$.
Its remainder coefficient approximates $5\pi/(4n^2)$ within $10^{-15}$.
The actual $10^{-8}$ comparison addition therefore preserves every
feasible list.
In \texttt{supdual}, all proposed entries are checked by
$|v_a|\le 10^4$, with a negated comparison that also rejects NaN. They
are then exact stored test coefficients. The four right-hand-side
errors against $-nc$ are below $3\cdot 10^{-11}$ each. The signed term
consequently has input error at most
\[
2\sum_{a=1}^4 |v_a|\,3\cdot 10^{-11}\le 2.4\cdot 10^{-6}<2.5\cdot 10^{-6}.
\]
Its absolute term sum is below $10^7$, so its arithmetic error is below
$10^{-7}$. Each of the nine row-dot residuals has error below
$4\cdot 10^{-9}$: its four input contributions are at most
$4\cdot 10^4\varepsilon_{\exp}$, and its short dot costs less than
$10^{-10}$. Each pair norm has error below $10^{-10}$ by
\eqref{round:norm}. The remainder coefficient uses $1.251$ instead of
$1.25$, which leaves positive slack after its own rounding. The final
absolute scale is below $2\cdot 10^7$. Summing these contributions shows
that the ideal support expression is at most the computed proposal plus
$10^{-4}$. The actual addition is $0.01$. If no usable proposal improves
the initial $10^{10}$ bound, that initial bound already exceeds the
trivial support.
Integer differences lie in $[-111,111]$, including those stored in
\texttt{int8\_t}; all their squared/absolute sums are exact and small.
The comparison radius is
\[
\texttt{RN2}=\operatorname{fl}((\operatorname{fl}(.032\cdot112))^2)
-.0001
\]
with the source's indicated operation ordering. The total possible
error before the subtraction is below $10^{-12}$, so this is smaller
than the ideal $(.032\cdot112)^2$. This reasoning requires the
\texttt{pow(x,2)} calls in \texttt{covout} to be lowered to binary64
multiplication, as checked for the fresh build documented in
. It is not a claim about an arbitrary library
implementation of \texttt{pow}. The other square-power value $T$,
bucket reach, and integer cuts only select reuse attempts; missing an
attempt does not discard the current bin pair. A newly inserted midpoint
needs radius only $1.5/112<.032$. Finally, the nine division-rounded
seed phases differ from their rational midpoints in norm by less than
$10^{-14}$, which is absorbed when the A seed radius is increased to
$.0320001$.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Frames, comparison blocks, and A exclusion tests.}
\label{round:frames}
% ---------------------------------------------------------------
All ideal quantities in this subsection use exact unit exponentials at
the specified centers, and the same coordinate indices selected by the
code. We never need the selected coordinate to maximize the \emph{ideal}
projection score.
For \texttt{frame}, the four-term sum $S$ has error below
$2.5\cdot 10^{-13}$; its squared modulus has error below
$2.2\cdot 10^{-12}$, and its modulus below $2.6\cdot 10^{-13}$ by
\eqref{round:norm}. If the stored $m=\sqrt{4-|S|}$ passes $m>1.25$, the
ideal $m>1.24$. Using $|\sqrt a-\sqrt b|=|a-b|/(\sqrt a+\sqrt b)$ at
these lower bounds gives errors below $3\cdot 10^{-13}$ in $m$ and
$\mathrm{top}=\sqrt{4+|S|}$. In particular $16-|S|^2>6.1$; its error is
below $2.3\cdot 10^{-12}$, so its reciprocal has error below
$7\cdot 10^{-14}$. The two pseudoinverse numerators have modulus at most
eight and error below $5.1\cdot 10^{-13}$. Division/scaling then gives
\begin{equation}\label{round:frame-g}
|\Delta G_{ij}|<10^{-12},\qquad |G_{ij}^{\rm ideal}|<.81,
\end{equation}
the second inequality following from the lower singular bound.
The first constant basis column is exact. The next column has numerator
error below $1.3\cdot 10^{-13}$ and denominator
$\sqrt{4-|S|^2/4}>1.5$ with error below $2\cdot 10^{-13}$. Its ideal
entries have modulus at most one, so the quotient error is below
$4\cdot 10^{-13}$. For each of the remaining two columns, use the actual
selected coordinate projection $q$. In the first such projection the
earlier two products give entry error below $9\cdot 10^{-13}$, and its
norm error is below $2\cdot 10^{-12}$. The accepted computed norm is
greater than $.449999$, so the ideal projection is nonzero and
normalization gives entry error below $7\cdot 10^{-12}$. In the last
projection the extra column-product error gives entry error below
$1.6\cdot 10^{-11}$ and norm error below $3.5\cdot 10^{-11}$. The same
denominator check gives entry error below $1.2\cdot 10^{-10}$. This
argument inductively constructs the ideal orthonormal complements with
the selected coordinates. The two-product kernel projector has entry
error below $3\cdot 10^{-10}$; multiplying two complement entries and
$\sqrt6$ gives a $W$ entry error below $8\cdot 10^{-10}$.
The following table records the rest of the propagation, with the source
array names. Its errors include ordinary operation rounding.
\begin{center}
\small
\begin{tabular}{>{\raggedright\arraybackslash}p{.22\linewidth}p{.23\linewidth}p{.45\linewidth}}
\hline
Quantity & Absolute error per entry & Ideal scale used\\
\hline
\texttt{base} & $10^{-11}$ & two factors of modulus at most 2 times G\\
\texttt{fixed} & $2\cdot 10^{-9}$ & projector times a factor at most 1.62\\
\texttt{bm}, \texttt{cm} & $2\cdot 10^{-12}$ & unit entry times G\\
\texttt{uc.U}, angular derivatives & $2\cdot 10^{-13}$ & products of unit phases/trig entries\\
\texttt{K}, \texttt{D} & $4\cdot 10^{-9}$ & four W products and base\\
\texttt{mag} & $4\cdot 10^{-8}$ & $|D|<4.6$\\
$\sqrt{\texttt{mag}}$ & $5\cdot 10^{-9}$ & norm Lipschitz\\
\texttt{jac} & $10^{-8}$ & ideal modulus at most 6\\
flat real J & $1\cdot 10^{-7}$ & dot of D and jac\\
flat c & $5\cdot 10^{-8}$ & $(|D|^2-1)/2$\\
pre-equation c & $4\cdot 10^{-13}$ & short unit-exponential sums\\
pre-equation coefficients & $6\cdot 10^{-14}$ & unit exponentials\\
\hline
\end{tabular}
\end{center}
\medskip\noindent\textbf{The horizontal-projector identity..}
An entry of the exact double-centering projector is the product of two
centering-matrix entries. The computed \texttt{PI} differs by less than
$10^{-15}$ per entry. In the matrix $L_P$, each row has only 12 nonzero
entries, with error below $1.3\cdot 10^{-13}$ from the edge products.
Forming \texttt{LP} gives error below $5\cdot 10^{-12}$ per entry
against $L_P\Pi$; the ideal entry magnitude is at most $8.4$. In the
residual $\Pi-TT^*-ML_P\Pi$ each computed entry then has error below
$10^{-9}$. The propagated $LP$ error is at most
$30\cdot 2\cdot 5\cdot 10^{-12}=3\cdot 10^{-10}$; the remaining PI,
$T$-product, and short accumulation errors fit in the remainder of the
budget, with absolute sum below 560. The 1296 entry errors have
Frobenius total below $3.6\cdot 10^{-8}$. Norm accumulation is below
$10^{-8}$ even with residual norm 20161. Thus the actual
$2\cdot 10^{-6}$ addition in \texttt{id} suffices.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Root certificates and normal constants.}
\label{round:root}
% ---------------------------------------------------------------
The tangent operator has the certified bound
\begin{equation}\label{round:t-operator}
\|T\|_{\rm op}\le\sqrt{1+\mathtt{ort}}<1.00001,
\end{equation}
since $T$ is a genuine Root. The evaluation of $T^*T-I$ for the stored
$T$ (each $|T_{ja}|\le 1.01$) yields $\mathtt{ort}<10^{-5}$; the
evaluation of $(I-\Pi)T$ yields $\mathtt{hor}<10^{-5}$; the evaluation
of $F(P)$ yields $\mathtt{froot}<10^{-6}$.
For the auxiliary maps $M$: the source forms $M=H_VL_V^*$ with at most
$25$ columns of $H_V$ and $L_V$, each checked to magnitude at most $2$.
Absolute row sums of the column Grams of $H_V$ and $L_V$ give upper
bounds on $\|H_V\|$ and $\|L_V\|$; their product overestimates
$\|M\|$. Adding the $2.1\cdot 10^{-11}$ Frobenius bound on the
difference $H_VL_V^*-M$ (where $M$ denotes the accumulated stored
product) gives a valid $\mathtt{mnorm}$.
The two candidate choices for $C$ are $\sqrt 8\|M\|$ and
$\sqrt 6\|\mathcal S\|$; the source takes the smaller of upper bounds
for these two quantities. Each is a valid $\pi b$ in
Lemma~\ref{phase:remainder} by Lemma~\ref{phase:remainder-constants}.
The root evaluation thresholds are $\mathtt{val}<10^{-18}$ for the
overall pass, $\mathtt{ort}<10^{-5}$, $\mathtt{hor}<10^{-5}$,
$\mathtt{id}<10^{-4}$ after the $2\cdot 10^{-6}$ addition, and
$\mathtt{mnorm}\cdot\mathtt{froot}\le\mathtt{n0}$.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Candidate conversion and ellipsoidal coverage.}
\label{round:candidate}
% ---------------------------------------------------------------
These coverage estimates are needed only after a successful Candidate
preparation, at final unitary-grid levels $N_u\ge 64$.
Let $m_j=|D_j|$ denote an ideal comparison-block magnitude. The stored
magnitude test $\sqrt{\mathtt{mag}_j}>.15$ implies $m_j>.14999998$.
Normalizing $D$ by its computed magnitude has unit-direction error
below
\[
\frac{2\cdot 4\cdot 10^{-9}}{.14999998}
+\text{norm/division rounding}<5.5\cdot 10^{-8}.
\]
Candidate additionally checks a computed chord smaller than $10^{-7}$
between this normalized value and the exponential of its proposed
phase. With \eqref{round:exp-budget}, the true unit chord error is less
than $1.56\cdot 10^{-7}$. If a circular phase distance in cycles is
$v\in[0,1/2]$, its chord is $2\sin(\pi v)\ge 4v$. Across $16$ entries
the phase-center correction therefore has norm below
$1.56\cdot 10^{-7}$, which is below the $2\cdot 10^{-6}$ allowance used
in . The nine seed entries use exact stored
seed centers and an exact integer subtraction.
For the D-position coefficients $X$, the normalized-direction error
times the Jacobian bound six, plus the $10^{-8}$ Jacobian error and
short-dot rounding, is below $10^{-6}$ per entry. The other used $X$
entries are exact. The $\mathtt{mi}$ subtraction of $10^{-6}$ therefore
gives a positive magnitude lower bound greater than $.148$, and
$\mathtt{invmag}_j=1/(\widehat m_j-10^{-6})+10^{-9}$ is an upper bound
for $1/m_j$, below $6.8$.
The flat-sum error is below $10^{-6}$, so the $10^{-5}$ addition in
\texttt{slack}, followed by its additional $10^{-10}$, gives the
required nonnegative weighted chord-square upper bound. In the passed
branch $q<3.8$ and $1-q/4>.049$. The arithmetic in dividing
\texttt{slack} by \texttt{mi} has error below $10^{-12}$:
\texttt{slack}$<4.1$ and the denominator is above $.148$. The ratio
$(1-q/4)^{-1/2}$ has error below $10^{-12}$ before its $10^{-7}$
addition; the denominator lower bound bounds the derivative, and
$q/4$ is a binary scaling. The remaining positive products/roots and
the division by $(2\pi)^2$ in \texttt{d2} each have error below
$10^{-12}$ before their $10^{-8}$ additions. In particular the stored
quantities are upper bounds with
\begin{equation}\label{round:candidate-scales}
\mathtt{arc}<8.8,\qquad \mathtt{weighted}<10,
\qquad \mathtt{invmag}_j<6.8.
\end{equation}
\medskip\noindent\textbf{Tangent offsets and external errors..}
The computed \texttt{off} approximates $2\pi T^*\Delta$ within $10^{-9}$
in norm, by~\eqref{round:t-operator}, the projected-difference error,
and short-dot rounding. The error in the whole $X$ array is at most
$\sqrt{208}\cdot 10^{-6}$ in Frobenius norm; there are 13 columns and
16 affected positions. Multiplication by $T^*$ and the less than
$2\cdot 10^{-13}$ rounding of each dot give
\begin{equation}\label{round:ab-error}
\|\mathtt{AB}-T^*X_{\rm ideal}\|<1.5\cdot 10^{-5}.
\end{equation}
The parameter-vector norm is below
$\sqrt{.203^2+4(\pi/64)^2}<.226$, so substitution of stored AB costs
less than $3.4\cdot 10^{-6}$ radians. Formation of the corner offsets
adds less than $10^{-11}$, using coefficients of magnitude below 25.
The same bound applies to convex combinations of these corners.
For the second covariance, use stored \texttt{weighted} and
\texttt{invmag}, with rational $\beta=.8$. Each $D$-weighted dot has
absolute sum less than $16(1.00001)^2\cdot 6.8<110$. Its
two-product-per-term accumulation has error below $5\cdot 10^{-13}$ by
\eqref{round:dot} with one extra factor. The multiplying coefficient is
below $225$. Coefficient translation and the subsequent arithmetic,
even at absolute intermediate scale $26000$, give total entry error
below $10^{-9}$. Both $G$ arrays are explicitly stored symmetrically.
The added $10^{-7}I$ dominates $k$ times either entry-error bound for
$k\le 10$. Hence the exact stored $G$ is positive semidefinite and
dominates the required covariance.
% ---------------------------------------------------------------
\medskip\noindent\textbf{The stored spectral upper bound.}
\label{round:ell}
% ---------------------------------------------------------------
\begin{lemma}\label{round:spectral}
For an exact stored symmetric $k\times k$ matrix $G$ with $k\le 10$ and
$|G_{ij}|\le 10$, \texttt{spectralUpper} returns an upper bound for
$\|G\|$.
\end{lemma}
\begin{proof}
An absolute row sum has error below $3\cdot 10^{-13}$. The code adds
$10^{-10}$ to their maximum. For its other bound, put
$m=\max_{ij}|G_{ij}|$ and $A=\|G^2\|_F$. Because one diagonal sum of
squares contains a term $m^2$, $A\ge m^2$. Then
$|\operatorname{fl}(\sqrt{\mathtt{fro4}})-A|<4\cdot 10^{-13}A+10^{-139}$
by dot-product and norm rounding. Taking $\sqrt A$ costs less than
$4\cdot 10^{-13}\sqrt A+10^{-69}$, and since $\sqrt A\le 100$ the total
possible downward error is below $10^{-10}$. The code adds
$2\cdot 10^{-7}$. The exact $\sqrt{\|G^2\|_F}$ bounds $\|G\|$. Both
returned candidates are upper bounds, so their minimum is also an upper
bound, and its value is nonnegative and at most $100.1$.
\end{proof}
The \texttt{ell} guard checks symmetry and $|G_{ij}|\le 10$, and each
offset component has magnitude at most ten. Its invalid return value
$100$ cannot certify a cell: after conversion to cycles it exceeds $10$,
so the associated bootstrap returns its nonsuccess sentinel. A different
valid bound may still succeed. For valid data the initial triangle bound
is $\sqrt{\mathtt{mx}}+\sqrt{\mathtt{gu}}+10^{-8}$, with $\mathtt{mx}$
the maximum squared offset norm. Its norm error is below $10^{-13}$, and
the triangle bound is smaller than $43$.
For each proposed resolvent parameter $\gamma$: $0<\gamma<418$, and
\[
\mathtt{gap}=\operatorname{fl}(\gamma-\mathtt{gu}-10^{-9})
\le\gamma-\lambda_{\max}(G).
\]
Only gaps at least $10^{-6}$ are used. Taking $\max$ with zero and over
all offsets gives \texttt{high}. The ellipsoid formula uses
$\sqrt{\gamma(1+\mathtt{high})}+10^{-6}$; the added $10^{-6}$ and the
whole-candidate error bound $10^{-12}$ dominate at the used scales.
% (bibliography moved to end of file)
% ---------------------------------------------------------------
\medskip\noindent\textbf{The base grid and the outside test.}
\label{round:base-grid}
% ---------------------------------------------------------------
Use the stored positive width $w_B$ as an exact grid parameter. The
allowed range is $10^{-4}\le w_B\le 0.05$. The stored grid radius
$E=\sqrt{1+\mathtt{ort}_P}\,R+10^{-7}$ exceeds the true tangent radius
by more than $0.999\cdot 10^{-7}$. For
$L=\left\lceil E/w_B+0.500001\right\rceil$, the quotient scale is at
most $1261$ and its error is below $10^{-12}$. The $0.000001$ excess
over the half-cell term is much larger, so the integer base range
covers $[-E,E]$ in every coordinate. In particular $L\le 1262$.
Exact base centers have magnitude below $0.203$ per coordinate. Their
computed products have error below $10^{-16}$. Halving a stored width
is exact here; after at most six subdivision levels, each computed
center differs from the intended exact dyadic child center by less
than $10^{-15}$ per coordinate. The tangent-restriction identity
$\|T^*\delta\|\le\|(\Pi T)^*\delta\|+\|(I-\Pi)T\|\,\|\delta\|$ is used
to reduce the parent-distance test to a single comparison against $E_P$,
without loss elsewhere.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Single-vector filters and vector-ball construction.}
\label{round:single}
% ---------------------------------------------------------------
The C stage validates its Root and $0\le r\le 0.03$. At every bin
resolution $M\ge 6$, the halfwidth, Euclidean bin radius, one-norm
bound, and infinity-norm bound are formed as
\[
\begin{aligned}
a&=0.5/M+10^{-12},&
X&=\sqrt{5-1/6}\,a+10^{-12},\\
x_1&=5a+10^{-12},&
x_\infty&=1.5a+10^{-12}.
\end{aligned}
\]
Their operations have error below $10^{-14}$, and the $10^{-12}$
additions are upward. They satisfy $X<0.19$, $x_1<0.42$, and
$x_\infty<0.126$. The stored tangent and normal bounds
\[
r_t=\sqrt{1+\mathtt{ort}}\,r+10^{-9},\qquad
r_n=\mathtt{n0}+\mathtt{id}\,r+Cr^2+10^{-9}
\]
are upper bounds, with $r_t<0.030001$ and $r_n<0.019$. The
square-root/product and positive-polynomial errors are below
$10^{-13}$. The constant $K=1.00001$ exceeds the operator bound
regardless of literal translation; the Root residual gives
\[
\|e(P)/\sqrt6\|
\le\sqrt{1+(2\pi\sqrt2/6)\,10^{-6}}<1.00000075<K.
\]
Storage and elementary arithmetic are bounded directly: the phase
magnitude is below $4.0$; the six-term $S$ error is below $6\cdot
10^{-14}$; $|\,|S|/\sqrt 6-4\,|<10^{-13}$. The arrays $J$ and $B$
have errors below $3\cdot 10^{-11}$ and $4\cdot 10^{-11}$. The
rounded constants $a,X,I_M,L_M,1/\sqrt6,(2\pi)^2$ each carry unit
relative error $u$, and the two required products and the sum carry
at most $3u$ relative. Each bin-face primitive check has constant
$2\cdot 10^{-7}$ above
$r_n\|J\|+(r_t/4)$, using $|x_j|\le a$; its threshold is $10^{-6}$.
The combined bound uses the same $2\cdot 10^{-7}$. The nested
squared-magnitude crossing band is $2\cdot 10^{-7}$; the squared
increment is $7\cdot 10^{-7}$.
\medskip\noindent\textbf{The directional single-vector support..}
The heuristic solve is accepted only with finite coordinates and a
maximum magnitude between $10^{-100}$ and $10^{100}$. After division
by that maximum, regard the resulting stored $|y_i|\le 1$ as the test
coefficients. Define exact stored arrays
\[
\begin{aligned}
a_y&=\sum_i y_i c_i,& h_a&=\sum_i y_iB_{ia},\\
z_j&=\sum_i y_iJ_{ij},&
\widetilde z_j&=\sum_i y_iD^{\rm re}_{ij},\\
N_J&=\sum_{i,j}y_i^2J_{ij}^2,& N_S&=\sum_i y_i^2|S_i|^2,\\
M_S&=\max_i y_i^2|S_i|^2,& m_y&=\max_i|y_i|,\\
\quad M_z&=\max_j\bigl(z_j^2+\widetilde z_j^2\bigr).&
\end{aligned}
\]
The support compared with $|a_y|$ is the stored evaluation of
\begin{align}
&a\sum_{j>0}|z_j|+r_t\|h\|+r_n\sqrt{N_J}\nonumber\\
&\quad+\pi\Bigl[
m_y(KX+r)^2+
\min\bigl\{\sqrt{M_z}\,X,\ K\sqrt{N_S}\,x_\infty\bigr\}X
+2\frac{\sqrt{N_S}}{\sqrt6}Xr
+\frac{\sqrt{M_S}}{\sqrt6}r^2 \Bigr].
\label{round:c-single-support}
\end{align}
This is the necessary support bound of the Exclusion stage.
Its coefficient errors against ideal arrays satisfy
\[
|\Delta a_y|<7\cdot 10^{-12},\qquad
|\Delta z_j|,|\Delta\widetilde z_j|<7\cdot 10^{-12},\qquad
|\Delta h_a|<4.3\cdot 10^{-11}.
\]
The first two use six input errors and a short dot; the last uses six
$B$ errors and the same dot estimate. The two bracketed products cost
less than $10^{-14}$ each; the short absolute sums and the
normalization by $m_y\ge 10^{-100}$ cost less than $10^{-15}$.
Squares and square-root evaluation at the used scales cost less than
$10^{-15}$. Three times either dot error is below $3\cdot 10^{-11}$.
Thus the certified support contribution is within $10^{-11}$ of the
ideal expression inside $\pi$, and the whole error against the ideal
support is below $10^{-10}$. The code adds $2\cdot 10^{-7}$.
% ---------------------------------------------------------------
\medskip\noindent\textbf{The joint coefficient support.}
\label{round:joint}
% ---------------------------------------------------------------
The joint test is separated over the three tangent pieces and has at
most $6+6+2=14$ rows. The proposed direction is checked finite and
normalized as before, so its stored entries have magnitude at most one.
The constant pairing has error below $1.5\cdot 10^{-11}$. Each
accumulated shared-Hadamard coefficient has error below $9\cdot 10^{-11}$:
it receives at most twelve $B$ input errors of $7\cdot 10^{-12}$ plus
short-dot rounding. Each vector-phase coefficient, including its edge
contribution, has error below $9\cdot 10^{-12}$: six single-vector
entries and at most two edge entries each contribute at most $12u$.
Combining, the ideal support is at most the computed support plus
$3\cdot 10^{-10}$. The code adds $3\cdot 10^{-7}$.
The final comparison allowance of $3\cdot 10^{-7}$ is monotone: NaN and
infinity are rejected by the direction checks, so only finite proposals
enter. Enlarging the proposed direction or its stored coefficients can
only raise the computed support, so the comparison remains a rigorous
upper bound at the recorded threshold.
% ---------------------------------------------------------------
\medskip\noindent\textbf{Finite paths and NaN handling.}
\label{round:summary}
% ---------------------------------------------------------------
The preceding estimates apply to bounded finite certificate
evaluations. Exceptional heuristic values cannot bypass those
preconditions.
\begin{enumerate}
\item In the seed cover, all direct lists come from finite integer
midpoints and \texttt{expair}. The dual solver's entries must
pass the negated bounded comparison $|v_a|\le 10^4$; NaN and
infinity fail it. A failed proposal leaves the large safe
support bound. Gap/hash/reach cuts select only possible reuse
attempts.
\item For A, a failed frame normalization makes setup fail and the
entry point records an unresolved seed ball. It does not
certify the seed ball empty. The passed frames, fixed radius,
and dyadic unitaries give finite bounded block coefficients.
The potentially ill-conditioned solves only propose directions;
the final \texttt{trial::certificate} checks every direction
entry and its maximum before use.
\item The snapping iteration calls \texttt{gen}, which rejects
nonfinite or out-of-range phases. Root certification requires
$\mathtt{val}<10^{-18}$, checks proposed eigenvalues finite,
checks HV/LV before using them, and checks $T$ before its
certificate products. Those bounded arrays make the products
finite. The optional tensor change of basis is range-checked;
an unusable one cannot lower the previously valid $C$. The
final conjunction of strict Root thresholds also rejects NaN.
\item In B, failed cells are retained as bad parents; they do not
enlarge the covered region. The marker file is replaced only
after a successful atlas update, and the final count is the
terminal cover size. A failed temporary write records an
unresolved parent.
\item The cluster output is checked as in
at $|v_j|\le 2$ and $0\le s<.1$. The two edges share the
distance bound $r_e=s_u+s_v$. If the direct tests reject an
edge, its weight is zero with infinity cost. Otherwise the
stored $D_E$ and $q_E$ are finite. The joint coefficient
guard from then rejects any failed
proposal and keeps the corresponding previous bound. A failed
proposal cannot certify a pair excluded.
\end{enumerate}
\end{document}