% =====================================================================












usetikzlibrary{arrows.meta,positioning}
input{glyphtounicode}
input{glyphtounicode-cmex}
pdfgentounicode=1
















makeatletter
newcommand{theoremspace}{if@nobreakelseNeedspace{3baselineskip}fi}
makeatother








emergencystretch=1em



% =====================================================================
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 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.
 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 complex Hadamard matrix, and
 the question becomes: which H admit two further unbiased bases?
 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
 4times4 block of unknown phases.
 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.
 Complete and refine. For a seed ball, the 4times4
 block must satisfy Hadamard equations that compress to a single
 Vin 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.
 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.
 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.


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.

% =====================================================================
medskipnoindentIntroduction.
% =====================================================================

Two orthonormal bases B and C of C^d are
mutually unbiased (MUB) if
[
 |< b,c>|^2=frac 1d  for every binB, cinC .
]
Here < u,v>=sum_ju_jbarv_j.
For dge 2 let N(d) be the maximum size of a pairwise mutually
unbiased family of orthonormal bases in C^d. The main result is:

[Computer-assisted]
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.


The upper bound is the hard direction. The lower bound (three bases
exist) is an explicit computation, .

medskipnoindentMUBs, 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, HH^*=6I .
]
Such an H is called a 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~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.

medskipnoindentHistorical 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).

