#!/usr/bin/env python3 r"""lbl3se-w5-seed.py -- AMFlow-FREE derivation of the kite-family w=5 boundary Laurent seed (the single transport seed used by the LBL3SE/LBL3KP/LBL3VP blog evaluators, previously pinned to numeric/amflow_kite_bnd_w5_out.json). FAMILY (QED 2-loop electron-self-energy "kite", equal masses, m^2 = 1): loops k1,k2; external P with P^2 = w; d = 4-2*eps D1 = k1^2-1 D2 = k2^2-1 D3 = (P-k1-k2)^2-1 D4 = (P-k1)^2 D5 = (P-k2)^2 (3 massive electron lines, 2 massless photon lines.) 9 masters J = (J1..J9), AMFlow-bare normalization d^dk/(i pi^{d/2}) per loop, Feynman +i0 (physical sheet = w+i0), indices: J1=[1,1,0,0,0] J2=[1,1,0,1,0] J3=[1,1,1,0,0] J4=[1,1,2,0,0] J5=[1,1,3,0,0] J6=[0,0,1,1,1] J7=[0,0,2,1,1] J8=[1,1,0,1,1] J9=[1,1,1,1,1] SELF-CONTAINED DEFINITION (no AMFlow number anywhere in the computation): [seed] w=0 is the p^2=0 VACUUM point. Every master is analytic there (no cut reaches w=0: all Cutkosky cuts cross a massive line; thresholds sit at w=1 and w=9 only). The vacuum values are classical closed forms: T = -Gamma(-1+eps) (1-loop tadpole, m^2=1) J1(0)=J2(0)=J8(0) = T^2 J3(0) = U111 (equal-mass 2-loop vacuum sunset, series below) J4(0) = U112 = (d-3)/3 * U111 [vacuum IBP] J5(0) = U113 = ((d-8)*U112 + 2*Gamma(eps)^2)/6 [vacuum IBP] J6(0) = W = -Gamma(eps)Gamma(1-eps)^2 Gamma(2eps-1)/Gamma(2-eps) J7(0) = W2 = (d-3)*W J9(0) = U111 - 2*Umm0 + W, Umm0 = Gamma(eps)^2/((1-eps)(1-2eps)) Derivations (all from scratch, verified live): * U111 = int int 1/((k1^2-1)(k2^2-1)((k1+k2)^2-1)). Chain: inner equal-mass bubble in Feynman-parameter form B2(q^2) = G(e) Int_0^1 dx [1 - x(1-x) q^2]^{-e}, then U111 = Int d^dq B2(q^2)/(q^2-1); the radial q-integral, split as 1/(1+u) = 1/u - 1/(u(1+u)) so each piece sits INSIDE its Mellin convergence strip (the naive one-piece Euler formula drops a residue term -- verified failure mode), gives exactly U111 = -(G(e)/G(2-e)) [ B(1-e,2e-1)B(e,e) - B(1-e,2e) * Int_0^1 2F1(e,1-e;1+e;1-x(1-x)) dx ]. The z->1-z Gauss connection formula (c-a-b = e; the b=c term collapses to (1-y)^{-e}) and Int_0^1 [x(1-x)]^s dx = B(s+1,s+1) turn the x-integral into two 4^{-n}-convergent sums (this is the known 2F1(...;1/4)-class form of the equal-mass vacuum sunset, Davydychev-Tausk / Davydychev-Kalmykov): Int_0^1 2F1 dx = e*B(e,e)*S1 + [G(1+e)G(-e)/(G(e)G(1-e))]*S2 S1 = 2F1(e,1;3/2;1/4), S2 = sum_n [(2e)_n/(1+e)_n] B(n+1+e,n+1+e) Verified vs direct simplex quadrature (e=0.55, 0.7; 15 digits). Cross-checked at runtime against the independent FJJ/Martin eps^0 closed form (vendored): U111[-2]=3/2, U111[-1]=9/2-3g, U111[0]=(xi+21+3*zeta2)/2 - 9g + 3g^2, xi(1,1,1) = -6*sqrt(3)*Cl2(2pi/3). * U112, U113: vacuum IBP 0 = (D-2a-c)U(a,b,c) - 2a U(a+1,b,c) - c[U(a-1,b,c+1)-U(a,b-1,c+1)+U(a,b,c+1)] at (1,1,1),(1,1,2),(2,1,1), with U(2,0,2)=Gamma(eps)^2, U(0,1,3)=U(1,0,3) (k1<->k2). * W (sunset masses 1,0,0): massless bubble = c_G*(-l^2)^{-eps} then a radial Euler integral -> pure Gamma. W2 = dW/dm^2 = (1-2eps)W. * Umm0 (sunset masses 1,1,0): same chain; the x-integral collapses by Gauss 2F1(a,b;c;1) to Gamma's: Umm0 = G(e)^2/((1-e)(1-2e)). * J9(0): partial fractions 1/((k^2-1)k^2) = 1/(k^2-1) - 1/k^2 on both photon/electron pairs at P=0. J9(0) is ALSO fixed by ker A_{-1} membership (row 9) -- both expressions are computed and compared. [transport] The banked exact rational 9x9 connection dJ/dw = A(w,d) J (lbl3se-kite-de.json, Kira IBP, 43-sample reconstruction, max validation diff 0; singularities w in {0,1,9}) is eps-expanded to eps^4 and solved as a Laurent-block system on the 5 slots eps^{-2..+2}: * Frobenius at the regular-singular point w=0. Res_{w=0} A has eigenvalues {0 x4, -1+eps x4, -2+2eps} -- NO positive-integer and NO fractional eigenvalues, so the bounded (=physical, analytic) solution is the unique power-series branch with J(0) in ker A_{-1}; the recursion (n - A_{-1}) c_n = sum A_{n-1-m} c_m is invertible for all n>=1. Both facts are ASSERTED numerically at runtime, as is A_{-1} J(0) = 0. (Same threshold-vacuum anchor argument as tools/rowscripts/eval_row23.py.) * Taylor march from w1=2/5 along w1 -> w1+iH -> 5+iH -> 5 (H=2), i.e. through the UPPER half plane: Feynman +i0 = w+i0 branch. Per-step order NORD ~ 1.7*dps+20 with SAFETY=0.25 (per-step truncation ~ SAFETY^NORD; TOOL_CHANGELOG 2026-07-03 entry). OUTPUT eps-graded seed vector J_i(w=5), orders eps^{-2..+2}, at --dps N. The vendored lbl3se-w5-derived.json here and lbl3kp-w5-derived.json in ../lbl3kp (byte-copy) were emitted by this script (--dps 80 --json); the page evaluators consume THEM as the transport seed -- the AMFlow JSONs below are held-out cross-checks only. GATE (held out -- read ONLY for the final comparison, never used upstream): /phys_lbl3se/numeric/amflow_kite_bnd_w5_out.json (AMFlow solve_integrals, goal_digits=140; blog copy lbl3se-kite-boundary-w5.json). USAGE python3 lbl3se-w5-seed.py --dps 120 # seed + gate digits table python3 lbl3se-w5-seed.py --dps 120 --check # + rerun at dps+60, diff python3 lbl3se-w5-seed.py --dps 60 --json out.json python3 lbl3se-w5-seed.py --dps 30 --mutate # NEGATIVE CONTROL: corrupt one vacuum constant by rel 1e-25 -> gate must collapse to ~25d (mp.dps is set inside main AFTER argparse -- the tools/dps_lint.py rule.) CHANGELOG 2026-07-05 (axis3 wave) WP1 certified-bound treatment (pilot pattern, /wiring/WIRING_LOG.md items 12-13): NFROB and NORD are STARTING guesses only. The Frobenius sum and every Taylor-march step now carry a runtime-certified trailing-window geometric tail bound (Frobenius: ratio r = w1 = 2/5, radius 1 to the nearest other singularity; march: r = |h|/dist <= SAFETY, singularities {0,1,9} verified at load); on failure the SAME recurrences are continued exactly (N x1.5) to cap 8*N0, RuntimeError at cap. Per-step march bounds accumulate and print as [cert] lines. The held-out AMFlow gate is now a RAISING gate (bar = min(dps,128)-12, calibrated: measured worst ~ dps-0.1). Guards: FROB_GUARD/MARCH_GUARD below; calibration in /axis3_wave/lbl3se-wp3/ROW_REPORT.md. Deps: python3 + mpmath + the sibling kite_de_exact_parse.py (it parses and expands the shipped exact rational connection in exact rational arithmetic; no sympy; all numerics are mpmath). Wall time is measured and printed. """ import argparse import hashlib import importlib.util import json import os import sys import time import mpmath as mp HERE = os.path.dirname(os.path.abspath(__file__)) DE_JSON = os.path.join(HERE, 'lbl3se-kite-de.json') # sha256 of the pinned sibling module (computed by the producer that cut this file, never typed); a byte change # or a missing file is REFUSED by name (exit 3) before any value is printed SCRIPT = "lbl3se-w5-seed" PINS = { "kite_de_exact_parse.py": "0cc967d6362b63963b2e7ddd3e41f9edc1a4a0cab9cd2d7ebd7aef45aa83d7c1", } def _pinned_path(name): """Path of a sibling file; a PINNED file (the exact parser module) is refused on any byte change (exit 3, the recorded and recomputed sha256 named).""" path = os.path.join(HERE, name) want = PINS.get(name) if want is not None: if not os.path.exists(path): sys.stderr.write(f"{SCRIPT} REFUSED: pinned file {name} is missing (recorded sha256 {want})\n") raise SystemExit(3) have = hashlib.sha256(open(path, "rb").read()).hexdigest() if have != want: pos = next((k + 1 for k, (x, y) in enumerate(zip(want, have)) if x != y), 0) sys.stderr.write(f"{SCRIPT} REFUSED: {name} integrity pin mismatch (recorded {want}, recomputed " f"{have}; first differing hex position {pos} of 64, 1-based) -- the shipped file " f"was altered\n") raise SystemExit(3) return path def _import_pinned_module(name): """Import a sibling module through the pin check (a byte change is REFUSED, exit 3).""" spec = importlib.util.spec_from_file_location(name[:-3], _pinned_path(name)) mod = importlib.util.module_from_spec(spec) spec.loader.exec_module(mod) return mod kx = _import_pinned_module("kite_de_exact_parse.py") # the exact parser of the connection strings (no sympy) GATE_BANKED = ('/phys_lbl3se/numeric/' 'amflow_kite_bnd_w5_out.json') GATE_LOCAL = os.path.join(HERE, 'lbl3se-kite-boundary-w5.json') MASTERS = [(1,1,0,0,0),(1,1,0,1,0),(1,1,1,0,0),(1,1,2,0,0),(1,1,3,0,0), (0,0,1,1,1),(0,0,2,1,1),(1,1,0,1,1),(1,1,1,1,1)] NM = 9 K0 = -2 # lowest eps order kept NS = 5 # slots: eps^{-2..+2} NAEPS = 5 # A(w,eps) expanded to eps^{0..4} SING = (0, 1, 9) # singularities of A (verified at load); exact ints SAFETY = 0.25 # exact binary float (per-step truncation ~ SAFETY^NORD) # axis3 wave (2026-07-05): certified-bound guards. Closed-form N's are # STARTING guesses; refine-until-bound loops enforce these tolerances. TAIL_WINDOW = 8 # trailing-term window of the geometric envelope FROB_GUARD = 6 # Frobenius tail tol = 10^-(dps+FROB_GUARD); calibrated # 2026-07-05: measured bound = 10^-(dps+8.4) (dps 40/80), # design w1^N0 = 10^-(dps+9.95) x scale ~39 -> ~10^2.4 # headroom, uniform in dps; guard 8 left only 1.4x MARCH_GUARD = 12 # per-step march tail tol = 10^-(dps+MARCH_GUARD) NCAP_FACT = 8 # series-order cap = NCAP_FACT * starting N GATE_MARGIN = 12 # held-out gate RAISE bar = min(dps, 128) - GATE_MARGIN # ============================ eps-Laurent series ============================= # Ser: sum_j c[j] eps^{k0+j}; dense window, mpf/mpc coefficients. # (Vendored pattern: /solve_row23/main/vacuum_seed.py.) class Ser: __slots__ = ("k0", "c") def __init__(self, k0, c): self.k0 = k0; self.c = list(c) @staticmethod def const(v, n): return Ser(0, [mp.mpf(v)] + [mp.mpf(0)]*(n-1)) def __mul__(self, o): if isinstance(o, Ser): n = min(len(self.c), len(o.c)) out = [mp.mpf(0)]*n for i in range(n): ai = self.c[i] if ai == 0: continue for j in range(n-i): out[i+j] += ai*o.c[j] return Ser(self.k0+o.k0, out) return Ser(self.k0, [x*o for x in self.c]) __rmul__ = __mul__ def __add__(self, o): if not isinstance(o, Ser): o = Ser.const(o, len(self.c)) k0 = min(self.k0, o.k0) n = min(self.k0+len(self.c), o.k0+len(o.c)) - k0 out = [mp.mpf(0)]*n for i, v in enumerate(self.c): j = self.k0 - k0 + i if j < n: out[j] += v for i, v in enumerate(o.c): j = o.k0 - k0 + i if j < n: out[j] += v return Ser(k0, out) __radd__ = __add__ def __sub__(self, o): return self + (o*(-1) if isinstance(o, Ser) else -o) def __neg__(self): return self*(-1) def inv(self): # normalize away exact leading zeros (e.g. poly factors like 2*eps) k0, cc = self.k0, self.c while cc and cc[0] == 0: k0 += 1; cc = cc[1:] self = Ser(k0, cc) a0 = self.c[0] n = len(self.c); out = [mp.mpf(0)]*n; out[0] = 1/a0 for k in range(1, n): s = mp.mpf(0) for j in range(1, k+1): s += self.c[j]*out[k-j] out[k] = -s/a0 return Ser(-self.k0, out) def coeff(self, k): j = k - self.k0 return self.c[j] if 0 <= j < len(self.c) else mp.mpf(0) def window(self, k0, n): """re-window to [k0, k0+n).""" return Ser(k0, [self.coeff(k0+j) for j in range(n)]) def ser_exp(s): assert s.k0 >= 0 n = len(s.c) full = [mp.mpf(0)]*n for i, v in enumerate(s.c): if s.k0 + i < n: full[s.k0+i] = v out = [mp.mpf(0)]*n out[0] = mp.exp(full[0]) if full[0] != 0 else mp.mpf(1) g = full[:]; g[0] = mp.mpf(0) for k in range(1, n): s_ = mp.mpf(0) for j in range(1, k+1): s_ += j*g[j]*out[k-j] out[k] = s_/k return Ser(0, out) def gamma_ser(x0, coef, n): """Gamma(x0 + coef*eps) as Ser; poles at nonpositive-integer x0 handled by Gamma(z) = Gamma(z+m)/prod_{j0).""" tol = mp.mpf(10)**(-(dps+10)) G = gamma_ser # S1 = sum (e)_k/((3/2)_k 4^k) S1 = Ser(0, [mp.mpf(0)]*n); term = Ser.const(1, n) k = 0 while True: S1 = S1 + term term = term * lin(k, 1, n) * (1/(4*(mp.mpf(3)/2 + k))) k += 1 if max(abs(x) for x in term.c) < tol and k > 4: break if k > 200000: raise RuntimeError("S1 no convergence") # S2 = sum [(2e)_n/(1+e)_n] B(n+1+e, n+1+e) S2 = Ser(0, [mp.mpf(0)]*n) poch = Ser.const(1, n) # (2e)_k/(1+e)_k k = 0 while True: B = G(k+1, 1, n)*G(k+1, 1, n)*G(2*k+2, 2, n).inv() t = poch * B S2 = S2 + t poch = poch * lin(k, 2, n) * lin(1+k, 1, n).inv() k += 1 if max(abs(x) for x in t.c) < tol and k > 4: break if k > 200000: raise RuntimeError("S2 no convergence") Bee = G(0,1,n)*G(0,1,n)*G(0,2,n).inv() # B(e,e) B1e_2em1 = G(1,-1,n)*G(-1,2,n)*G(0,1,n).inv() # B(1-e,2e-1) B1e_2e = G(1,-1,n)*G(0,2,n)*G(1,1,n).inv() # B(1-e,2e) eBee = Ser(1, [mp.mpf(1)] + [mp.mpf(0)]*(n-1)) * Bee cf2 = G(1,1,n)*G(0,-1,n)*(G(0,1,n)*G(1,-1,n)).inv() intF = eBee*S1 + cf2*S2 bracket = B1e_2em1*Bee - B1e_2e*intF return -(G(0,1,n)*G(2,-1,n).inv()) * bracket def vacuum_seed(dps): """The 9 vacuum values J_i(w=0) as Ser windows [K0, K0+NS).""" n = NS + 4 # working slots; re-windowed at the end G = gamma_ser T = -G(-1, 1, n) # tadpole T2 = T*T U111 = U111_ser(n, dps) # --- FJJ/Martin eps^0 anchor (independent literature closed form) --- g = mp.euler xi = -6*mp.sqrt(3)*mp.clsin(2, 2*mp.pi/3) # xi(1,1,1), Cl2(2pi/3) fjj = [mp.mpf(3)/2, mp.mpf(9)/2 - 3*g, (xi + 21 + 3*mp.zeta(2))/2 - 9*g + 3*g*g] err = max(abs(U111.coeff(-2+j) - fjj[j]) for j in range(3)) assert err < mp.mpf(10)**(-dps+8), f"U111 vs FJJ anchor: {err}" U112 = U111 * lin(1, -2, n) * (mp.mpf(1)/3) # (d-3)/3 * U111 Ge2 = G(0, 1, n)**2 if False else G(0, 1, n)*G(0, 1, n) U113 = (lin(-4, -2, n)*U112 + 2*Ge2) * (mp.mpf(1)/6) W = -(G(0,1,n) * G(1,-1,n)*G(1,-1,n) * G(-1,2,n) * G(2,-1,n).inv()) W2 = lin(1, -2, n) * W Umm0 = Ge2 * (lin(1,-1,n)*lin(1,-2,n)).inv() J9 = U111 - 2*Umm0 + W seed = [T2, T2, U111, U112, U113, W, W2, T2, J9] # slot format: state[t][i] = coeff of eps^{K0+t} in J_i(0) slots = [[seed[i].coeff(K0+t) for i in range(NM)] for t in range(NS)] return slots, {'U111': U111, 'xi_err': err} # ===================== [transport] exact rational connection ================= def load_connection(): """Parse the banked exact rational A(w,d), substitute d=4-2eps, and return per-eps-order sparse rational entries: comps[j] = list of (i, jc, num_coeffs_asc, den_coeffs_asc) j=0..NAEPS-1 with exact-rational -> mpf coefficients (exact at current dps). The parse and the eps-expansion are the sibling kite_de_exact_parse.py's (exact rational arithmetic, the same integer lists in lowest terms the earlier sympy parse produced; no sympy); the singular points are proved to lie in {0,1,9} by exact division; a string outside the grammar is REFUSED by name (exit 3).""" J = json.load(open(DE_JSON)) assert [tuple(m) for m in J['masters']] == list(MASTERS) def wpoly_mpf(coeffs): return [mp.mpf(c.numerator)/mp.mpf(c.denominator) for c in coeffs] try: lists = kx.graded_lists(J['A'], range(NM), NAEPS - 1) except kx.ParseError as ex: sys.stderr.write(f"{SCRIPT} REFUSED: {os.path.basename(DE_JSON)}: {ex}\n") raise SystemExit(3) try: sing = kx.singular_certificate(lists) except kx.ParseError as ex: raise AssertionError(f"unexpected singularities: {ex}") comps = [[] for _ in range(NAEPS)] for (i, jc, j), (num, den) in lists.items(): comps[j].append((i, jc, wpoly_mpf(num), wpoly_mpf(den))) assert sing <= {0, 1, 9}, f"unexpected singularities {sing}" return comps def _shift_asc(coefs, x0): """Taylor coeffs about x0 of poly with ASCENDING coeffs (Horner shift).""" p = [coefs[-1]] for c in reversed(coefs[:-1]): pn = [c] + [0]*len(p) for m_, a in enumerate(p): pn[m_] += a*x0 pn[m_+1] += a p = pn return p def a_taylor(comps, w0, N): """Aser[j][m] = sparse list of (i, jc, coeff) : Taylor coefficient m about w0 of the eps^j component of A. Requires w0 not a root of any den.""" out = [] for j in range(NAEPS): ser_m = [[] for _ in range(N+1)] for (i, jc, nc, dc) in comps[j]: ph = _shift_asc(nc, w0) qh = _shift_asc(dc, w0) q0 = qh[0] r = [1/q0] for m_ in range(1, N+1): s = 0 for l in range(1, min(m_, len(qh)-1)+1): s += qh[l]*r[m_-l] r.append(-s/q0) for m_ in range(N+1): s = 0 for l in range(min(m_, len(ph)-1)+1): s += ph[l]*r[m_-l] if s != 0: ser_m[m_].append((i, jc, s)) out.append(ser_m) return out def a_laurent0(comps, N): """Laurent matrices of A about w=0: returns (Am1, Alist) where Am1[j] = dense 9x9 (mpf) eps^j component of Res_{w=0}A, and Alist[j][m] = sparse Taylor part (m=0..N).""" Am1 = [[[mp.mpf(0)]*NM for _ in range(NM)] for _ in range(NAEPS)] Alist = [] for j in range(NAEPS): ser_m = [[] for _ in range(N+1)] for (i, jc, nc, dc) in comps[j]: # den = w^p * reg(w); p<=1 (verified for this family) p = 0 dcl = list(dc) while abs(dcl[0]) == 0: dcl = dcl[1:]; p += 1 assert p <= 1 q0 = dcl[0] r = [1/q0] for m_ in range(1, N+2): s = 0 for l in range(1, min(m_, len(dcl)-1)+1): s += dcl[l]*r[m_-l] r.append(-s/q0) # entry = num(w)/(w^p reg(w)) : coefficients of w^{m-p} full = [] for m_ in range(N+2): s = 0 for l in range(min(m_, len(nc)-1)+1): s += nc[l]*r[m_-l] full.append(s) if p == 1: if full[0] != 0: Am1[j][i][jc] = full[0] coefs = full[1:] else: coefs = full[:N+1] for m_ in range(N+1): if coefs[m_] != 0: ser_m[m_].append((i, jc, coefs[m_])) Alist.append(ser_m) return Am1, Alist # state = list of NS slots, each a length-9 vector (slot s <-> eps^{K0+s}) def matvec_slots(Asparse, state, out=None): """out[t] += sum_j Asparse[j] . state[t-j]""" if out is None: out = [[mp.mpf(0)]*NM for _ in range(NS)] for j in range(NAEPS): Aj = Asparse[j] if not Aj: continue for t in range(j, NS): st = state[t-j] ot = out[t] for (i, jc, a) in Aj: ot[i] += a*st[jc] return out def frobenius(Am1, Alist, seed, w1, N, comps=None, dps=None): """Power-series branch at w=0: c_0 = seed, (n - A_{-1}) c_n = rhs_n. Returns state at w1 (real). axis3 wave: N is a STARTING guess. The trailing-TAIL_WINDOW geometric tail bound (ratio r = w1/1, radius 1 = distance from 0 to the nearest other verified singularity) must beat 10^-(dps+FROB_GUARD), else the SAME recursion continues exactly at N*1.5 (Alist re-derived to the higher order -- identical lower coefficients) up to NCAP_FACT*N; RuntimeError at cap. Returns (acc, r0, ev, certified tail bound, N used).""" # residual check: A_{-1} c_0 = 0 res = matvec_slots(Am1_sparse(Am1), seed) r0 = max(abs(x) for sl in res for x in sl) # eigenvalues at eps=0: no positive integers (Frobenius solvability) E = mp.matrix([[Am1[0][i][j] for j in range(NM)] for i in range(NM)]) ev = mp.eig(E, left=False, right=False) for lam in ev: assert abs(mp.im(lam)) < 1e-10 and (mp.re(lam) < 0.5 or abs(lam) < 1e-10), \ f"resonant eigenvalue {lam}" c = [ [list(sl) for sl in seed] ] # c[m][slot][i] acc = [ [x for x in sl] for sl in seed ] # running value sum c_m w1^m wp = [mp.mpf(1)] # running power w1^n (list: mutable) Ineye = mp.eye(NM) def _extend(n_lo, n_hi, Alist): for n_ in range(n_lo, n_hi+1): rhs = [[mp.mpf(0)]*NM for _ in range(NS)] for m_ in range(n_): k = n_ - 1 - m_ for j in range(NAEPS): lst = Alist[j][k] if k < len(Alist[j]) else [] if not lst: continue for t in range(j, NS): st = c[m_][t-j] ot = rhs[t] for (i, jc, a) in lst: ot[i] += a*st[jc] # solve (n - A_{-1}) x = rhs slotwise: L0 = n I - Am1[0] L0 = n_*Ineye - mp.matrix([[Am1[0][i][j] for j in range(NM)] for i in range(NM)]) xs = [] for t in range(NS): b = [rhs[t][i] for i in range(NM)] for j in range(1, NAEPS): if t-j < 0: break Aj = Am1[j] xprev = xs[t-j] for i in range(NM): row = Aj[i] s = mp.mpf(0) for jc in range(NM): if row[jc] != 0: s += row[jc]*xprev[jc] b[i] += s sol = mp.lu_solve(L0, mp.matrix(b)) xs.append([sol[i] for i in range(NM)]) c.append(xs) wp[0] *= w1 for t in range(NS): at = acc[t]; xt = xs[t] for i in range(NM): at[i] += xt[i]*wp[0] _extend(1, N, Alist) r = w1 / 1 # certified envelope ratio (radius 1) fac = r / (1 - r) def _bound(): mx = mp.mpf(0) pw = wp[0] / w1**(TAIL_WINDOW - 1) for m_ in range(len(c) - TAIL_WINDOW, len(c)): for sl in c[m_]: for x in sl: ax = abs(x)*pw if ax > mx: mx = ax pw *= w1 return mx*fac def _frob_finite(b, n): # 2026-09-06 (Q20f): fail CLOSED on a non-finite tail bound BEFORE the # continue-test below -- `bound >= tol` is False for a NaN by IEEE # ordering, so a NaN bound would leave the loop as 'converged' (the # same fail-open form the evaluator's refine loop carried until # 2026-09-06). The bound is the trailing-window max |c| x r/(1-r), so # a non-finite bound names a non-finite Frobenius coefficient or # envelope factor at this order. if not mp.isfinite(b): raise RuntimeError( f"[cert] Frobenius at w1={mp.nstr(mp.mpf(w1), 8)}: the certified " f"tail bound (trailing-window max |c| x r/(1-r)) is not finite " f"(nan/inf) at N={n} -- fail closed") tol = mp.mpf(10)**(-(dps + FROB_GUARD)) if dps else None ncap = NCAP_FACT*N bound = _bound() _frob_finite(bound, N) if tol is not None: while bound >= tol: if N >= ncap: raise RuntimeError( f"[cert] Frobenius at w1={mp.nstr(mp.mpf(w1), 8)}: certified " f"tail bound {mp.nstr(bound, 4)} >= tol {mp.nstr(tol, 4)} " f"at N={N} (cap {ncap}) -- refine-until-bound FAILED, STOP") N2 = min(ncap, int(1.5*N) + 1) _, Alist = a_laurent0(comps, N2) # same recursion, higher order _extend(N+1, N2, Alist) N = N2 bound = _bound() _frob_finite(bound, N) return acc, r0, ev, bound, N def Am1_sparse(Am1): out = [] for j in range(NAEPS): lst = [] for i in range(NM): for jc in range(NM): if Am1[j][i][jc] != 0: lst.append((i, jc, Am1[j][i][jc])) out.append(lst) return out def _step_recursion(C, Aser, m_lo, m_hi): """Continue the state-coefficient recursion for orders m_lo..m_hi-1 (EXACT continuation: same recurrence, same lower coefficients).""" for m_ in range(m_lo, m_hi): nxt = [[mp.mpc(0)]*NM for _ in range(NS)] for j in range(NAEPS): Aj = Aser[j] for l in range(m_+1): lst = Aj[l] if not lst: continue Cm = C[m_-l] for t in range(j, NS): st = Cm[t-j] ot = nxt[t] for (i, jc, a) in lst: ot[i] += a*st[jc] inv = mp.mpf(1)/(m_+1) C.append([[x*inv for x in sl] for sl in nxt]) def taylor_step(comps, state, w0, h, NORD, dmin=None, tol=None, ncap=None): """One local-Taylor step of the eps-graded system; returns state at w0+h. axis3 wave: NORD is a STARTING guess; the trailing-TAIL_WINDOW geometric tail bound (ratio r = |h|/dmin <= SAFETY, singularities verified {0,1,9}) must beat tol, else the SAME recurrence continues exactly at N*1.5 up to ncap (RuntimeError at cap). Returns (state, bound, N used).""" Aser = a_taylor(comps, w0, NORD) C = [ [list(sl) for sl in state] ] # C[m][slot][i] _step_recursion(C, Aser, 0, NORD) bound = None if tol is not None: absh = abs(h) r = absh/dmin fac = r/(1 - r) while True: mx = mp.mpf(0) pw = absh**(NORD - TAIL_WINDOW + 1) for m_ in range(NORD - TAIL_WINDOW + 1, NORD + 1): for sl in C[m_]: for x in sl: ax = abs(x)*pw if ax > mx: mx = ax pw *= absh bound = mx*fac if bound < tol: break if NORD >= ncap: raise RuntimeError( f"[cert] march step at w0={mp.nstr(mp.mpc(w0), 10)}, " f"h={mp.nstr(mp.mpc(h), 10)}: certified tail bound " f"{mp.nstr(bound, 4)} >= tol {mp.nstr(tol, 4)} at " f"N={NORD} (cap {ncap}) -- refine-until-bound FAILED, STOP") N2 = min(ncap, int(1.5*NORD) + 1) Aser = a_taylor(comps, w0, N2) _step_recursion(C, Aser, NORD, N2) NORD = N2 out = [[mp.mpc(0)]*NM for _ in range(NS)] hp = mp.mpc(1) for m_ in range(NORD+1): Cm = C[m_] for t in range(NS): ot = out[t]; ct = Cm[t] for i in range(NM): ot[i] += ct[i]*hp hp *= h return out, bound, NORD def march(comps, state, z0, targets, NORD, verbose=False, dps=None): """axis3 wave: every step certified (tol 10^-(dps+MARCH_GUARD), cap NCAP_FACT*NORD); returns (state, nstep, accumulated bound, worst bound, escalations, max N used).""" nstep = 0 tol = mp.mpf(10)**(-((dps or mp.mp.dps) + MARCH_GUARD)) ncap = NCAP_FACT*NORD err = mp.mpf(0) worst = mp.mpf(0) ngrow = 0 maxn = NORD for zt in targets: while abs(zt - z0) > mp.mpf(10)**(-mp.mp.dps): dist = min(abs(z0 - s) for s in SING) hmax = SAFETY*dist dz = zt - z0 h = dz if abs(dz) <= hmax else dz/abs(dz)*hmax state, bnd, nu = taylor_step(comps, state, z0, h, NORD, dmin=dist, tol=tol, ncap=ncap) err += bnd if bnd > worst: worst = bnd if nu > NORD: ngrow += 1 if nu > maxn: maxn = nu z0 = z0 + h nstep += 1 if verbose: print(f" step {nstep:3d} -> {mp.nstr(z0, 8)}") return state, nstep, err, worst, ngrow, maxn # ================================ gate ======================================= def parse_arb(s): s = str(s).strip() if s.startswith('['): return mp.mpf(s[1:].split('+/-')[0].strip()) return mp.mpf(s) def load_gate(): path = GATE_BANKED if os.path.exists(GATE_BANKED) else GATE_LOCAL J = json.load(open(path)) out = {} for r in J['result']: idx = tuple(r['integral']['indices']) out[idx] = {c['order']: mp.mpc(parse_arb(c['value']['re']), parse_arb(c['value']['im'])) for c in r['coefficients']} return out, path def run(dps, verbose=True, mutate=False): mp.mp.dps = dps + 15 t0 = time.time() seed, seedinfo = vacuum_seed(dps) if mutate: # NEGATIVE CONTROL: corrupt one vacuum constant (U111 eps^0 slot of J3) # by a relative 1e-25 and confirm the held-out gate COLLAPSES to ~25d. bump = seed[2][2]*mp.mpf(10)**(-25) seed[2][2] += bump print(f"[MUTATE] J3(0) eps^0 += {mp.nstr(bump, 3)} (rel 1e-25) -- " f"gate must collapse to ~25d") if verbose: print(f"[seed] vacuum closed forms built; FJJ eps^0 anchor err " f"{mp.nstr(seedinfo['xi_err'], 3)}") comps = load_connection() if verbose: print(f"[DE] exact rational connection loaded, eps^0..{NAEPS-1}, " f"sing {{0,1,9}} ({time.time()-t0:.1f}s)") w1 = mp.mpf(2)/5 NFROB = int(dps/mp.log10(1/w1)) + 25 Am1, Alist = a_laurent0(comps, NFROB) acc, r0, ev, fbound, nfrob_used = frobenius(Am1, Alist, seed, w1, NFROB, comps=comps, dps=dps) if verbose: evs = sorted(set(round(float(mp.re(l)), 6) for l in ev)) print(f"[frob] Res_0 A eigenvalues(eps=0) ~ {evs}; |A_-1 c_0| = " f"{mp.nstr(r0, 3)}; N={nfrob_used}, tail ~ {mp.nstr(fbound, 3)} " f"({time.time()-t0:.1f}s)") print(f"[cert] Frobenius certified tail bound {mp.nstr(fbound, 3)} < tol " f"1e-{dps+FROB_GUARD} (r=2/5, window {TAIL_WINDOW}; N start {NFROB}, " f"used {nfrob_used}, cap {NCAP_FACT*NFROB})") if mutate: print(f"[MUTATE] ker-A_-1 residual now {mp.nstr(r0, 3)} (assert waived)") else: assert r0 < mp.mpf(10)**(-dps+10), "seed not in ker A_{-1}" state = [[mp.mpc(x) for x in sl] for sl in acc] NORD = int(mp.ceil(mp.mpf('1.7')*(dps+15))) + 20 H = mp.mpf(2) state, nstep, merr, mworst, mgrow, mmaxn = march( comps, state, mp.mpc(w1), [mp.mpc(w1, H), mp.mpc(5, H), mp.mpc(5)], NORD, dps=dps) if verbose: print(f"[march] {nstep} Taylor steps (NORD={NORD}, SAFETY={SAFETY}), " f"path w1 -> w1+{H}i -> 5+{H}i -> 5 ({time.time()-t0:.1f}s)") print(f"[cert] march: {nstep} steps certified, per-step tol " f"1e-{dps+MARCH_GUARD}; worst step bound {mp.nstr(mworst, 3)}; " f"accumulated {mp.nstr(merr, 3)}; escalations {mgrow} " f"(N start {NORD}, max used {mmaxn}, cap {NCAP_FACT*NORD})") return state, time.time()-t0 def main(): ap = argparse.ArgumentParser(description="AMFlow-free kite w=5 boundary seed") ap.add_argument('--dps', type=int, default=60) ap.add_argument('--check', action='store_true', help="rerun at dps+60 and diff (two-precision rule)") ap.add_argument('--json', default=None, help="write seed vector JSON here") ap.add_argument('--mutate', action='store_true', help="negative control: corrupt one vacuum constant " "(rel 1e-25) and show the gate collapse") ap.add_argument('-q', '--quiet', action='store_true') args = ap.parse_args() dps = args.dps state, wall = run(dps, verbose=not args.quiet, mutate=args.mutate) mp.mp.dps = dps + 15 gate, gpath = load_gate() print(f"\n=== eps-graded seed vector at w=5 (dps {dps}) and HELD-OUT gate " f"({os.path.basename(gpath)}) ===") worst = mp.inf digs = {} for i, idx in enumerate(MASTERS): for t in range(NS): o = K0 + t if o not in gate[idx]: continue mine = state[t][i] ref = gate[idx][o] scale = max(abs(ref), mp.mpf(1)) d_ = -mp.log10(abs(mine-ref)/scale) if mine != ref else mp.mpf(dps) digs[(idx, o)] = d_ worst = min(worst, d_) line = " ".join(f"e^{K0+t}:{mp.nstr(digs[(idx, K0+t)], 5)}d" for t in range(NS) if (idx, K0+t) in digs) print(f" J{i+1} {str(idx):18s} {line}") print(f" WORST gate agreement: {mp.nstr(worst, 6)} digits " f"(oracle goal_digits=140)") print(f" wall {wall:.1f}s") bar = min(dps, 128) - GATE_MARGIN print(f" [cert] held-out gate: worst {mp.nstr(worst, 5)} d >= bar {bar} d " f"(min(dps,128)-{GATE_MARGIN}) -- RAISES on fail") if not worst >= bar: raise RuntimeError(f"[cert] held-out AMFlow w=5 gate: worst " f"{mp.nstr(worst, 6)} d < bar {bar} d at dps {dps} " "-- seed not certified, STOP") if args.json: out = {'w': 5, 'dps': dps, 'masters': [list(m) for m in MASTERS], 'orders': list(range(K0, K0+NS)), 'values': [[[mp.nstr(mp.re(state[t][i]), dps), mp.nstr(mp.im(state[t][i]), dps)] for t in range(NS)] for i in range(NM)], 'provenance': 'lbl3se-w5-seed.py vacuum->w=5 transport, AMFlow-free'} json.dump(out, open(args.json, 'w'), indent=1) print(f" seed vector written to {args.json}") if args.check: print(f"\n[check] rerun at dps {dps+60} ...") state2, wall2 = run(dps+60, verbose=False) mp.mp.dps = dps + 90 stab = mp.inf for i in range(NM): for t in range(NS): a, b = state[t][i], state2[t][i] sc = max(abs(b), mp.mpf(1)) if a != b: stab = min(stab, -mp.log10(abs(a-b)/sc)) print(f"[check] two-precision stability: {mp.nstr(stab, 6)} digits " f"(wall {wall2:.1f}s)") if __name__ == '__main__': main()