{smallhttps://www.pnas.org/doi/10.1073/pnas.46.4.570}

{smallhttps://doi.org/10.1088/0305-4470/14/12/019}

{smallhttps://www.sciencedirect.com/science/article/pii/0003491689903229}

{smallhttps://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.

{smallhttps://arxiv.org/pdf/quant-ph/0309120}

{smallhttps://www.gerhardzauner.at/documents/gzdiplthd.pdf}

{smallhttps://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.

{smallhttps://arxiv.org/pdf/quant-ph/0406175}

{smallhttps://arxiv.org/pdf/math.OA/0609076}

{smallhttps://arxiv.org/pdf/0901.4051}

{smallhttps://arxiv.org/pdf/1103.1025}

The closest antecedent is the certified discretization of Jaming,
Matolcsi, M'ora, Sz"ollH 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.

{smallhttps://arxiv.org/pdf/0902.0882}

{smallhttps://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.

{smallhttps://arxiv.org/pdf/2503.14752}

{smallhttps://arxiv.org/pdf/2608.18053}

medskipnoindentThe strategy at a glance.

The proof has five layers, each of which strictly refines the previous
one and preserves every exact candidate solution.
 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).
 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 4times 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 2times 2
 matrix whose polar part is a unitary Vin 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~sec:phase-cover).
 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.
 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.
 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.


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 necessary condition with a certified margin; and no
step in the chain ever relies on a floating-point center being exact.

medskipnoindentThe 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.

[Three mutually unbiased bases in dimension six]

There exist three pairwise mutually unbiased orthonormal bases in
C^6.


Proof. 
Identify C^6 with C^2 ⊗ C^3 in the standard product
coordinates. Define
[
A_0=I_2,
A_1=frac 1{sqrt2}1&11&-1,
A_2=frac 1{sqrt2}1&1 i&-i.
]
These are unitary, and for rne 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 121+i&1-i1-i&1+i.
]
Let omega=e^{2pi i/3} and set
[
F=frac 1{sqrt3}(omega^{jk})_{j,k=0}^2,
D=diag(1,omega,omega),
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 13sum_{x=0}^2omega^{x^2+(k-j)x}.
]
If t=k-j is read modulo 3, the sum is 1+2omega 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 rne s.

For r=0,1,2, let the columns of U_r=A_r ⊗  B_r form the rth
basis. Each U_r is unitary. For rne s we have
U_r^*U_s=(A_r^*A_s) ⊗ (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.


medskipnoindentA Lie-algebra consequence.

An orthonormal basis B of C^6 determines its Cartan
subalgebra mathfrak h_frakh_{B} of traceless matrices that are
diagonal in that basis. It is spanned by the matrices bb^*-I/6 for
binB. The Killing form on sl_6(C) is
K(X,Y)=12tr(XY); two subalgebras are
Killing-orthogonal when this form vanishes on every pair of
elements. For unit vectors b,c,
[
 tr((bb^*-I/6)(cc^*-I/6))
 =|< b,c>|^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.

[Adjoint-closed Cartan subalgebras]

Under the arithmetic and compiler conditions of
the certified arithmetic, using the complete execution recorded above
and certified in Theorem~thm:main, the maximum
number of pairwise Killing-orthogonal Cartan subalgebras of
sl_6(C) that are all closed under the standard adjoint
X|-> X^*=overlineoverline X^X^{,T} is three.


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~thm:main
gives the upper bound under its stated contract, and
Lemma~lem:three-mub-six gives attainment.


{smallhttps://arxiv.org/pdf/quant-ph/0506089}

% =====================================================================
medskipnoindentEquivalences and an exhaustive seed cover.

% =====================================================================

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"ollH 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.

{smallhttps://arxiv.org/pdf/0902.0882}

{smallhttps://arxiv.org/pdf/1201.0640}

Throughout, phases are in cycles and e(t)=exp(2pi i t)
entrywise. Matrix indices run 0,...,5, entries are stored in
row-major order with (i,j) in position 6i+j, and the Hermitian inner
product is < u,v>=sum_ju_jbarv_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.

% ---------------------------------------------------------------
medskipnoindentEquivalences and phase gauges.
% ---------------------------------------------------------------


A matrix HinC^{6times 6} is Hadamard if every entry has
modulus one and HH^*=6I. It is extendible if the standard
basis and the columns of H/sqrt6 belong to a family of four
mutually unbiased orthonormal bases.


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).

[Equivalences]
Extendibility is invariant under row or column permutations,
multiplication of rows or columns by unit scalars, complex conjugation,
and transposition.


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.


For a real 6times 6 matrix z define double centering by

 (Pi z)_{ij}
 =z_{ij}-frac 16sum_b z_{ib}-frac 16sum_a z_{aj}
       +frac 1{36}sum_{a,b}z_{ab}.

Its range is the horizontal subspace of matrices whose every
row and column sum is zero. We also use unwrapped dephasing

 D_0(P)_{ij}=P_{ij}-P_{i0}-P_{0j}+P_{00}.


[Horizontal phase representatives]
Pi is the orthogonal projection onto the horizontal subspace; in
particular |Pi z|_Fle |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.


Proof. 
Equation~(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~seed:equivalences. Integer additions leave
e(P) unchanged.


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.

% ---------------------------------------------------------------
medskipnoindentA symmetry arrangement for nine phases.
% ---------------------------------------------------------------

The 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.

[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

 x=(0,x_1,x_2,x_3,x_4,x_5),
 0le x_1le...le x_5le 1,
 g_x=1-x_5, x_1le x_5-x_4 .

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

 y=(0,xi,y_2,y_3,y_4,y_5), 0le y_2le y_3le y_4le y_5le 1 ,

where xi is the selected entry of x. These lists satisfy

 sum_{j=0}^5 e(x_j)=0,
 sum_{j=0}^5 e(y_j)=0, g_yge g_x .



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~(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, x'_i=x_5-x_{5-i}(1le ile 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~(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_yge g_x.


Set N=112 and w=1/N. Use closed bins

 I_b=[b/N,(b+1)/N], u_b=(b+tfrac12)/N, 0le b<N .

They cover [0,1] including both endpoints. Sorted exact phases admit
nondecreasing bin indices, so write the five indices for
x_1,...,x_5 as
b^X=(a,b,c,d,l) with 0le ale ble cle dle l<N. The orientation
condition implies the necessary integral restriction

 frac aNle x_1le x_5-x_4le(l+1-d)/(N),
  ale l-d+1 .


For this index tuple choose sin{0,...,4} as the first minimizer
of |2b^X_s+1-N|, the rule used in 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~(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
0le b^Y_1le b^Y_2le b^Y_3le 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.

[Stability of circular gaps]
If each point of one circular multiset is matched to a point of another
at circular distance ledelta, their gaps differ by le 2delta.


Proof. 
Take an empty open arc of length L in the first multiset. If L>2delta,
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-2delta.
If Lle 2delta the same bound follows from nonnegativity. Apply to a
largest gap and interchange the two multisets. No preservation of sorted
order is used.


[Necessary midpoint gap filters]
For every bin pair containing lists as in
Lemma~seed:arrangement, the midpoint gaps (with the fixed zero)
satisfy

 g_{x,mid}le 1-(l-1)/(N),
 g_{y,mid}ge 1-(l+2)/(N).



Proof. 
Each varying phase is within w/2 of its midpoint and the initial zero
is unchanged. Lemma~seed:gap-stability bounds each gap change by
w. Since x_5in[l/N,(l+1)/N] we have
1-(l+1)/Nle g_xle 1-l/N. Combine with g_yge g_x.


For sorted midpoint bin indices q_0,...,q_4 the gap used by
mk is
[
 maxBigl{(q_0+tfrac12)/(N),
  (q_1-q_0)/(N),...,(q_4-q_3)/(N),
 1-(q_4+tfrac12)/(N)Bigr}.
]
For a Y tuple it sorts a copy solely for this computation, preserving
the stored shared-first ordering.

% ---------------------------------------------------------------
medskipnoindentVanishing-sum filters and the shared displacement.
% ---------------------------------------------------------------

For either of the two lists x or y, let u=(u_0,...,u_4) be its
five varying midpoints. Define real two-vector c and real 2times 5
matrix J by

 c=frac{1+sum_{j=0}^4 e(u_j)}{2pi i}, J_j=e(u_j),

realified in both definitions.

[Necessary sum test]
If 1+sum_{j=0}^4 e(u_j+Delta_j)=0 and |Delta_j|le w/2 then

 |c+JDelta|le(5pi w^2)/(4).

In particular for every vinR^2,

 |v· c|le(5pi w^2)/(4)|v|
                  +frac w2sum_{j=0}^4|v· J_j|.



Proof. 
The second derivative of e(t) has modulus 4pi^2, so Taylor's formula
with integral remainder gives
[
 |(e(u_j+Delta_j)-e(u_j))/(2pi i)-e(u_j)Delta_j|
 lepiDelta_j^2 .
]
Sum the five remainders and use the vanishing-sum hypothesis to
obtain~(seed:taylor-residual). Dot with v, apply the triangle
inequality and |Delta_j|le w/2 to get~(seed:sum-support).


The function sumok tests~(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=NDelta, write

 E=(5pi)/(4N), |Jh+Nc|le E, |h_j|letfrac12 .


[Shared-coordinate intervals]
Let t be the position of a shared entry in either list. Every feasible
displacement satisfies, for every vinR^2,

 |Nv· c+(v· J_t)h_t|
 le E|v|+frac 12sum_{jne t}|v· J_j|.

If a=v· J_tne 0 and the right side is r, this is the interval
constraint

 -(Nv· c)/(a)-frac r{|a|}
 le h_tle
 -(Nv· c)/(a)+frac r{|a|}.

The intervals obtained for the shared entry from X and Y must
intersect.


Proof. 
Separate the shared term in~(seed:grid-error), take its dot
product with v, and bound each remaining coordinate by 1/2. Solving
the resulting absolute-value inequality gives~(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.


The function 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 seeds::build now has an unambiguous
meaning. It visits every weakly increasing five-bin X tuple and
applies~(seed:orientation-filter), the first inequality
in~(seed:gap-inequalities), and sumok. Separately, for
every possible first bin, it visits all weakly increasing four-bin Y
tails and applies sumok. The program 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).

% ---------------------------------------------------------------
medskipnoindentCertified reuse of seed centers.
% ---------------------------------------------------------------

The nine seed coordinates, before moving the selected column, are

 q=(x_1,x_2,x_3,x_4,x_5,y_2,y_3,y_4,y_5)inR^9 .

For a bin pair let uinR^9 be its rational midpoint and put
h=N(q-u). Thus hin[-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~(seed:residual-data). Define the real 4times 9 matrix J_*
by its action

 J_*h=
 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
 ,
 p_*=-Nc_X -Nc_Y.

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.

[Support bounds for center reuse]
Let u' be a previously retained midpoint with the same selected
position s, and let k=N(u-u')inZ^9. Every feasible
displacement in the current bin pair satisfies

 N^2|q-u'|^2=|k+h|^2
 le|k|^2+frac 94+2k· h .

Two upper bounds for the last term are

 2k· h&lesum_{j=0}^8|k_j|,
 
 2k· h&le
 sum_{j=0}^8|k_j-(J_*^Tv)_j|
 +2v· p_*+2E(|v_X|+|v_Y|),
  vinR^4,
 

where v_X,v_Y are the two real two-vector blocks of v.


Proof. 
Expand the square in~(seed:distance-bound) and use
|h|^2le 9/4. The coordinate bounds give~(seed:box-support).
For the other bound write
[
 2k· h
 =2(k-J_*^Tv)· h+2v· p_*+2v· 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|).


The function supdual proposes vectors v and evaluates the
right side of~(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.25pi/N by the larger
1.251pi/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, testBucket first
tries~(seed:box-support) and then, if applicable, the dual bound.
It returns success only when the upper bound in~(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~(seed:nine-coordinates), though its displacement enters both
equations in~(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.

% ---------------------------------------------------------------
medskipnoindentThe exhaustive seed cover proposition.
% ---------------------------------------------------------------

[Exhaustive seed cover]
For a completed execution of the specified 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~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.


Proof. 
Given an exact Hadamard matrix, first apply
Lemma~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~seed:gap-filters, seed:sum-test, and
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~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
[
 (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.


% =====================================================================
medskipnoindentCompletion over a fixed seed ball.

% =====================================================================

This section reduces "complete the seed phase ball to a full dephased
Hadamard matrix" to a covering of the 2times 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.

% ---------------------------------------------------------------
medskipnoindentSeed 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

 E_0=1&11&a,
 B_0=1&1&1&1b_0&b_1&b_2&b_3,
 (C_0)_{i,:}=(1,c_i)  (0le i<4).

A feasible completion is a unit-modulus matrix
[
 H'=E'&B' C'&D',
 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},...,s_{b_3},s_{c_0},...,s_{c_3})inR^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

 (7,8,9,10,11,13,19,25,31).

Let r_s be the enlarged seed radius in cycles and R_s=2pi 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

 e_1^2+(b_1')^2+(c_1')^2le|s|^2le R_s^2.


[Fixed frames]
Let pinC^4 have unit-modulus entries, S=sum_j p_j, and
M=({bf 1} p). The squared singular values of M are 4pm|S|. If
|S|<4 then

 M(M^*M)^{-1}
 =frac 1{16-|S|^2}(4{bf 1}-poverlineoverline SS,,;4p-{bf 1}S).

Under the same condition there is an orthonormal 4times 2 matrix K
with 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.


Proof. 
M^*M=4&S overlineoverline S&4,S&4, whose
eigenvalues and inverse give the singular-value formula
and~(completion:inverse-frame). The first two Gram--Schmidt
vectors span 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.


Apply Lemma~completion:frames with p=(overlineoverline b_j)_jb_j)_j and with
p=(c_i)_i. Throughout the seed ball keep the resulting frames fixed
and write

 G_B=B_0^*(B_0B_0^*)^{-1},
 G_C=C_0(C_0^*C_0)^{-1},
 ranK_B=ker B_0,
 ranK_C=ker C_0^*.

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

 P_K=K_BK_B^*, P_B=I-P_K=G_BB_0,
 P_C=I-K_CK_C^*=G_CC_0^*.

Also

 A_0=|E_0|_{rm op}=√(2+|1+a|),
 A_*=min(2,A_0+R_s).

Both A_0 and |E'|_{rm op} are at most A_*: indeed A_0le 2,
every entry of E' has modulus one so |E'|_{rm op}le |E'|_F=2,
and |E'|_{rm op}le A_0+e_1le 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.

% ---------------------------------------------------------------
medskipnoindentA 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.

[Polar comparison]
For a feasible completion let Y=K_C^*D'K_B. Then


 6I-Y^*Y&=K_B^*B'^*B'K_B+K_B^*D'^*P_CD'K_Bsucceq 0,

 tr(6I-Y^*Y)
 &le(b_1')^2+frac{(|E'|_{rm op}b_1'+sqrt6,c_1')^2}{m_C^2}le q,

where

 q=R_s^2(1+(A_*^2+6)/(m_C^2)).

If q<6 the polar factor V=Y(Y^*Y)^{-1/2} lies in U(2) and

 |Y-sqrt6 V|_Fle
 p:=frac q{sqrt6+√(6-q)}.



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
(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}lesqrt6, giving the first inequality
in~(completion:deficit-bound). Cauchy--Schwarz with
(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
√(6-ell_j) and its polar factor is unitary. Moreover
[
 sqrt6-√(6-ell_j)
 =frac{ell_j}{sqrt6+√(6-ell_j)}
 lefrac{ell_j}{sqrt6+√(6-q)}.
]
Taking the Euclidean norm of these two differences and using
√(ell_1^2+ell_2^2)le ell_1+ell_2le q
proves~(completion:polar-bound).


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
Vin U(2).

% ---------------------------------------------------------------
medskipnoindentThe unitary cells and their Taylor bounds.
% ---------------------------------------------------------------

[A closed covering of U(2)]
Every Vin U(2) has a representation

 V=e(alpha)
 
 e(beta)cos(2pi t)&e(gamma)sin(2pi t)
 -e(-gamma)sin(2pi t)&e(-beta)cos(2pi t)
 ,

with 0lealphaletfrac12, 0le tletfrac14,
0lebeta,gammale 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.


Proof. 
Choose a square root of det V having phase alphain[0,1/2]
in cycles. Then e(-alpha)Vin SU(2) has the form
z&w-overlineoverline w&w&overlineoverline zz
with |z|^2+|w|^2=1. Choose tin[0,1/4] with
|z|=cos(2pi t), |w|=sin(2pi 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.


Fix one cell and denote its midpoint matrix by U. If V is represented
inside the cell, let

 theta=(theta_alpha,theta_t,theta_beta,theta_gamma),
 |theta_j|le h:=pi/N_u

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
(completion:unitary-parameterization) uses cycle coordinates.

[Uniform cell Taylor estimates]

For the cell and displacement in~(completion:unitary-cell),

 |V-U|_Flesqrt6,h,
 |V-U-dU[theta]|_Flefrac{√(42)}2,h^2.



Proof. 
Follow the straight segment in the four angle coordinates, writing its
matrix as U(tau), 0letaule 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^2cos^2 T+w^2sin^2 T)le 6h^2 .
]
For the two diagonal entries, the squared magnitudes of the real
factors in the second derivative have sum

 &(((x+z)^2+y^2)^2+((x-z)^2+y^2)^2)cos^2 T
 &=2((x^2+z^2+y^2)^2+4x^2z^2)cos^2 T
 le 26h^4cos^2 T,

and the imaginary factors contribute
4y^2((x+z)^2+(x-z)^2)sin^2 Tle 16h^4sin^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^2le 42h^4. Integrating the first derivative and
then using the integral Taylor remainder with weight 1-tau proves
(completion:unitary-taylor-bounds).


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 uc::init and
uc::initDeriv.

% ---------------------------------------------------------------
medskipnoindentComparison with a cell center.
% ---------------------------------------------------------------

For the fixed seed and unitary cell define

 D=D_{rm span}+D_{rm ker},
 D_{rm span}=-C_0E_0^*G_B^*,
 D_{rm ker}=sqrt6,K_CUK_B^*, Z=D'-D .

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.

[Orthogonal 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


 ZP_B&=(-Delta C E_0^*-C'Delta E^*-D'Delta B^*)G_B^*,

 P_CZP_K&=G_C(-E'^*Delta B-Delta C^*D')P_K,

 (I-P_C)ZP_K&=K_C(Y-sqrt6,U)K_B^* .

Consequently every feasible completion associated to this unitary cell
satisfies

 |Z|_Fle z:=
 √(((R_s(sqrt6+A_*))/(m))^2+(6h+p)^2),.



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~(completion:piece-B).
The identity C'^*D'=-E'^*B', together with B_0P_K=0 and
P_CDP_K=0, gives~(completion:piece-C). The definitions of Y,
K_B, K_C give~(completion:piece-K).

Put r_{EB}=√(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

 |ZP_B|_F&lefrac{sqrt6,r_{EB}+A_0c_1'}{m_B},
 |P_CZP_K|_F&le(A_*b_1'+sqrt6,c_1')/(m_C)
 lefrac{A_*r_{EB}+sqrt6,c_1'}m .

The vector of these two norms is bounded componentwise by
[
 frac 1msqrt6&A_* A_*&sqrt6
 r_{EB} c_1'.
]
The displayed matrix has operator norm sqrt6+A_*, and its input norm
is at most R_s by~(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~(completion:z-bound).


[A row-block operator bound]
For every cell midpoint Uin U(2),

 |[C_0 D]|_{rm op}le
 L:=max(sqrt6, |C_0|_{rm op}√(1+(A_0/m_B)^2)).

Also |DP_K|_{rm op}=sqrt6 and

 |D|_{rm op}lemax(sqrt6, frac{|C_0|_{rm op}A_0}{m_B}).



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

 [C_0 D][C_0 D]^*
 =C_0[I+E_0^*(B_0B_0^*)^{-1}E_0]C_0^*+6K_CK_C^* .

The first summand has range in ranC_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~(completion:L-definition). Omitting the
C_0C_0^* contribution gives~(completion:D-operator-bound).
Finally DP_K=D_{rm ker}, whose two nonzero singular values are
sqrt6.


% ---------------------------------------------------------------
medskipnoindentA 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

 
 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^*.
 

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.

[Linear remainder bound]

Every feasible completion associated to the cell satisfies

 |Z-J(s,theta)|_Fle E_r,
 
 E_r=sqrt{e_{rm off}^2+
       (frac{√(6· 42)}2,h^2+p)^2},

where

 e_{rm off}=
 (R_s^2(L+A_0+1))/(2m)+(zR_s)/(m).



Proof. 
The first two terms of~(completion:linear-map) lie in the spaces
of the first two projections of
Lemma~completion:comparison-pieces; the last term lies in the
third. Bound the remainders in each.

First replace C',D' in~(completion:piece-B) by C_0,D. The
additional numerator terms are
-Delta C,Delta E^*-ZDelta B^*. Replacing E',D'
in~(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√((b_1')^2+(c_1')^2)}mle(R_s^2)/(2m).
]
The product terms involving Z have combined norm at most
[
 frac{|Z|_F√((b_1')^2+(c_1')^2)}mlefrac{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(sum_{j=1}^9s_j^4)^{1/2}
 lefrac 12|s|^2lefrac 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~completion:L-bound bounds the norms of their two linear
output pieces componentwise by
[
 frac 1mL&A_0 A_0&sqrt6
 x y.
]
For the second component use |DP_K|_{rm op}=sqrt6. Since
Lgesqrt6 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~(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+√(6· 42),h^2/2. This space is orthogonal to the first
two, so Pythagoras
proves~(completion:linear-remainder-bound).


The entrywise form for assembling~(completion:linear-map) is: for
0le i,j,k<4 the shared seed coefficient is
[
 -c_ioverlineoverline{ii a(G_B)_{j1}},
]
the b_k coefficient is
[
 -[(G_C)_{i0}+(G_C)_{i1}overlineoverline a]a]
       , i b_k(P_K)_{kj}
 -D_{ik}overlineoverline{ii b_k(G_B)_{j1}},
]
and the c_k coefficient is
[
 -mathbf 1_bf1_{i=k},i c_k
       (overlineoverline{(G_B)_(G_B)_{j0}}+overlineoverline{a(G_B)_a(G_B)_{j1}})
 -(G_C)_{i1}i c_kbar(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,floor( l/2)}
                 overlineoverline{(K_B)_(K_B)_{j,lbmod 2}}, 0le l<4,
]
and contract with each angular derivative of U. These are the fixed
and moving parts of block::setup, makeD, and
makeJ.

% ---------------------------------------------------------------
medskipnoindentNecessary flatness and separating inequalities.
% ---------------------------------------------------------------

We have bounded the total block displacement Z=D'-D by z and its
remainder after the linear comparison J(s,theta) by E_r.
The next step uses these bounds inside the necessary equations
|D'_{ij}|=1.

[Norm exclusion]
A feasible completion associated to a cell must satisfy

 sum_{i,j=0}^3(|D_{ij}|-1)^2le z^2 .



Proof. 
Every D'_{ij} has modulus one, and the reverse triangle inequality
gives ||D_{ij}|-1|le |D'_{ij}-D_{ij}|. Sum squares
and use~(completion:z-bound).


For the second test define two complex constants and their linear parts
by

 c_X&=(1+a+textstylesum_j b_j)/i,&
 J_Xs&=a s_a+textstylesum_j b_j s_{b_j},
 c_Y&=(1+a+textstylesum_i c_i)/i,&
 J_Ys&=a s_a+textstylesum_i c_i s_{c_i}.

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

 |rho_X|letfrac12 R_s^2, |rho_Y|letfrac12 R_s^2.

Take real and imaginary parts of these two equations and append the
following sixteen real constants and linear forms:

 c_{ij}=frac{|D_{ij}|^2-1}{2},
 J_{ij}(s,theta)
 =Re(overlineoverline{D_D_{ij}},J_{ij}(s,theta)).

Denote the resulting real constant vector by c_{rm all}inR^{20} and
the real linear matrix by J_{rm all}inR^{20times 13}, with the
first nine columns corresponding to s and the last four to theta.

[Directional exclusion]
Let yinR^{20} be arbitrary, split as two pairs y_X,y_YinR^2 and
a flat part y_DinR^{16}, and write g=J_{rm all}^*y=(g_s,g_theta)
inR^9timesR^4. Every feasible completion of the cell satisfies

 
 y· c_{rm all}le{}&
 R_s|g_s|+h|g_theta|_1+E_r|Dcirc y_D|_F
 &+frac 12 z^2max(0,max_{i,j}(-y_{D,ij}))
 +frac 12 R_s^2(|y_X|+|y_Y|).
 

Here circ is entrywise multiplication, with the sixteen flat
coordinates arranged as a 4times 4 matrix.


Proof. 
Write R=Z-J(s,theta), so |R|_Fle E_r. Expanding the
exact flatness equation gives
[
 0=frac{|D'_{ij}|^2-1}{2}
   =c_{ij}+J_{ij}(s,theta)
    +Re(overlineoverline{D_D_{ij}}R_{ij})
    +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· 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|Dcirc y_D|_F; the seed remainders contribute the last term
by~(completion:seed-sum-remainders).

The quadratic term has a fixed sign before weighting. After it is moved
to the right its contribution is
[
 -frac 12sum_{i,j}y_{D,ij}|Z_{ij}|^2
 lefrac 12max(0,max_{i,j}(-y_{D,ij}))
               sum_{i,j}|Z_{ij}|^2 .
]
Use |Z|_Fle z. This proves the remaining term and the stated
inequality.


[Dual certificates]
A cell contains no feasible completion if~(completion:norm-test)
fails, or if some real direction y violates
(completion:directional-test). The choice of y may be heuristic;
only its final inequality must be certified. Positive normalization
of y preserves the certificate.


Proof. 
Both inequalities are necessary for every feasible completion in the
cell. All terms of~(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.


The implementation stacks exactly these twenty equations in
eq. The routine trial::pass proposes directions by
regularized linear solves and weight updates; it excludes only when
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~(completion:directional-test). No accuracy or
convergence assertion about the proposing solves is needed. The used
Stage A path is Asearch::testA; the auxiliary
Usearch::test and 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.

% ---------------------------------------------------------------
medskipnoindentFinite recursion and unresolved cells.
% ---------------------------------------------------------------

[Completion coverage conditional on successful
termination]
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.


Proof. 
Given any feasible completion, Lemma~completion:polar produces a
unitary V. Lemma~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.


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.

% =====================================================================
medskipnoindentCertified phase balls and their atlas.

% =====================================================================

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 6times 6 matrices are identified with
R^{36} in row-major order; a displacement is horizontal if it
lies in ranPi. All vector norms are Euclidean, phase
matrix norms are Frobenius, and stars on real matrices mean transpose.

% ---------------------------------------------------------------
medskipnoindentResidual and Root certificates.

% ---------------------------------------------------------------

For a real phase matrix P define F(P)inC^{15} by

 F_{ij}(P)=frac 1{2pi i}sum_{g=0}^5e(P_{ig}-P_{jg}),
  0le j<ile 5,

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}toR^{30} in complex
notation is

 (L_Px)_{ij}=sum_{g=0}^5e(P_{ig}-P_{jg})(x_{ig}-x_{jg}).

The factor 2pi i from differentiation cancels the denominator
in~(phase:residual). These are the residual and full derivative in
buildH::gen.

[Root certificate]
A Root consists of a stored real center P, a stored real
36times k matrix T, and nonnegative constants
(C), n0, id, ort, hor,
mnorm, froot. The source stores column a of T
as T[a] and uses 0le kle 10. Write

 n_P(s)=n0+id,s+C,s^2, sge 0 .

The data are certified if the basic bounds

 |T^*T-I|lemathtt{ort},
 |(I-Pi)T|lemathtt{hor},
 |F(P)|lemathtt{froot}

hold and in addition, for every horizontal delta with F(P+delta)=0
and eta=T^*delta, s=|delta|,

|(I-TT^*)delta|&le n_P(s), 
s^2&le(1+ort)|eta|^2+
                  |(I-TT^*)delta|^2, 

and, for every such delta and every horizontal v, with
h'=delta-v,

 |(I-TT^*)h'|
 le[id+C(|delta|+|v|)]|h'|
       +mnorm,|F(P+v)| .

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.


For a displacement delta of an exact solution, given s_0ge|delta|
and t_0ge|T^*delta|, nonnegativity of the coefficients of n_P
gives
[
 |delta|le
 (n_P(s_0)^2+(1+ort)t_0^2)^{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.


The columns of T are numerical directions only. Neither exact
horizontality nor exact orthogonality is assumed; in particular
|T|lesqrt{1+ort}. The 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.


% ---------------------------------------------------------------
medskipnoindentConstructing Root certificates.
% ---------------------------------------------------------------

The construction uses an auxiliary real map M:R^{30}toR^{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.

[Quadratic and difference remainders]

There is a real 30-vector R_P(x) with
F(P+x)=F(P)+L_Px+R_P(x). Put
a_{ij}(x)=|x_i-x_j|^2 (where x_i is row i of x). Then

|(R_P(x))_{ij}|&lepi a_{ij}(x), 
|(R_P(x)-calR_P(x)-R_P(v))_{ij}|
&lepisqrt{a_{ij}(x-v)}
      (sqrt{a_{ij}(x)}+sqrt{a_{ij}(v)}).


For a real unit yinR^{36} let w_{ij}(y) be the Euclidean norm of
the two edge coordinates of M^*y. If

 sum_{i>j}w_{ij}(y)a_{ij}(x)le b|x|^2
 quadfor every horizontal x and every unit y,

then for horizontal x,v,

 |MR_P(x)|&lepi b|x|^2,
 
 |M(R_P(x)-calR_P(x)-R_P(v))|
 &lepi b,|x-v|(|x|+|v|).
 



Proof. 
For real t put r(t)=(e(t)-1-2pi i t)/(2pi i). Its derivative is
e(t)-1, of modulus at most 2pi|t|, so |r(t)|lepi t^2 by
integrating from zero. On the segment from v to x, convexity of
absolute value gives
[
 |r(x)-r(v)|
 le 2pi|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~(phase:edge-difference).

Pairing MR_P(x) with a real unit y and using the
two-coordinate Cauchy--Schwarz inequality on each edge
proves~(phase:matrix-remainder). For the difference, apply
separately
[
 sum_{i>j}w_{ij}
     sqrt{a_{ij}(x-v)a_{ij}(x)}
 le
 (sum_{i>j}w_{ij}a_{ij}(x-v))^{1/2}
 (sum_{i>j}w_{ij}a_{ij}(x))^{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.


[Two computable remainder constants]

Two sufficient choices in~(phase:weighted-constant) are
b=sqrt 8,|M| and b=sqrt 6,|S|, where the
30times 216 real matrix S has row
[
 S_{a,:}
   =M_{:,a}^* ⊗ ({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.


Proof. 
For horizontal x with row vectors summing to zero, let s=|x| and
u_i=|x_i|^2. Expansion gives

 sum_{i>j}a_{ij}(x)&=6s^2,nonumber
 sum_{i>j}a_{ij}(x)^2
   &=6sum_i u_i^2+s^4+2sum_{i,j}(x_i· x_j)^2 .
 

The cross terms use sum_{jne i}x_i· x_j=-u_i. Also
u_ile 5(s^2-u_i) and sum_{i,j}(x_i· x_j)^2le s^4. Thus
sum_i u_i^2le (5/6)s^4, and the second sum
in~(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 S to
y ⊗  x_{:,g} and sum squared norms over g: for unit y this
gives sum w_{ij}(y)^2a_{ij}(x)le|S|^2s^2. Another
Cauchy--Schwarz, now using sum a_{ij}=6s^2, proves the claim. The row
Gram of S is computable without forming its 216 columns:

 (S S^*)_calS^*)_{ab}
  =(({bf e}_i-{bf e}_j)·({bf e}_{i'}-{bf e}_{j'}))
           (M_{:,a}· M_{:,b}),

with the edge labels repeated for real and imaginary parts. This is the
sign-factor formula in buildH::patch.


[Sufficient conditions for a Root certificate]

Let P,T and the nonnegative constants be as in
Definition~phase:root-data. Suppose~(phase:root-basic) holds
and an auxiliary stored real 36times 30 matrix M satisfies

 |M|lemathtt{mnorm},
 mnorm,frootlemathtt{n0},
 |Pi-TT^*-ML_PPi|lemathtt{id},

and C is an upper bound for a remainder constant pi b in
Lemma~phase:remainder. Then these data form a certified Root.


Proof. 
Let delta be horizontal with F(P+delta)=0 and put eta=T^*delta.
Write E=Pi-TT^*-ML_PPi. Taylor's formula at the exact solution gives
[
 (I-TT^*)delta=Edelta-MF(P)-MR_P(
]
The certificate and Lemma~phase:remainder
prove~(phase:root-normal). For the squared inequality (no
projector assumption on TT^*),
[
 |delta-Teta|^2
   =|delta|^2-2|eta|^2+eta^*(T^*T)eta
   ge|delta|^2-(1+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(R_P(-R_P(v)).
]
Apply~(phase:matrix-difference). This
proves~(phase:root-difference), including the residual term at the
possibly inexact P+v.


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.

[Certificates from the direct checks]

A successful call to buildH::patch under the arithmetic
conditions of the certified arithmetic
certifies~(phase:root-basic),
(phase:projector-certificate), and a valid remainder constant C
in the sense of Definition~phase:root-data.


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 mnorm also includes the
difference between the accumulated stored M and the exact product
H_VL_V^*. The required additions and the 2.1· 10^{-11}
Frobenius bound on this difference are proved in the certified
arithmetic. Thus mnorm bounds the actual stored map, and the
first choice in Lemma~phase:remainder-constants
supplies a certified C.

For the optional tensor improvement the row Gram
(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}|^2le 1/(1-beta). A certified absolute row
sum bound r for U^*(S S^*)UcalS^*)U gives
|S|^2le r/(1-beta). This uses congruence by an
invertible matrix, not assumed orthogonal diagonalization. The second
choice in Lemma~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~(phase:projector-certificate) is evaluated, not just its
proposed horizontal form. The padded Frobenius residual gives
id. Corresponding checks on T^*T-I, (I-Pi)T and F(P)
give the basic bounds, and the stored n0 is an upper bound
for mnorm,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~phase:root-bounds applies to
every successful call.


% ---------------------------------------------------------------
medskipnoindentTransport and alignment.
% ---------------------------------------------------------------

Callbacks in variants select a real 6times 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=PhiPhi^*=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,

 F(V)=F(Phi P), L_{Phi P}(Phidelta)=L_Pdelta,
  F(V+Phidelta)=F(P+delta)=0.

Hence a Root transported by Phi remains a valid Root certificate with
the same constants; the source writes such a transport through
Root::transport. This is the map P|->Phi P used below.

After transport, replace P by its dephasing D_0(P). Dephasing is a
row and column rephasing; by Lemma~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/(2pi)+ε_X with the further split of u
into a seed part s and an angular part d; then take a compact
sub-ball.

[Alignment]
For real centers P,Q and integer matrix N,

 d=Pi(D_0(P)-D_0(Q)-N)
 ⇒
 e(Q+d) is phase-equivalent to e(P).

Under that alignment every exact solution displacement delta in a
completion cell takes the form

 delta=Delta+Pi(ε+frac u{2pi}),
 Delta=Pi(P_0-D_0(P)-N),

with the same equivalence of delta-parameterized matrices.


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 kerPi 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~seed:equivalences.


The standard dephase and 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.

% ---------------------------------------------------------------
medskipnoindentDisplacement 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, 1le ale 4, and
Z=D'-D has |Z|_Fle z and |Z-J(s,theta)|_Fle E_r as
in . Let I_s=(7,8,9,10,11,13,19,25,31)calI_s=(7,8,9,10,11,13,19,25,31)
denote the nine seed indices, and I_D the sixteen lower-right
indices (i,j) with 2le i,jle 5. The first row and column have zero
raw displacement.

For a retained cell we want a bound on the exact horizontal delta
in~(phase:delta). Let beta be a fixed parameter; the actual
codes use beta=.45 or beta=.8.

[Arc and weighted-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_jin[-pi,pi] from D_j^circ to D'_j. Let slack be a
nonnegative upper bound for z^2-sum_j(m_j-1)^2, let mi be
a positive lower bound for min_jm_j, and choose
q_dge slack/mi, q_d<3.8. Then valid upper bounds
are

 |d|lemathtt{arc},
 (sum_jm_jd_j^2)^{1/2}lemathtt{weighted},

provided
arcge√(q_d)/√(1-q_d/4) and
weightedgesqrt{slack}/√(1-q_d/4).


Proof. 
Since D'_j has modulus one,

 |D'_j-D_j|^2=(m_j-1)^2+m_j|D'_j-D_j^circ|^2 .

Summing proves sum_jm_j|D'_j-D_j^circ|^2lemathtt{slack}, and
hence sum_j|D'_j-D_j^circ|^2le q_d. In particular no individual
chord reaches 2. If c_j=|D'_j-D_j^circ|,
|d_j|=2arcsin(c_j/2). On 0le cle√(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.


The source's Candidate::prep only attempts this argument when
its direction and magnitude checks pass. It uses the stored flat
sum and sets

 slack=max(0,z^2-sum+10^{-5})+10^{-10},
 q_d=slack/mi+10^{-8}.

The subtraction 10^{-6} from each computed magnitude gives mi
and the upper reciprocal
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~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.

[Stored center and support representation]

Let P_0 be the stored Candidate phase matrix. On a successful
preparation path there is a correction ε, in cycles, with

 |ε|le 2· 10^{-6},

such that, up to entrywise integer choices, the actual feasible phase
matrix has the form

 P_0+ε+frac u{2pi}.

Here the raw angular matrix u is s at I_s, is d at
I_D, and is zero elsewhere, so in particular

 |u|lesqrt{R_s^2+arc^2}.

Let X:R^{13}toR^{36} give s at the seed positions and, at the
lower-right positions, give

 (X(s,theta))_j
       =Im(D_j^circbarJ(s,theta)_j).

Then

 u=X(s,theta)+e_D+e_{sin},
 |e_D|le E_r,
 |e_{sin}|_1lemathtt{arc}^3/6,

where both errors are supported on I_D.


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~completion:linear-remainder,
[
 D'_j=D_j+J(s,theta)_j+R_j, |R|le E_r,
]
and Lemma~phase:arc-bounds bounds |d_j|lemathtt{arc}. Divide
D'_j by m_j: the real part is cos d_j, the imaginary part is
sin d_j. Thus
[
 Im(D_j^circbarJ(s,theta)_j)
 =m_j^{-1}Im(D_j^circbar(D'_j-D_j))
 +O(E_r)= d_j+O(E_r)+O(arc^3).
]
The O(E_r) contribution is at most E_r per component and the
O(arc^3) sine expansion error is at most 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~(phase:raw-displacement).
The stored center P_0 differs from an exact solution by
ε+u/(2pi) up to integer entries. Dephasing an exact phase matrix
changes P_0 by less than 2· 10^{-6} cycles per entry, as proved
in the certified arithmetic.


[A direct displacement bound]

Either of the following is an upper bound for Delta· u:

 S_1={}&R_s|(X^*Delta)_s|
          +h|(X^*Delta)_theta|_1
          +E_r|Delta_D|
          +|Delta_D|_∞arc^3/6, 
 S_2={}&R_s|Delta_s|
          +weighted
                  (sum_{jinI_D}/m_j)^{1/2}.
 

With S_Delta=min(S_1,S_2), a valid upper bound is

 |delta|le
 (|Delta|^2+
       frac{R_s^2+arc^2}{(2pi)^2}
       +(S_Delta)/(pi))^{1/2}
       +|ε|.



Proof. 
Equation~(phase:raw-support), the seed ball, and the angle box
give~(phase:support-linear). Direct weighted Cauchy--Schwarz
using~(phase:arc-weighted) gives~(phase:support-weighted).
Since Delta is horizontal, Delta·Pi u=Delta· u.
Expand the squared norm of Delta+Pi u/(2pi), contract by Pi, and
add |Piε|le|ε| by the triangle inequality.


% ---------------------------------------------------------------
medskipnoindentEllipsoid 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,
2pi T^*Delta+hT^*X_theta[-1,1]^4 is the convex hull of the 16
offsets

 U_sigma=2pi T^*Delta+
                hsum_{a=1}^{4}sigma_a T^*X_{theta,a},
  sigma_ain{-1,1}.


[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

 Gsucceq(1+beta)r^2AA^*+(1+1/beta)q^2BB^*.



Proof. 
For each real test vector z, the support is at most
r|A^*z|+q|B^*z|, which by the inequality
(x+y)^2le(1+beta)x^2+(1+1/beta)y^2 is at
most (z^*Gz)^{1/2}. Diagonalize G: the vectors satisfying
|z· x|le(z^*Gz)^{1/2} for every z are exactly G^{1/2} times
the unit ball, including when G is singular.


By Lemma~phase:ellipsoid-containment with beta=.45, the sum of
the seed and e_D contributions to T^*u lies in an ellipsoid whenever

 G_{rm lin}succeq
    1.45R_s^2(T^*X_s)(T^*X_s)^*
       +(1+1/.45)E_r^2T_D^*T_D.

The remaining sine error has norm at most
(max_{jinI_D}|T_{j,:}|)arc^3/6.
There is also a bound using the raw seed and weighted angle balls with a
single offset U_0=2pi T^*Delta:

 G_{rm dir}succeq
    1.8R_s^2T_s^*T_s+
       (1+1/.8)weighted^2
                     T_D^*diag(1/m_j)T_D.


[Offset ellipsoid norm bound]

Let Gsucceq 0 and UinR^k. For every
gamma>lambda_{max}(G) and |w|le 1,

 |U+G^{1/2}w|^2
   legamma(1+U^*(gamma I-G)^{-1}U).

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.


Proof. 
Put A=gamma I-Gsucc 0. Expanding and completing the square,
[
 |U+G^{1/2}w|^2
 legamma+|U|^2+
       2w^*G^{1/2}U-w^*Aw
 legamma+|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.


Let E(G;{U_sigma}) denote any certified upper bound from
Lemma~phase:ellipsoid-norm. Two valid tangent bounds are

 |T^*delta|
 le{}&
 frac{E(G_{rm lin};{U_sigma})+
      max_{jinI_D}|T_{j,:}|,arc^3/6+
      horsqrt{R_s^2+arc^2}}{2pi}
       +sqrt{1+ort},|ε|,
 
 |T^*delta|
 le{}&
 frac{E(G_{rm dir};{U_0})+
      horsqrt{R_s^2+arc^2}}{2pi}
       +sqrt{1+ort},|ε|.
 


Numerical realization. In covercompute,
off is the computed 2pi T^*Delta; AB is the
computed T^*X; U holds the 16 corners. The covariance uses
the stored AB coefficients and charges their ideal-coefficient
error separately outside the ellipsoid.  gives
|widehatwidehat A-T^*XA-T^*X|<1.5· 10^{-5},
|(s,theta)|<.226, so this substitution costs less than
3.4· 10^{-6} radians. For both covariance choices a diagonal
10^{-7}I is added. The inverses in the second choice are certified
upper invmag values.  proves PSD
domination of the exact covariance expressions in the stored-coefficient
model.

% ---------------------------------------------------------------
medskipnoindentBootstrap and local coverage.
% ---------------------------------------------------------------

[Finite improvement of a certified radius]

Suppose s_0ge|delta| and t_0ge|T^*delta| for every feasible
solution in a completion cell, using the particular
alignment~(phase:delta). Starting from any such upper s, a
further upper is

 (n_P(s)^2+(1+ort)t_0^2)^{1/2}.

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.


Proof. 
Lemma~phase:root-data, nonnegativity of the coefficients, and
monotonicity of n_P give
|delta|^2le n_P(|delta|)^2+(1+ort)|T^*delta|^2
le n_P(s)^2+(1+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.


% ---------------------------------------------------------------
medskipnoindentAtlas merge of certified balls.
% ---------------------------------------------------------------

For a recorded phase ball, the Atlas stores the Root (P,T), the
radius r, and the certified constants. Two entries are 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.

[Atlas merge]
Suppose two Root records (P,r_P), (Q,r_Q) satisfy
[
 |(P+delta_P)-(Q+delta_Q)|le d_{PQ}
 quadwhenever 
 |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
S(P,r_P)calS(P,r_P)cupS(Q,r_Q)calS(Q,r_Q) is contained in
[
 S(P, max{r_P, d_{PQ}+r_Q+allowance}).
]


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 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
S in Definition~refine:solution-set is closed under this
equivalence, so it contains the aligned representative. Inflate to
r_Pvee d_{PQ}+r_Q and add the outward allowance.


The implementation uses mapH for Phi.
Its two index arrays are genuine permutations. In variants, the
leading two row indices are distinct and 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,col[j]+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, 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.

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. Atlas::add grows the selected
radius to at least the incoming radius plus this distance, with another
10^{-8} added and a separate 2· 10^{-8}-margin comparison to the
cap. It remains valid when an alternative map found by 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.

[Output of a successful A stage]

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.


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~phase:atlas-merge). Boundaries cause no exception:
overlapping assignments are harmless.


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· 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 topK=4 in an array of
size~8.

% =====================================================================
medskipnoindentRefining the phase-ball cover.

% =====================================================================

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; V=im is the
horizontal subspace.

For a certified center P and radius Rge 0 put

 S(P,R)=calS(P,R)=
 bigl{[e(P+delta)]:deltainV, |delta|le R,
       F(P+delta)=0bigr},

where brackets denote the equivalences of
Lemma~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.

% ---------------------------------------------------------------
medskipnoindentCertified relations between centers.
% ---------------------------------------------------------------

For a Root at P denote its stored tangent matrix T:R^ktoR^{36}
and put
 
 o_P=ort_P, alpha_P=id_P, gamma_P=mathtt C_P,
 mu_P=mnorm_P, ε_P=froot_P,
 n_P(s)=n0_P+alpha_Ps+gamma_Ps^2.

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}toR^{36}. In particular

 |T^*T-I|le o_P,
 |S^*S-I|le o_Q,
 |F(Q)|leε_Q.



Put A=I-TT^*. For every xinR^{36}, with zeta=T^*x,

 |T|le√(1+o_P),
 |A|le 1+o_P,
 |x|^2le(1+o_P)|zeta|^2+|Ax|^2.

These statements do not require x or the columns of T to be
horizontal.


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+|Tzeta|^2
ge|x|^2-(1+o_P)|zeta|^2, using
|Tzeta|^2ge(1-o_P)|zeta|^2.


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

 d_q=Pi(D_0Q-D_0P-N),
 eta_q=T^*d_q, n_q=|A d_q| .

Because D_0X-XinkerPi, Q-P-d_qin N+kerPi; vectors in kerPi
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
deltainV with |delta|le R. Set

 eta=T^*delta,
 h'=delta-d_q, d=|h'|
 ,
 xi=T^*h'=eta-eta_q ,
 nu=Ah'.

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.

[Relation estimates]
The following three basic bounds hold.

 d&lemin{R+|d_q|, √((1+o_P)R^2+n_q^2)},
 
 |nu|&lemin{nu_0,(1+o_P)d,lambda d+e_0},
 
 |S^*h'|&lebeta t+sigma|nu|.
 

Here

 nu_0=n_P(R)+n_q,
 lambda=alpha_P+gamma_P(R+|d_q|),
 e_0=mu_Pε_Q,


 beta=√((1+o_P)(1+o_Q)),
 sigma=|S_h|,(1+o_P)^{1/2},
 t=|eta-eta_q|,

and S_h=Pi S; and for the next step the child Root's normal and
tangent bounds apply to h':

 d^2le n_Q(d)^2+(1+o_Q)|S^*h'|^2,
  n_Q(d)=n0_Q+alpha_Qd+gamma_Qd^2.




Proof. 
The triangle inequality gives dle R+|d_q|. The identity
|x|^2le(1+o_P)|T^*x|^2+|(I-TT^*)x|^2 from
Lemma~refine:approximate-projection, applied to h', gives the
second inequality in~(refine:parent-distance).

For~(refine:normal-three): nu=A(delta-d_q), so
|nu|le|Adelta|+|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
(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ε_Q, giving
|Ah'|lelambda d+e_0.

For~(refine:child-tangent), first use horizontality of h' to
write S^*h'=S_h^*h'. The identities h'=Txi+nu and
T^*nu=(I-T^*T)xi yield
[
 S^*h'=S_h^*Txi+S_h^*Anu+S_h^*T(I-T^*T)xi.
]
Put beta:=√((1+o_Q)(1+o_P))ge|S_h^*T| and
sigma:=|S_h|√(1+o_P)ge|S_h|,|A|, using
|A|le 1le√(1+o_P). Then
|S_h^*Txi|lebeta|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 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~(refine:child-distance), apply the squared inequality at the
child Root to the horizontal exact-solution displacement h' of Q:
d^2=|h'|^2le n_Q(d)^2+(1+o_Q)|S^*h'|^2.


% ---------------------------------------------------------------
medskipnoindentUniform iteration without contraction.
% ---------------------------------------------------------------

The relation estimates yield a useful upper iteration without requiring
either contraction or a fixed point.


 
 z&=√(1+o_P),t,
 v(b)&=min{v_0, (1+o_P)b, lambda b+e_0},
 psi_P(b)&=√(z^2+v(b)^2),
 psi_Q(b)&=
 √((1+o_Q)(beta t+sigma v(b))^2+n_Q(b)^2).
 

Start at

 b_0=min{D, √(z^2+v_0^2)},

and replace b_j by

 b_{j+1}=min{b_j, psi_P(b_j), psi_Q(b_j)}.



Every iterate in~(refine:iteration-start)--(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.


Proof. 
Both candidates in~(refine:iteration-start) are upper bounds by
Proposition~refine:relation-estimates. If dle b_j,
(refine:normal-three) and nonnegativity of all coefficients give
|nu|le v(b_j). Thus~(refine:parent-distance) gives
dlepsi_P(b_j); and~(refine:child-tangent)--(refine:child-distance),
together with n_Q(d)le n_Q(b_j), give dlepsi_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.


% ---------------------------------------------------------------
medskipnoindentThe 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,ε_Q. For every tangent threshold
tge 0, let widehatwidehat D(t)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
(refine:iteration-step) at
b_0(t)=min{D, √((1+o_P)t^2+v_0^2)}; by
Lemma~refine:iteration-certificate, the final iterate is a valid
upper. The source stores it back through the
alignment~(refine:alignment).


 |eta-eta_q|le t  and  widehatwidehat D(t)D(t)le a
  ⇒
 [e(P+delta)]inS(Q,a)calS(Q,a)
 quadfor every parent solution delta.


For a queried sphere radius 0le ale .126, the source's init
searches for a usable tangent distance rho as follows. Set ell=0,
u=1.01a. If widehatwidehat D(0)>a-10^D(0)>a-10^{-7}, report no relation. Otherwise
make 25 midpoint trials, each with m=(ell+u)/2: if
widehatwidehat D(m)<a-2D(m)<a-2· 10^{-7}, replace ell by m; otherwise
replace u by m. Finally set

 rho=max(0,ell-10^{-8}),

and report a usable relation only if rho>10^{-7}.


Every usable returned threshold certifies
[
 |eta-eta_q|lerho
 ⇒
 [e(P+delta)]inS(Q,a)calS(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 widehatwidehat D.D.


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 widehatwidehat D(D(ell)<a-2· 10^{-7}, and
rho=ell-10^{-8}<ell. Take any parent solution delta with
|eta-eta_q|lerholeell. By construction widehatwidehat D(D(ell) is
an upper bound for |h'| at this solution, so
|h'|<a-2· 10^{-7}<a. Apply~(refine:eval-semantics). No
comparison between widehatwidehat DD at distinct arguments is needed.


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 widehatwidehat DD 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.

% ---------------------------------------------------------------
medskipnoindentExhaustive grid in tangent coordinates.
% ---------------------------------------------------------------

For every parent solution in the parent record,

 |eta|=|T^*delta|le√(1+o_P),R=:E_P.

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~(refine:tangent-range),
and chooses
[
 L=leftlceil E/w+0.500001rightrceil, 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 uinR^k and hge 0 write

 B(u,h)=prod_{a=1}^k[u_a-h,u_a+h].

The base boxes have u_a=b_aw, -Lle b_ale L, and h=w/2. Their
union contains [-E_P,E_P]^k, so every tangent vector
in~(refine:tangent-range) belongs to at least one base box. All
intervals are closed; grid boundaries are included.

medskipnoindentBlocker induction.

A 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.


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 A be the
final atlas of children and B the list of input parents
whose refinement failed. Then

 
 ∪_{(P,R),input}S(P,R)calS(P,R)
 &⊆
 ∪_{(Q,r)inA} S(Q,r)calS(Q,r)
 & cup
 ∪_{(P,R)inB} S(P,R).calS(P,R).
 

In particular, a completed run with no bad parents replaces the input
cover by a child cover. A nonempty bad list remains an obligation.


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~phase:arc-bounds
and~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~(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 B. Any children produced while attempting it
remain valid stored spheres and are reused without assuming the failed
parent has been covered.


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.

medskipnoindentBatch 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

 r_iletau,
 S(Q_i,r_i)⊆calS(Q_i,r_i)⊆S(Q_i,tau),

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~(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 0le Rle 0.126, and the child cap and width satisfy
[
 0.0001letaule 0.126, 0.0001le wle .05, 0le kle 10.
]
Normal constants are checked against
[
 gamma_P<20,
 n0_P,ε_P<10^{-6},
 alpha_P<10^{-4},
 o_P,hor_P<10^{-5},
 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 Lle 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.

% =====================================================================
medskipnoindentExcluding two further mutually unbiased bases.

% =====================================================================

The final stage takes a certified Root (P,T) with 0le rle .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(2pi i t). A real matrix in R^{36} is
identified with row-major indexing. The Root constants omega,n_0,
iota,Gamma denote ort,n0,id,C
as in the arithmetic bounds above. Its certified data give

 |T^*T-I|leomega,
 |F(P)|le f_0<10^{-6}, 0leomega<10^{-5},

and, for every horizontal delta with |delta|le r such that
e(P+delta) is Hadamard,

 eta=T^*delta, nu=delta-Teta,
 |eta|le r_H:=√(1+omega),r,
 |nu|le n:=n_0+iota r+Gamma r^2.


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,...,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=2pisqrt2,|F(P)| and~(exclude:root-data) give

 |A/sqrt6|_{op}
 lesqrt{1+frac{2pisqrt2}6f_0}<K, K:=1.00001,

where A=e(P).

% ---------------------------------------------------------------
medskipnoindentClosed phase bins.
% ---------------------------------------------------------------

At resolution M the five variable phases are covered by closed
intervals [b/M,(b+1)/M], 0le 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(sum_j d_j){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)
ith entry of the row over inner products. Set

 c_i=(|S_i|^2-1)/(4pi),
 J_{ij}=-Im(S_ibarq_{ji}),
 B_{ik}=sum_j J_{ij}T_{(j,i),k},

where q_{ji}=e(P_{ij}-P_{i0}) (so sum_jq_{ji}=S_i; and note
sum_jJ_{ij}=0). For yinR^6 put
[
 G_y:=max_j|sum_i y_iS_ibarq_{ji}|,
 
 N_y:=(sum_{i,j}y_i^2J_{ij}^2)^{1/2}.
]

[Direct displacement tests]

If |x|le X, |x|_1le L, and |delta|le r, then

 |S'-S|le 2pi(KX+r),
 
 |S'_i-S_i|le 2pi(L/sqrt6+r).
 

Consequently unbiasedness requires

 max{0,1-tau}^2le |S_i|^2le(1+tau)^2,
 
 sum_i(|S_i|-1)^2le[2pi(KX+r)]^2,

with tau=2pi(L/sqrt6+r).


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 2pi KX. For the second, the ith component has
magnitude at most (2pi/sqrt6)sum_j|delta_{ji}|
le 2pi|delta_{:,i}|, and their combined norm is at most 2pi r.
The scalar bound is the same with the first change at most
2pi L/sqrt6. Apply ||S_i|-|S'_i||le|S_i-S'_i| and use
|S'_i|=1.


[Weighted Taylor estimate]

Suppose |x|le X and |x|_∞le I. For every yinR^6,
unbiasedness implies

 |y·(c+Jx-Beta)|leE(y;X,I),calE(y;X,I),

where

 E(y;X,I):=calE(y;X,I):={}&nN_y+piBigl[
   |y|_∞(KX+r)^2
   +min{G_yX^2, K|Scirc y|IX}notag
 &hspace{28mm}
   +frac 2{sqrt6}|Scirc y|Xr
   +frac 1{sqrt6}|Scirc y|_∞ r^2Bigr].
 

For a bin, this gives the necessary support inequality

 |y· c|
 le asum_{j>0}|(J^*y)_j|+r_H|B^*y|
       +E(y;X_M,I_M).calE(y;X_M,I_M).




Proof. 
Write Delta_i=S'_i-S_i. The expansion of
(|S'_i|^2-1)/(4pi) has constant term c_i, linear term
sum_jJ_{ij}(x_j-delta_{ji}), and the squared increment
|Delta_i|^2/(4pi). After multiplication by y_i and summation the
latter is bounded by pi|y|_∞(KX+r)^2 via
(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|Scirc y|IX, and the minimum is valid. The
Hadamard-delta linear term is Beta; its remainder with nu is
bounded by nN_y. The cross term and pure delta term give the
remaining |Scirc y|Xr and |Scirc y|_∞ r^2 contributions.
Subtracting the bin's linear box term yields the necessary support
bound~(exclude:bin-support).


The bin routine Layer::test applies~(exclude:direct-tests),
then the scalar test~(exclude:scalar-linear)
(see below), and finally~(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

 |c_i|le asum_{j>0}|J_{ij}|+r|J_{i,:}|+rho,

where

 rho=minBigl{&
 pi[(KX+r)^2+
   S_{max}(KIX+(2Xr+r^2)/sqrt6)],notag
 &pi[(L/sqrt6+r)^2+
   S_{max}(X+r)^2/sqrt6]Bigr}.
 

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~(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.

% ---------------------------------------------------------------
medskipnoindentFrom 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 12to 24to 48to 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 vinR^6 and sge 0, define the closed vector ball

 B_{vec}(v,s)=
 Bigl{(lambda e(v+x))/(sqrt6):
       |lambda|=1, xinR^6, sum_jx_j=0, |x|le sBigr}.


[Certified clustering]

Suppose each retained bin midpoint u is assigned to a center v and
an integer vector minZ^6. Let
[
 h=u-v-m-frac16(sum_j(u_j-v_j-m_j)){bf 1}.
]
If the assigned radius is at least |h|+X_M, then
B_{vec}(v,s) covers the whole bin up to global
phase. The final assignment in cluster has this property.


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.


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 0le 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~exclude:cluster-cover. In particular, the argument does not
assume that the selected wrapped lift minimizes torus distance.

[Distinct labels]

Two vectors in the same ball B_{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.


Proof. 
After removing global phases, their unscaled inner-product sum has the
form sum_je(z_j) with |z|le 2s<.2. Therefore
[
 |sum_je(z_j)|
 ge 6-sum_j|e(z_j)-1|
 ge 6-2pisqrt6,|z|
 >6-frac25pisqrt6>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.



% ---------------------------------------------------------------
medskipnoindentBin geometry and edge tests.
% ---------------------------------------------------------------

[Bin geometry]

Fix a bin with midpoint uinR^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 jge1. 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 ,

 |x|le X_M, |x|_1le L_M, |x|_∞le I_M.



Proof. 
sum_jx_j=0 is immediate from the definition.
Write d=barbar d+(d-d+(d-barbar d)d) with barbar d=(1/6)d=(1/6)sum_jd_j. Then
|x|^2=|d-barbar dd{bf 1}|^2=|d|^2-6barbar d^2d^2le|d|^2le5a^2,
so X_Mgesqrt5,a suffices. Each |d_j-barbar d|d|le 2a, so
L_Mge10a and I_Mge2a also suffice. 
accounts for the 10^{-12} margin.


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

 Q_j=e(u_j-v_j), S=sum_jQ_j,
 z=x_u-x_v, r_e=s_u+s_v, |z|le r_e.

The actual unscaled inner product is, up to conjugation and global
phase, S'=sum_jQ_je(z_j). Since |S'-S|le 2pisqrt6,r_e,
necessary direct tests for orthogonality (O) and
unbiasedness (U) are

 (|S|)/(2pisqrt6)le r_e (O),
 (||S|/sqrt6-1|)/(2pi)le r_e (U).



 (D_{O}{:,j}=
 -ImQ_j ReQ_j,
 
 (D_{U}
 frac{ReS,(D_{O}{1j}
       +ImS,(D_{O}{2j}}6.


[Edge Taylor estimates]

For any actual pair from two distinct closed vector balls, with the
data of~(exclude:edge-data), and for
S_{O},
S_{U}|S|^2-6)/(4pi),

 |Re(S'-S)-(D_{O}|
 &le E_{O}:=4pisqrt6,r_e^2, 
 |(|S'|^2-6)/(4pi)-S_{U}-D_{U}|
 &le E_{U}:=
   (2pimax(|S|,sqrt6))/(sqrt6),r_e^2 +4pi r_e^2.
 

Therefore the necessary support conditions

 |S_{O}|
 &le r_e|D_{O}|_F+E_{O},
 
 |S_{U}|
 &le r_e|D_{U}|+E_{U}
 

hold for every actual pair.


Proof. 
Write e(z_j)-1=2pi iz_j+ε_j with
|ε_j|le2pi^2z_j^2 by the integral Taylor remainder.
Multiply by Q_j=a_j+ib_j and separate real and imaginary parts:
Re(Q_j(e(z_j)-1))approx-2pi b_j z_j=(D_{O}
and Im(Q_j(e(z_j)-1))approx2pi a_jz_j up to a
residual of modulus at most 4pi^2 r_e^2 per coordinate. Summing over
j and using |z|le r_e proves~(exclude:orth-residual).

For unbiasedness, expand (|S'|^2-6)/(4pi) at S. The constant is
S_{U}, the linear term is D_{U} (the division by
4pi normalizes the derivative exactly), and the quadratic remainder
is at most (2pi|S|/sqrt6)r_e^2+4pi^2r_e^2; substituting the
residual of the orthogonality calculation
gives~(exclude:unbiased-residual). Cauchy--Schwarz on
|z|le r_e gives the support inequalities.



% ---------------------------------------------------------------
medskipnoindentThe joint test with one shared tangent variable.
% ---------------------------------------------------------------

For each ball t=0,1, form c^t,J^t,B^t,S^t from
(exclude:single-sums) at its actual center. Write its radius as s_t and put

 E_t^{ball}(y):=
 n(sum_{i,j}y_i^2(J^t_{ij})^2)^{1/2}
 +piBigl[
   |y|_∞(Ks_t+r)^2
   +s_t^2G_y^t
   +frac 2{sqrt6}|S^tcirc y|s_tr
   +frac 1{sqrt6}|S^tcirc y|_∞ r^2Bigr].

This is exactly the weighted estimate~(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 Ein{O,U} write ell=2 or
1 for the dimension of D_E, and let q_EinR^ell be the
constant part of the linearized edge equation:

 q_{O} := (ReS, ImS) inR^2,
 
 q_{U} := (|S|^2-6)/(4pi) inR^1,

with D_E and E_E as in~(exclude:edge-derivatives)
and~(exclude:orth-residual)--(exclude:unbiased-residual).

[Joint necessary support]

Suppose one vector from each of the two balls is unbiased to the
coordinate basis and to the same second basis, whose phase matrix
is e(P+delta)/sqrt6, and the two vectors have relation
Ein{O,U}. For every y_0,y_1inR^6 and
y_einR^ell,

 &|y_0· c^0+y_1· c^1+y_e· q_E|notag
 &quadle
 sum_{t=0}^1Bigl[
   s_tbigl|(J^t)^*y_t+(-1)^tD_E^*y_ebigr|
   +E_t^{ball}(y_t)Bigr]notag
 &+r_Hbigl|(B^0)^*y_0+(B^1)^*y_1bigr|
   +E_E|y_e|.
 



Proof. 
For each single vector the linearized equation is
c^t+J^tx_t-B^teta with weighted residual bounded by
E_t^{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.


This is exactly the comparison implemented by graph::joint.
It first calls single::support with a zero vector-displacement
coefficient array for each ball, and only then adds the edge
derivatives. Thus the maximum coefficient in~(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.

% ---------------------------------------------------------------
medskipnoindentExact graph exclusion.
% ---------------------------------------------------------------

Build one vertex for each retained ball. Join two distinct vertices by
an O edge iff their orthogonality type has not been excluded
by the direct and joint tests above; similarly for 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~exclude:cluster-cover and
Lemma~exclude:distinct-labels, two additional mutually unbiased
orthonormal bases would give disjoint vertex sets A,B satisfying

 |A|=|B|=6, A,B are O-cliques,calO-cliques,
 Atimes B⊆ E(U).
 

The following finite search tests this necessary configuration.

[Complete graph search]

The mode-1 search of graph::exist, started with both pools
equal to the full vertex set and first-clique target size six, returns
false only if~(exclude:graph-configuration) is impossible. A
work-cap exhaustion never gives a negative certificate.


Proof. 
First consider clique on a candidate set P and target size
m. For m=0 it succeeds, and for m=1 it tests nonemptiness. For
mge2 it repeatedly selects the least remaining vertex i, removes
it from P, and searches Pcap N_{O} 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~(exclude:graph-configuration).
Interchange their names if necessary so that
a=min(Acup B) belongs to A. At the root, 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 Acup 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 clique on the first-pool target
six, and for each resulting candidate A, on the second-pool target
six. When clique succeeds, its exact enumeration is used to
test Atimes A'⊆ E(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.


[Local exclusion]

Suppose
 the retained bins cover the closed phase ball of the Root in the
      sense of Lemma~exclude:cluster-cover;
 every retained bin has either been excluded by
      (exclude:direct-tests), (exclude:scalar-linear), or
      (exclude:bin-support), or been accepted and covered by the
      assigned ball;
 the O-calO- and U-graphscalU-graphs on the retained balls,
      after excluding any pair by the direct or joint tests conditional
      on eta, admit no
      configuration~(exclude:graph-configuration) as decided by a
      completed graph::exist search.

Then the phase ball contains no two additional mutually unbiased
orthonormal bases.


Proof. 
Suppose two such bases existed. Each of their twelve vectors is
covered by some retained ball (Lemma~exclude:cluster-cover).
Lemma~exclude:distinct-labels puts them in twelve distinct balls.
Within each basis the six vectors are pairwise orthogonal, and every
intra-pair survived its O direct and joint tests because the
tests are necessary conditions; thus they are joined by O
edges. Across the two bases every pair is unbiased, and similarly
they survive the U tests, so the cross-pairs are mathcal
UcalU edges. This produces a configuration~(exclude:graph-configuration),
contradicting the third hypothesis.


The conclusion of Proposition~exclude:local-certificate is
conditional on a completed negative result for the individual Record;
no collection of unexecuted tests is itself such a result.

% ---------------------------------------------------------------
medskipnoindentImplementation margins and invocation.

% ---------------------------------------------------------------

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

 widehatwidehat r_H=√(1+r_H=√(1+omega),r+10^{-9},
 widehatwidehat n=n_0+n=n_0+iota r+Gamma r^2+10^{-9}.

The upper squared-magnitude threshold in~(exclude:direct-tests) is
increased by 2·10^{-7}, and its lower threshold, when tau<1, is
decreased by that amount. The combined squared-distance threshold, the
right side of~(exclude:scalar-linear), and the right side
of~(exclude:bin-support) each receive a 2· 10^{-7}
comparison allowance.

For a pair, the radius sum s_u+s_v is increased by 10^{-10}. The
distance tests~(exclude:edge-distance) permit a further
2· 10^{-7} on that radius. Both edge residual bounds
E_{O} and E_{U} are increased by 2· 10^{-8};
the direct comparisons~(exclude:orth-support)
and~(exclude:unbiased-support) receive another 2· 10^{-7}.
The right side of the joint
comparison~(exclude:joint-inequality) receives 3· 10^{-7}.
Cluster assignment radii include the 10^{-9} allowance described in
Lemma~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 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
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~exclude:local-certificate is conditional on a
completed negative result for the individual Record; no collection of
unexecuted tests is itself such a result.


% =====================================================================
medskipnoindentFinite execution and the recorded certificate.

% =====================================================================

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.

% ---------------------------------------------------------------
medskipnoindentPrograms, records, and machine arithmetic.

% ---------------------------------------------------------------

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
FLT_EVAL_METHOD=0, 32-bit int, 64-bit
size_t, 64-bit uint64_t equal to
unsigned long, eight-bit bytes, and MXCSR 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.

medskipnoindentHow 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.


% ---------------------------------------------------------------
medskipnoindentControl flow and terminal conditions.

% ---------------------------------------------------------------

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:
 successful builds, a normally completed cover, and a complete
      center file;
 four A exits 0, all assigned centers processed, terminal
      failedBalls=0, and four empty bad text files;
 both B exits 0, all input parents processed, agreement of
      logged and saved output counts, and zero-count B bad files;
 every initial C batch processed exactly once with a consistent
      exit, done count, bad count, and saved record count;
 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.

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 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.


% ---------------------------------------------------------------
medskipnoindentThe recorded execution.

% ---------------------------------------------------------------

The complete run consists of 100 program invocations:
 one pass over the seed cover;
 four shards of the completion stage;
 two refinement passes, a mid one and a fine one;
 88 batches of the initial exclusion pass;
 one fallback refinement pass;
 four final exclusion batches.

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.


% ---------------------------------------------------------------
medskipnoindentThe finite-execution certificate.

% ---------------------------------------------------------------

[Finite-execution certificate]

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.


Proof. 
Suppose four such bases exist. Normalize the first to the standard
basis and write the second as H/sqrt6.

Seeds. Proposition~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.

Stages A and B. Proposition~phase:A-cover therefore places
an equivalent exact solution in a coarse Record. Each equivalence
preserves extendibility by two further bases
(Lemma~seed:equivalences). On that Record the refinement stage
replaces the ball by a finite atlas of children together with any
unresolved parent; Proposition~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.

Stage C. Proposition~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.


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.

% =====================================================================
medskipnoindentConclusion of the proof.

% =====================================================================

The proof of Theorem~thm:main is now complete. The
overall argument is the following chain of strictly finer covers, each
preserving every exact candidate solution.
 By Lemma~seed:equivalences, every complex Hadamard matrix,
      up to allowed symmetries, has its nine-phase seed inside one of
      finitely many closed balls (Proposition~seed:cover).
 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~phase:A-cover).
 Each Phase ball is subdivided into smaller children, and the
      original balls are covered by new children together with any
      unresolved bad parents (Proposition~refine:cover-preservation).
 Each small phase ball is ruled out by the C-stage graph test on
      independent vector balls and pairwise tests
      (Proposition~exclude:local-certificate).


The recorded execution of the pipeline
meets every terminal condition listed above, so the
finite-execution certificate Theorem~execute:finite-theorem
applies and yields the negative upper bound. The lower bound is
Lemma~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

% =====================================================================
medskipnoindentCertified arithmetic for the recorded execution.

% =====================================================================

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.

% ---------------------------------------------------------------
medskipnoindentPrimitive-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 -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 environment sets and checks nearest rounding and
checks the binary64 properties of 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 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 nearbyint and 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 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· 10^{-16} when compared with a displayed decimal.

Put

 u=2· 10^{-16}, tau=10^{-300},
 gamma_m=2.1· 10^{-16}(m+1).

For a correctly rounded basic operation with finite exact result z at
the scales below,

 |fl(z)-z|le u|z|+tau .

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.

[Accumulation and norm bounds]
For mle 2000, a dot product on exact stored inputs, formed as
separate products and successive additions, satisfies

 |fl!(sum_{i=1}^m a_i b_i)
                  -sum_{i=1}^m a_i b_i|
 le gamma_msum_{i=1}^m|a_i b_i|+4(m+1)tau .

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,

 |fl(|a|)-|a||
 le 2.2· 10^{-16}(m+2)|a|+10^{-140}.

If the computed components approximate other components, their
Euclidean difference is added to this bound.


Proof. 
Iterating |z_k|le (1+u)|z_{k-1}|+u+btau 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)^msum_{i=1}^m|a_ib_i| plus a O(mtau) term. The bound
gamma_m=2.1· 10^{-16}(m+1) is a conservative choice: for
mle 2000 one has (1+u)^m-1le mule 2m· 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.


% ---------------------------------------------------------------
medskipnoindentPi and exponential evaluations.

% ---------------------------------------------------------------

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· 10^{-15} is an ample
absolute error bound for each.

The source evaluates cosine through degree 18 and sine through degree
19. In the jth 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· 10^{-15} for jle 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_jle E_{j-1}/2+4tau, E_0=0. Thus E_j<8tau; 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· 10^{-21} and
1.25· 10^{-22}, respectively; both are below 10^{-16} even in
the worst scale considered here. Thus

 |fl(e^{iy}) - e^{iy}| le 10^{-14}

for all real y with |y|le 20 and no overflow or underflow.

For z=a+bi with |a|,|b|le20, the pair
(fl(cos b),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]; 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.


% ---------------------------------------------------------------
medskipnoindentDephasing and projection.

% ---------------------------------------------------------------

For |p_j|le 4, the three operations in unwrapped dephasing
p_{ij}-p_{i0}-p_{0j}+p_{00} have entry error below 8· 10^{-15};
their intermediate magnitudes are at most 8, 12, and 16.
Subtracting the integer selected by 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· 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· 10^{-15}, and its length-36 grand sum below
1.4· 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 Atlas::find's returned √(s)+10^{-7} is an upper
distance for the selected ideal projected lift. A merge uses radius
r+distance+10^{-8} and checks a cap with 2· 10^{-8}
slack. Its sum-rounding error is below 10^{-14}. The same arithmetic
applies to 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.

% ---------------------------------------------------------------
medskipnoindentSeed filtering and seed-ball reuse.

% ---------------------------------------------------------------

Here n=112 and the stored w approximates 1/n. Midpoint-gap
calculations contain only small integer arithmetic and a few
products/subtractions; (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
sumok has component error below 2· 10^{-13}. Indeed its
five exponential errors, after division by 2pi, total less than
4.8· 10^{-14} in complex modulus, and the divisions, additions,
and constant error cost less than 3· 10^{-15}. The sharper
componentwise bound 1.1· 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 chk cost together less than 10^{-11}.
Its remainder coefficient approximates 5pi/(4n^2) within 10^{-15}.
The actual 10^{-8} comparison addition therefore preserves every
feasible list.

In 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· 10^{-11} each. The signed term
consequently has input error at most
[
 2sum_{a=1}^4 |v_a|,3· 10^{-11}le 2.4· 10^{-6}<2.5· 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· 10^{-9}: its four input contributions are at most
4· 10^4ε_{exp}, and its short dot costs less than
10^{-10}. Each pair norm has error below 10^{-10} by
(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· 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
int8_t; all their squared/absolute sums are exact and small.
The comparison radius is
[
 RN2=fl((fl(.032·112))^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·112)^2. This reasoning requires the
pow(x,2) calls in 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 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.

% ---------------------------------------------------------------
medskipnoindentFrames, comparison blocks, and A exclusion tests.

% ---------------------------------------------------------------

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 ideal
projection score.

For frame, the four-term sum S has error below
2.5· 10^{-13}; its squared modulus has error below
2.2· 10^{-12}, and its modulus below 2.6· 10^{-13} by
(round:norm). If the stored m=√(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· 10^{-13} in m and
top=√(4+|S|). In particular 16-|S|^2>6.1; its error is
below 2.3· 10^{-12}, so its reciprocal has error below
7· 10^{-14}. The two pseudoinverse numerators have modulus at most
eight and error below 5.1· 10^{-13}. Division/scaling then gives

 |Delta G_{ij}|<10^{-12}, |G_{ij}^{rm ideal}|<.81,

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· 10^{-13} and denominator
√(4-|S|^2/4)>1.5 with error below 2· 10^{-13}. Its ideal
entries have modulus at most one, so the quotient error is below
4· 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· 10^{-13}, and its
norm error is below 2· 10^{-12}. The accepted computed norm is
greater than .449999, so the ideal projection is nonzero and
normalization gives entry error below 7· 10^{-12}. In the last
projection the extra column-product error gives entry error below
1.6· 10^{-11} and norm error below 3.5· 10^{-11}. The same
denominator check gives entry error below 1.2· 10^{-10}. This
argument inductively constructs the ideal orthonormal complements with
the selected coordinates. The two-product kernel projector has entry
error below 3· 10^{-10}; multiplying two complement entries and
sqrt6 gives a W entry error below 8· 10^{-10}.

The following table records the rest of the propagation, with the source
array names. Its errors include ordinary operation rounding.


small
{>{raggedrightarraybackslash}p{.22linewidth}p{.23linewidth}p{.45linewidth}}

Quantity & Absolute error per entry & Ideal scale used

base & 10^{-11} & two factors of modulus at most 2 times G
fixed & 2· 10^{-9} & projector times a factor at most 1.62
bm, cm & 2· 10^{-12} & unit entry times G
uc.U, angular derivatives & 2· 10^{-13} & products of unit phases/trig entries
K, D & 4· 10^{-9} & four W products and base
mag & 4· 10^{-8} & |D|<4.6
sqrt{mag} & 5· 10^{-9} & norm Lipschitz
jac & 10^{-8} & ideal modulus at most 6
flat real J & 1· 10^{-7} & dot of D and jac
flat c & 5· 10^{-8} & (|D|^2-1)/2
pre-equation c & 4· 10^{-13} & short unit-exponential sums
pre-equation coefficients & 6· 10^{-14} & unit exponentials




medskipnoindentThe horizontal-projector identity..
An entry of the exact double-centering projector is the product of two
centering-matrix entries. The computed 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· 10^{-13} from the edge products.
Forming LP gives error below 5· 10^{-12} per entry
against L_PPi; the ideal entry magnitude is at most 8.4. In the
residual Pi-TT^*-ML_PPi each computed entry then has error below
10^{-9}. The propagated LP error is at most
30· 2· 5· 10^{-12}=3· 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· 10^{-8}. Norm accumulation is below
10^{-8} even with residual norm 20161. Thus the actual
2· 10^{-6} addition in id suffices.


% ---------------------------------------------------------------
medskipnoindentRoot certificates and normal constants.

% ---------------------------------------------------------------

The tangent operator has the certified bound

 |T|_{rm op}lesqrt{1+ort}<1.00001,

since T is a genuine Root. The evaluation of T^*T-I for the stored
T (each |T_{ja}|le 1.01) yields ort<10^{-5}; the
evaluation of (I-Pi)T yields hor<10^{-5}; the evaluation
of F(P) yields 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· 10^{-11} Frobenius bound on the
difference H_VL_V^*-M (where M denotes the accumulated stored
product) gives a valid mnorm.

The two candidate choices for C are sqrt 8|M| and
sqrt 6|S|; the source takes the smaller of upper bounds
for these two quantities. Each is a valid pi b in
Lemma~phase:remainder by Lemma~phase:remainder-constants.
The root evaluation thresholds are val<10^{-18} for the
overall pass, ort<10^{-5}, hor<10^{-5},
id<10^{-4} after the 2· 10^{-6} addition, and
mnorm·mathtt{froot}lemathtt{n0}.

% ---------------------------------------------------------------
medskipnoindentCandidate conversion and ellipsoidal coverage.

% ---------------------------------------------------------------

These coverage estimates are needed only after a successful Candidate
preparation, at final unitary-grid levels N_uge 64.

Let m_j=|D_j| denote an ideal comparison-block magnitude. The stored
magnitude test sqrt{mag_j}>.15 implies m_j>.14999998.
Normalizing D by its computed magnitude has unit-direction error
below
[
 frac{2· 4· 10^{-9}}{.14999998}
     +norm/division rounding<5.5· 10^{-8}.
]
Candidate additionally checks a computed chord smaller than 10^{-7}
between this normalized value and the exponential of its proposed
phase. With (round:exp-budget), the true unit chord error is less
than 1.56· 10^{-7}. If a circular phase distance in cycles is
vin[0,1/2], its chord is 2sin(pi v)ge 4v. Across 16 entries
the phase-center correction therefore has norm below
1.56· 10^{-7}, which is below the 2· 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 mi subtraction of 10^{-6} therefore
gives a positive magnitude lower bound greater than .148, and
invmag_j=1/(widehatwidehat m_j-10^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
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
slack by mi has error below 10^{-12}:
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 (2pi)^2 in d2 each have error below
10^{-12} before their 10^{-8} additions. In particular the stored
quantities are upper bounds with

 arc<8.8, weighted<10,
  invmag_j<6.8.


medskipnoindentTangent offsets and external errors..
The computed off approximates 2pi T^*Delta within 10^{-9}
in norm, by~(round:t-operator), the projected-difference error,
and short-dot rounding. The error in the whole X array is at most
√(208)· 10^{-6} in Frobenius norm; there are 13 columns and
16 affected positions. Multiplication by T^* and the less than
2· 10^{-13} rounding of each dot give

 |AB-T^*X_{rm ideal}|<1.5· 10^{-5}.

The parameter-vector norm is below
√(.203^2+4(pi/64)^2)<.226, so substitution of stored AB costs
less than 3.4· 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 weighted and
invmag, with rational beta=.8. Each D-weighted dot has
absolute sum less than 16(1.00001)^2· 6.8<110. Its
two-product-per-term accumulation has error below 5· 10^{-13} by
(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
kle 10. Hence the exact stored G is positive semidefinite and
dominates the required covariance.

% ---------------------------------------------------------------
medskipnoindentThe stored spectral upper bound.

% ---------------------------------------------------------------


For an exact stored symmetric ktimes k matrix G with kle 10 and
|G_{ij}|le 10, spectralUpper returns an upper bound for
|G|.


Proof. 
An absolute row sum has error below 3· 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, Age m^2. Then
|fl(sqrt{fro4})-A|<4· 10^{-13}A+10^{-139}
by dot-product and norm rounding. Taking sqrt A costs less than
4· 10^{-13}sqrt A+10^{-69}, and since sqrt Ale 100 the total
possible downward error is below 10^{-10}. The code adds
2· 10^{-7}. The exact √(|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.


The 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{mx}+sqrt{gu}+10^{-8}, with 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
[
 gap=fl(gamma-gu-10^{-9})
       legamma-lambda_{max}(G).
]
Only gaps at least 10^{-6} are used. Taking max with zero and over
all offsets gives high. The ellipsoid formula uses
sqrt{gamma(1+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)


% ---------------------------------------------------------------
medskipnoindentThe base grid and the outside test.

% ---------------------------------------------------------------

Use the stored positive width w_B as an exact grid parameter. The
allowed range is 10^{-4}le w_Ble 0.05. The stored grid radius
E=sqrt{1+ort_P},R+10^{-7} exceeds the true tangent radius
by more than 0.999· 10^{-7}. For
L=leftlceil E/w_B+0.500001rightrceil, 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 Lle 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.

% ---------------------------------------------------------------
medskipnoindentSingle-vector filters and vector-ball construction.

% ---------------------------------------------------------------

The C stage validates its Root and 0le rle 0.03. At every bin
resolution Mge 6, the halfwidth, Euclidean bin radius, one-norm
bound, and infinity-norm bound are formed as
[
 
 a&=0.5/M+10^{-12},&
 X&=√(5-1/6),a+10^{-12},
 x_1&=5a+10^{-12},&
 x_∞&=1.5a+10^{-12}.
 
]
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_∞<0.126. The stored tangent and normal bounds
[
 r_t=sqrt{1+ort},r+10^{-9},
 r_n=n0+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|
 lesqrt{1+(2pisqrt2/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·
10^{-14}; |,|S|/sqrt 6-4,|<10^{-13}. The arrays J and B
have errors below 3· 10^{-11} and 4· 10^{-11}. The
rounded constants a,X,I_M,L_M,1/sqrt6,(2pi)^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· 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· 10^{-7}. The nested
squared-magnitude crossing band is 2· 10^{-7}; the squared
increment is 7· 10^{-7}.

medskipnoindentThe 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
[
 
 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|,
  M_z&=max_j(z_j^2+widetilde z_j^2).&
 
]
The support compared with |a_y| is the stored evaluation of

 &asum_{j>0}|z_j|+r_t|h|+r_n√(N_J)nonumber
 &+piBigl[
 m_y(KX+r)^2+
 minbigl{√(M_z),X, K√(N_S),x_∞bigr}X
 +2frac{√(N_S)}{sqrt6}Xr
 +frac{√(M_S)}{sqrt6}r^2 Bigr].
 

This is the necessary support bound of the Exclusion stage.
Its coefficient errors against ideal arrays satisfy
[
 |Delta a_y|<7· 10^{-12},
 |Delta z_j|,|Deltawidetilde z_j|<7· 10^{-12},
 |Delta h_a|<4.3· 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_yge 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· 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· 10^{-7}.

% ---------------------------------------------------------------
medskipnoindentThe joint coefficient support.

% ---------------------------------------------------------------

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· 10^{-11}. Each
accumulated shared-Hadamard coefficient has error below 9· 10^{-11}:
it receives at most twelve B input errors of 7· 10^{-12} plus
short-dot rounding. Each vector-phase coefficient, including its edge
contribution, has error below 9· 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· 10^{-10}. The code adds 3· 10^{-7}.

The final comparison allowance of 3· 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.


% ---------------------------------------------------------------
medskipnoindentFinite paths and NaN handling.

% ---------------------------------------------------------------

The preceding estimates apply to bounded finite certificate
evaluations. Exceptional heuristic values cannot bypass those
preconditions.
 In the seed cover, all direct lists come from finite integer
      midpoints and 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.
 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 trial::certificate checks every direction
      entry and its maximum before use.
 The snapping iteration calls gen, which rejects
      nonfinite or out-of-range phases. Root certification requires
      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.
 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.
 The cluster output is checked as in 
      at |v_j|le 2 and 0le 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.