#!/usr/bin/env python3 r"""kite-equal-evaluate.py -- the EQUAL-mass kite (row 13 of the paper's table of thirty integrals): standalone evaluator of the eps^0 top master J(1,1,1,1,1). This is the third configuration of the kite bundle. kite-evaluate.py is the UNEQUAL-mass kite (masses 1, 0, 1, 0, sqrt2; a non-modular curve, evaluated by transport) and kite-threshold-evaluate.py is the threshold configuration; the script here is the equal-mass kite, the one configuration with a finite modular closed form, and it evaluates that closed form directly. Two-loop kite, propagators D1=k1^2-1, D2=k2^2, D3=(k1-k2)^2-1, D4=(k1-p)^2, D5=(k2-p)^2-1 (masses m,0,m,0,m all equal, m^2=mu^2=1), measure int d^dk/(i pi^{d/2}) per loop, d=4-2eps, NO e^{gamma_E eps}. Kinematic variable t = p^2/m^2. The eps^0 top master is finite and equals the five-word modular Gamma_1(6) closed form printed in the equal-mass kite section of the paper (blind fit + PSLQ, all five coefficients clean): I_kite(t) = (1/(4t)) [ -(2 pi^2/3) G(1;t) - 8 G(0,1,1;t) + 4 G(1,0,1;t) - 108 Cl2(2pi/3) Ebar_{1;-1}(z3;1;-q) - 108 Ebar_{0,2;-2,0;2}(z3,z3;1,-1;-q) ], z3 = e^{2 pi i/3}, with q = e^{i pi tau}, tau = psi2/psi1 the geometric nome of the sunrise curve Dtilde = (t-1)^3 (t-9) (Adams-Bogner-Schweitzer-Weinzierl conventions, arXiv:1607.01571; Ebar/ELi defined in eqs. def_Ebar_weight_1, Ebar_integration there). Everything is computed AT RUNTIME from these definitions: AGM periods, eta-quotient hauptmodul q-series, ELi q-series, weight<=3 {0,1}-HPL closed forms. NOME BRANCH: the naive ellipk principal-branch nome is only a modular REPRESENTATIVE; away from the cusp-connected branch it can give the right Re and the WRONG Im (a lesson learned on a sibling elliptic integral). On the EUCLIDEAN domain of this script the ellipk nome IS the cusp-connected branch, and we PROVE it at runtime: the nome is independently recomputed by Newton path-tracking of the level-6 eta-quotient hauptmodul t(q_s)=9 q_s prod(...) (dictionary q_AW = -q_sunrise) from the t=0 cusp / t=-3 base point, and the two must agree to working precision or the script refuses. The physical region t>1 is reached by the elliptic Feynman-i0 prescription of Bogner, Schweitzer and Weinzierl -- the nome continued through the upper half t-plane, NOT a naive t+i0 per piece -- see PHYSICAL REGION below. DOMAIN: Euclidean t < 1, t != 0 (i.e. t < 0 or 0 < t < 1, both below the t=m^2 threshold) on the real-axis path, and the physical region t > 1, t != 9, on the continued path (--point t > 1: a complex value, see PHYSICAL REGION). Near the cusps t -> 0 (trivial zero of q), t -> 1 (from either side), t -> 9 and t -> -infinity the nome |q| -> 1; points with |q| >= QMAX=0.90 are refused explicitly, and the cusps t = 0, 1, 9 themselves by name. Precision is capped only by --dps (q-series lengths are auto-sized from measured coefficient growth + |q|). PHYSICAL REGION (t > 1; --point t > 1, added 2026-09-09): the closed form is continued with the elliptic Feynman-i0 prescription. The nome q_s(t + i0) is the served hauptmodul inverse followed along a polygon through the UPPER half t-plane from the Euclidean base point t_b = 1/2 (route H: Newton path-tracking of the served eta-quotient series), and independently the periods (psi1, psi2) are transported along the same polygon by Taylor steps of the Picard-Fuchs operator of the sunrise curve, t(t-1)(t-9) psi'' + (3t^2-20t+9) psi' + (t-3) psi = 0 (route P; the operator is verified at runtime on the served ellipk periods at three Euclidean points), q_P = exp(i pi psi2/psi1) with the served dictionary q_AW = -q_sunrise. Route H must agree with route P, and with itself at a second polygon height, to dps-6 digits, or the point is refused. The principal-branch ellipk nome at t + i delta is a DIFFERENT modular representative (it agrees with the continued nome to under one digit at every gated point); it is printed beside for the record and never used. The HPL pieces take their t + i0 boundary values, log(1-y) = log(t-1) - i pi, Li2(1-y) and Li3(1-y) real, Li2(t + i0) = pi^2/3 - log(t)^2/2 - Li2(1/t) + i pi log t; the served G(0,1,1) / G(1,0,1) evaluated at y = t + i 10^-(dps/2) must agree with these to about dps/2 digits (a runtime control). The words, the five coefficients and the N/2N certification are the served ones with the complex nome; the value is complex and prints as Re / Im. The imaginary part carries the Feynman +i0 sign: Im I_kite(t + i0) = -rho(t), rho the kite spectral density. The cusps t = 1 (threshold) and t = 9 (pseudo-threshold) are refused by name, as are points with |q_s(t + i0)| >= QMAX (a low-precision pass of route H sizes |q_s| first, so the refusal is cheap). The stretch 1 < t < 2 and the neighbourhoods of the cusps carry no independent reference here. MINKOWSKI REFERENCE VALUES (GATE-ONLY literals REF_MINKOWSKI / REF_RHO; the value path never reads them): at t = 4 and t = 16 the eps^0 top master from independent auxiliary-mass-flow runs at the same point in the same measure, the engine's own +i0 prescription (goal 60; the goal-40 twin agrees to 67.56 d at t = 4 and 67.37 d at t = 16, the certified count of each literal); at t = 12, 50, 100 the kite spectral density rho(w) = -Im I_kite(w + i0) computed by differential-equation transport and gated to 99.9 d against independent runs (Im only; Re has no independent reference at these points). Whenever a --point coincides with one of these five, the computed value is compared per component: digits = -log10 |computed - reference| / |reference|, capped at min(dps, certified), bar = cap - DIGIT_BAR_SLACK, and the SIGN of the imaginary part must agree (a wrong sign = the wrong i0 side, FAIL by name); below the bar the run exits 3, in the suite and in --point-only mode alike. At a physical point without a literal the value prints with its certified bound and no reference comparison. CONSTANTS: all classical -- pi, zeta(3), Cl2(2pi/3) = Im Li2(z3) -- computed at runtime; the five word coefficients {-2pi^2/3, -8, +4, -108 Cl2(2pi/3), -108} are the PSLQ-clean rationals of the blind fit. No interim / single-point-fitted constants anywhere in the construction. REFERENCE VALUES (independent evaluations, embedded as literals, GATE-ONLY -- the value path never reads them; two are blind-fit grid points and two were held out of the fit; two of them are printed in the paper, see PAPER DIGITS): t=-1 : 260-digit independent evaluation [blind-fit grid pt] t=1/2 : 260-digit independent evaluation [blind-fit grid pt] t=-3 : 160-digit auxiliary-mass-flow HELD-OUT run [NEVER in the fit] t=-10 : 130-digit auxiliary-mass-flow HELD-OUT run [NEVER in the fit] The two fit-grid points are still honest checks of THIS script (the coefficients are exact rationals, not per-point fits), but the headline is min(held-out). Each literal carries its provenance in the REF dict below: the auxiliary-mass-flow output file it was read from (name + sha256 prefix; the string is that run's eps^0 midpoint), the run's goal digits / eps order / thread count and the ball radius the run reported; the two 260-digit values are goal-250 midpoints each confirmed by a goal-280 twin run to the full 260 digits, the t=-3 value (goal 150) by a goal-120 twin to the twin's own length. PAPER DIGITS (checked on every default run): the paper's equal-mass kite section prints I_kite(-1) = -1.3317114... (its sample-value equation) and, at the held-out point, I_kite(-3) = -0.88550175...; the computed values must reproduce every printed digit (agreement to within one unit in the last printed decimal place) or the run exits 4. The same section quotes held-out agreement of 160 digits at t=-3 and 130 digits at t=-10, each count capped by the stored reference (161- and 131-digit strings; the last digit of each lies inside the ball its run reported); that claim is tested by --paper (dps 170), or at any --dps >= 162 for t=-3 and >= 132 for t=-10 (the requested precision must exceed the reference length, so that the count is capped by the reference and not by --dps), and is reported as "not testable at this dps" otherwise (the default dps 40 and --full's dps 100 included). The measured count is compared as printed (one decimal), the form the quoted count was read in. POSITIVE CONTROLS (computed at runtime at t=-1, independent of the references): 1607.01571 eq. integration_kernel_1 : (1/(i pi)) psi1^2/(W t) = 1 - 4 Ebar_{0;-1}(z3;-1;-q) with the sunrise Wronskian W = -12 pi i/(t(t-1)(t-9)) [Legendre + PF] eq. HPL_weight_1_2 depth 1 : G(1;y) = 3[Ebar_{1;0}(-1;1;-q) - Ebar_{1;0}(r6;1;-q)] eq. HPL_weight_1_2 depth 2 (2o=2) : G(1,1;y) = 9[Ebar_{0,1;-1,0;2}(...) 4-term combo] The depth-2 control exercises exactly the ELi code path used by the Ebar_{0,2;-2,0;2} word of the answer. NEGATIVE CONTROLS of the physical arm (each needs a --point t > 1 that carries a reference literal; each MUST exit 3): --sign-plant : the i0 sign flipped -- the polygon through the LOWER half plane and the conjugate boundary values give the complex conjugate; the reference gate FAILS by name on the sign of the imaginary part. --naive : the principal-branch ellipk nome at t + i delta fed to the same words (the naive per-piece t + i0); both components then disagree with the reference from the first digit and the gate FAILS by name. CERTIFICATION of every printed value (fail-closed): * nome branch gate (ellipk nome vs path-pinned hauptmodul nome); * N vs 2N doubling-agreement gate on each Ebar q-series word: the closed-form series length is a STARTING guess only; a word is accepted when |f(N)-f(2N)| < 10^-(dps+10), else N doubles by EXACT continuation of the same truncated double sum up to 16x the start, then the run refuses; * imaginary-part residual of the assembled real value; * the two agreement residuals propagated through the exact prefactors to a value-level l1 bound, printed on a "[certified]" line and compared (not asserted) against the tolerance. EXIT CODES (nothing prints PASS on a nonzero exit): 0 every check passed (reference gate, certified bounds, paper digits, positive controls; with --check also the two-precision rerun) 2 usage or domain refusal, named: t = 0, t = 1, t = 9 (cusps), |q| >= 0.90 (too close to a cusp, on either side of the threshold), incompatible flags, a control without a gated physical point, missing mpmath 3 reference gate or certification failure: a computed value agrees with its independent reference below the bar min(dps, reference digits) - 2, the imaginary part at a Minkowski reference point has the wrong sign, a certified l1 bound reaches the tolerance, a q-series or the nome Newton fails to converge, the nome branch gate or a positive control fails, the continued nome (route H vs route P, path independence), the Picard-Fuchs operator or the HPL boundary values are not certified -- the value cannot be certified and is not to be trusted 4 paper-digits mismatch: a printed digit of the paper is not reproduced, or (when testable) the quoted held-out digit counts are not reached 5 --check two-precision rerun unstable --mutate perturbs the -108 coefficient of the depth-2 word by a relative 1e-30 (a wrong closed form): every reference row then agrees to only ~30 digits, below the default bar of 38, and the run MUST exit 3. --sign-plant / --naive (physical arm, see NEGATIVE CONTROLS) MUST exit 3. REQUIREMENTS: python3 + mpmath ONLY. No network, no file reads, no computer-algebra system. mp.dps is set INSIDE main() after argparse (no module-level floats exist in this file). DEFAULT RUN (no flags): dps 40 -- positive controls, the four-point reference table with certified bounds, the paper-digits check. Single-thread mpmath; a few seconds on a laptop-class core (measured under ten seconds on one core of a heavily loaded build server). --full = the same suite at dps 100 (measured about twenty seconds under the same conditions); --paper = the same suite at dps 170, the tier that tests the paper's held-out digit counts (measured about half a minute on one core of a loaded shared host); --check reruns at dps+60. CLI: python3 kite-equal-evaluate.py # default: dps 40 suite python3 kite-equal-evaluate.py --full # the suite at dps 100 python3 kite-equal-evaluate.py --paper # the suite at dps 170: tests the paper's held-out digit counts python3 kite-equal-evaluate.py --dps 100 --check # two-precision rule: rerun at dps+60, diff python3 kite-equal-evaluate.py --point=-7/3 # add any Euclidean t<1, t!=0 to the table python3 kite-equal-evaluate.py --point=-7/3 --point-only # that point alone (batteries skipped) python3 kite-equal-evaluate.py --point 4 --point-only # PHYSICAL t > 1: I_kite(4 + i0) as Re / Im, gated against the vendored +i0 reference python3 kite-equal-evaluate.py --point 12 --point-only --dps 60 # Im gated against the spectral density; Re printed without reference python3 kite-equal-evaluate.py --point 2 --point-only # a physical point without a literal: value + certified bound, no reference python3 kite-equal-evaluate.py --point 4 --dps 60 # the physical point added to the default suite (a table row of two lines) python3 kite-equal-evaluate.py --mutate # control: MUST exit 3 python3 kite-equal-evaluate.py --point 4 --point-only --sign-plant # control (the i0 sign): MUST exit 3 python3 kite-equal-evaluate.py --point 4 --point-only --naive # control (the naive nome): MUST exit 3 (write --point=-7/3 with '=': a bare leading-dash value confuses argparse) PHYSICAL TIER WALLS (one process each, a shared 96-core host at loadavg 205.42, nice 10, wall clock by /usr/bin/time, --point-only): t = 4 at dps 30 6.7 s, dps 60 11.0 s, dps 90 17.7 s; at dps 60 t = 12 33.2 s, t = 16 29.1 s, t = 50 39.4 s, t = 100 59.4 s; peak resident set 28 MB. The cost grows with |q_s| (the q-series length) and with the number of polygon steps, not with t itself. CHANGELOG: 2026-07-05 -- certification hardening: (1) nome Newton RAISES on iteration exhaustion (was: silent fall-through returning a possibly non-converged q); path-tracking cap 12 -> NEWTON_TRACK_ITERS=20 (measured worst 8 used at dps<=100 over the reference set + 4 extras; need grows ~log2(dps) by quadratic convergence, so 20/80 cover any practical --dps). (2) Eb1/Eb2 answer words: the closed-form q-series length qseries_len() is a STARTING guess N0 only. Each word must pass an N vs 2N doubling-agreement gate |f(N)-f(2N)| < 10^-(dps+QGUARD=10); on failure N doubles by EXACT continuation of the same truncated double sum (identical terms jk<=N, plus more) up to NCAP_MULT=16 x N0, then RuntimeError naming word/t/|q|/residual/tol/N/cap. Accepted value is f(N) once gated, so healthy outputs are byte-identical to pre-edit; the per-point agreement residuals + the l1-propagated value-level bound are printed as "[certified]" lines. Calibration (instrumented run): worst residual over dps {30,60,100} x 8 points = 0.0 exactly (tail terms beyond N0 lie below the working-precision ulp), so the default path NEVER escalates and tol 10^-(dps+10) has unbounded measured headroom; a crippled N0 (mutation test) makes the residual finite and the gate escalate/raise. (3) _cert_line COMPARED, not asserted: the '[certified] ... l1 bound X < tol Y' relation was a HARDCODED '<' -- under an escalation sabotage at t=-7/3, dps=100 it printed '1.12e-110 < tol 1.0e-110', which is FALSE. The relation is now computed ('<' or '>='), and any bound >= tol is a gate failure: the full table still prints, then the run fails closed. Healthy outputs byte-identical. (4) gate on the reference table (loud = NONZERO EXIT, not just a printed table): every reference row must measure >= min(dps, ref digits) - DIGIT_BAR_SLACK(=2) digits; below-bar rows (e.g. a 1e-30 mutation of a reference string or of Cl2, which collapsed the table to ~29.6 d yet exited 0 before this change) now fail closed after the table prints. A '--check' FAIL verdict fails the same way. Slack 2 covers end-of-string rounding of the stored reference when dps ~ ref digits; healthy runs sit exactly AT the cap (digits are clipped there), so the bar has >= 2-digit measured headroom on every row. 2026-07-06 -- --point-only: plain point evaluation -- skips the always-on control battery (positive controls + the 4-point reference table) and evaluates ONLY the requested --point values. The VALUE PATH is untouched: nome branch gate, certified N/2N q-series agreement, the imaginary-residual check and the value-level l1 bound all still run per point, so the computed value is bit-identical to the same point in a battery run (printed here at full --dps; the [certified] line is byte-identical). Requires --point; mutually exclusive with --check. 2026-09-03 -- page release (this file): the research evaluator of row 13 ported into the kite bundle. Value path unchanged (every printed value line byte-identical to the research script at the same dps). Added: named exit codes 0/2/3/4/5 (the gate failures above used to surface as Python tracebacks, rc 1); the paper-digits check (exit 4) against the two printed values and the quoted held-out digit counts (70/71 then; 160/130 since 2026-09-05, below); --mutate control (relative 1e-30 on the -108 coefficient of the depth-2 word; MUST exit 3); --full (dps 100); default dps 100 -> 40 (the 1e-30 control is caught at 40: bar 38 vs ~30 d measured); line-buffered stdout with an immediate banner; positive controls now gated (each residual must lie below 10^-(dps+5)); public-register comments. 2026-09-05 -- held-out digit counts reprinted: the paper's quoted counts at t=-3 / t=-10 are 160 / 130 (were 70/71, the June blind fit's own residual at its working precision); the counts are those this script's gate prints at dps >= 162 / 132, each capped by the stored reference (161 / 131 digits; ball radii 2.64e-161 / 4.39e-131), 160.5 / 130.0 d measured. PAPER_HELDOUT_DIGITS follows; the claim is testable once the requested precision exceeds the reference length (dps >= count + 2) and the measured count is compared as printed (one decimal). NEW --paper = the suite at dps 170 (tests the claim); --full stays at dps 100 (there the claim is not testable; its help says so); default dps 40 untouched (the --mutate control is still caught there, exit 3). The REF dict now carries per-string provenance (output file name + sha256 prefix, goal digits, eps order, threads, ball radius; the higher-goal twin runs) and the t=-10 label says HELD-OUT like t=-3 (it never was in the fit; the table already printed HELD-OUT). Value path untouched: every value line at dps 40 / 100 / 170 / 300 byte-identical to the previous release (walls masked); only the t=-10 source label and the [paper] claim lines differ. 2026-09-09 -- the physical region t > 1 (--point t > 1; the arm of section 5b): the served closed form continued with the elliptic Feynman-i0 prescription -- the nome followed along t + i0 through the upper half plane by two independent routes (Newton path-tracking of the served hauptmodul series; Picard-Fuchs transport of the periods) that must agree, and at two polygon heights; the HPL pieces at their t + i0 boundary values (a runtime control against the served closed forms at t + i delta); the served words, coefficients and N/2N certification with the complex nome; the value printed as Re / Im. GATE-ONLY literals REF_MINKOWSKI (t = 4, 16: the eps^0 balls of independent auxiliary-mass-flow +i0 runs, goal 60, the goal-40 twin beside; certified 67 digits each) and REF_RHO (t = 12, 50, 100: the kite spectral density by differential-equation transport, 99 digits; Im only), compared per component with the served bar rule and a sign check on the imaginary part (exit 3). Refused by name: t = 1 (threshold cusp) and t = 9 (pseudo-threshold cusp), |q_s(t + i0)| >= QMAX (sized by a low-precision pass first). Controls --sign-plant / --naive (exit 3). The Euclidean path is untouched: every served function is byte-identical and every Euclidean tier prints what it printed before (stamps and walls masked); eval_kite still raises at t >= 1 and is reached only for t < 1 (the dispatch routes t > 1 to eval_kite_physical and refuses t = 1 by name before either). """ import argparse import sys import time sys.stdout.reconfigure(line_buffering=True) try: import mpmath as mp except ImportError: print("ERROR: this script needs mpmath (pip install mpmath)") sys.exit(2) EXIT_OK, EXIT_USAGE, EXIT_GATE, EXIT_PAPER, EXIT_CHECK = 0, 2, 3, 4, 5 QMAX = "0.90" # refuse |q| beyond this (cusp neighbourhoods t->1^-, 0, -inf) GUARD = 25 # working-precision guard digits on top of --dps QGUARD = 10 # N/2N agreement gate: Eb words must agree to 10^-(dps+QGUARD) NCAP_MULT = 16 # q-series escalation cap: N may grow to NCAP_MULT*N0, then raise NEWTON_TRACK_ITERS = 20 # tracking-step Newton cap (measured worst 8 used, dps<=100) DIGIT_BAR_SLACK = 2 # gate bar: reference rows must reach min(dps, ref digits) - this CONTROL_GUARD = 5 # positive-control residuals must lie below 10^-(dps+CONTROL_GUARD) DPS_DEFAULT = 40 # default run (measured: a few seconds on one core) DPS_FULL = 100 # --full DPS_PAPER = 170 # --paper: tests the paper's held-out digit counts (PAPER_HELDOUT_DIGITS) MUTATE = False # --mutate control: perturb the -108 depth-2 coefficient SIGN_PLANT = False # --sign-plant control (physical arm): the i0 sign flipped NAIVE = False # --naive control (physical arm): the principal-branch ellipk nome at t + i delta PHYS_T_BASE = "1/2" # Euclidean base point of the continuation (branch-gated there by the served path) PHYS_HEIGHTS = ("1.5", "2.5") # the two polygon heights in the upper half t-plane (path independence) PHYS_SING = (0, 1, 9) # singular fibres of q_s(t) on the real axis PHYS_PRESCAN_DPS = 20 # a cheap route-H pass sizes |q_s(t + i0)| before the working-precision run PHYS_QMAX_PATH = "0.95" # refuse by name when |q_s| on the continuation path reaches this (a cusp neighbourhood) # --------------------------------------------------------------------------- # Independent reference values (GATE ONLY). The 're' strings of the reference # set; the same values are quoted in the paper. # 'fitgrid': point was on the 12-point grid of the original blind coefficient # fit (coeffs since snapped to exact rationals by PSLQ); held-out points never. # 'prov': where the string was read from -- the auxiliary-mass-flow output # file (name + sha256 prefix) whose eps^0 midpoint it is, the run's goal # digits / eps order / thread count, and the ball radius the run reported # (the count a comparison can reach is capped by it); 'twin': a second run # at another goal that agrees with the string to the digits stated. # --------------------------------------------------------------------------- REF = { "-1": dict(fitgrid=True, src="independent evaluation, 260d", prov=dict(out="kite_tm1_g250_o8.json", sha256_16="374c93f1d0133e25", goal=250, eps_order=8, n_thread=8, ball_radius="7.88e-261"), twin=dict(out="kite_tm1_g280_o10.json", sha256_16="4e2ed5be446fa6c6", goal=280, eps_order=10, n_thread=8, ball_radius="2.37e-290", agrees_digits=260), val=( "-1.3317114414221095767985229828526649733648343335308848237269336649" "1108784462599027592272079304001792909391595391945255013541386706180" "6055294303664935019443519432524729642507788448614715717953981890094" "58336930190898228010435491752832143899855987475388371821396 92")), "1/2": dict(fitgrid=True, src="independent evaluation, 260d", prov=dict(out="kite_t1o2_g250_o8.json", sha256_16="92355e8e7d1a0c6f", goal=250, eps_order=8, n_thread=8, ball_radius="4.49e-261"), twin=dict(out="kite_t1o2_g280_o10.json", sha256_16="7ad1cd7ca184ca43", goal=280, eps_order=10, n_thread=8, ball_radius="3.87e-290", agrees_digits=260), val=( "-2.4590795483625630714549568899983825289812658026439695939882582373" "3313925725890241514582046918386357120438923712120858378987359390241" "6716346830685292158534682279185410348088271538953596923126694985423" "52570194352973493454320976845610105370303754626054826818363 18")), "-3": dict(fitgrid=False, src="AMFlow HELD-OUT run, 160d", prov=dict(out="kite_tm3_g150_o10.json", sha256_16="9eb91dbcff185d11", goal=150, eps_order=10, n_thread=8, ball_radius="2.64e-161"), twin=dict(out="kite_tm3_g120_o8.json", sha256_16="b7b82f19ccdcb18b", goal=120, eps_order=8, n_thread=8, ball_radius="1.16e-131", agrees_digits=130), val=( "-0.8855017520635901794137037145612123046455203550272134813868766996" "0394745765710810346613750186209565844353473166904539412738131863621" "15845371753734435023509289109")), "-10": dict(fitgrid=False, src="AMFlow HELD-OUT run, 130d", prov=dict(out="kite_tm10_g120_o8.json", sha256_16="ec4e89a3ccd7e111", goal=120, eps_order=8, n_thread=8, ball_radius="4.39e-131"), val=( "-0.4351420364817559093449603607174377351193979595611879541081485289" "2626591229233556610791762044513268240144194917543178943425962798 66")), } # Digits printed in the paper's equal-mass kite section (transcribed; the # trailing "..." of the paper is the truncation). Checked on every run that # evaluates these points: agreement to within one unit in the last printed # decimal place, else exit 4. PAPER_PRINTED = { "-1": ("-1.3317114", "sample-value equation, Euclidean point t=-1"), "-3": ("-0.88550175", "held-out point t=-3, sentence after that equation"), } # Held-out agreement the same section quotes: the five-word closed form vs # the independent oracle, 160 digits at t=-3 and 130 at t=-10, each count # capped by the oracle value's own precision (goal-150 / goal-120 runs: # 161- / 131-digit strings whose last digit lies inside the reported ball). # Testable only when dps - DIGIT_BAR_SLACK reaches the quoted count, i.e. # when the requested precision exceeds the reference length and the measured # count is capped by the reference, not by --dps: --paper, or --dps >= # 162 (t=-3) / 132 (t=-10). The count is compared as printed (one decimal). PAPER_HELDOUT_DIGITS = {"-3": 160, "-10": 130} # --------------------------------------------------------------------------- # Minkowski reference values (GATE ONLY; the physical arm, see PHYSICAL REGION # and MINKOWSKI REFERENCE VALUES in the header). REF_MINKOWSKI: at t = 4 and # t = 16 the eps^0 top master, real and imaginary part, from an independent # auxiliary-mass-flow run at goal 60 with the engine's own +i0 prescription # ('prov': output file name + sha256 prefix, goal / eps order / threads, the # two ball radii; 'twin': the goal-40 run at the same point and the digits to # which the two midpoints agree per component; 'certified' = the whole digits # of the smaller of the two, the count a comparison can reach). REF_RHO: at # w = 12, 50, 100 the kite spectral density rho(w) = -Im I_kite(w + i0) by # differential-equation transport, gated against independent runs to the # digits in 'prov' ('certified' = the whole digits of that agreement); Im only. # --------------------------------------------------------------------------- REF_MINKOWSKI = { "4": dict(src="AMFlow +i0 run, goal 60 (goal-40 twin beside)", prov=dict(out="kite_t4_g60.json", sha256_16="945f0bb541400251", goal=60, eps_order=8, n_thread=4, ball_radius_re="1.35e-111", ball_radius_im="2.89e-110"), twin=dict(out="kite_t4_g40.json", sha256_16="11c99eb8418e8eba", goal=40, eps_order=8, n_thread=4, ball_radius_re="1.68e-111", ball_radius_im="2.26e-110", agrees_digits_re=67.56273308734428, agrees_digits_im=72.40089419688749), certified=67, re=( "0.8061973750252974703948373709852941613333886499683351611454229889" "2038615082033855248840525187915744050688914700"), im=( "-2.177222842384698915978940960458172768489490190872105517655865043" "4626760802287908229289520550041673036849171975")), "16": dict(src="AMFlow +i0 run, goal 60 (goal-40 twin beside)", prov=dict(out="kite_t16_g60.json", sha256_16="aac92577b2f78f7d", goal=60, eps_order=8, n_thread=4, ball_radius_re="3.06e-111", ball_radius_im="2.29e-111"), twin=dict(out="kite_t16_g40.json", sha256_16="f4ac177b1795d717", goal=40, eps_order=8, n_thread=4, ball_radius_re="1.48e-111", ball_radius_im="3.99e-111", agrees_digits_re=67.37025979131192, agrees_digits_im=70.71294640053449), certified=67, re=( "0.5177886228653184647557313834150044434428991012059652005629456761" "6428978213270629529281747551027772184718952020"), im=( "-0.229337613164037052112195452065059725578426386712430870882468427" "80662091147347445471364300083218367540303370729")), } REF_RHO = { "12": dict(src="spectral density rho(w) by DE transport, gated 99.95 d against an independent AMFlow run", prov=dict(note_sha256_16="baa6b9903b63932f", agreement_digits="99.95", string_digits=99), certified=99, minus_rho=( "-0.391181609082472721402976623249127689606443541627541768278400581" "811628178499673193969248045186550608")), "50": dict(src="spectral density rho(w) by DE transport, gated 99.90 d against an independent AMFlow run", prov=dict(note_sha256_16="baa6b9903b63932f", agreement_digits="99.90", string_digits=98), certified=99, minus_rho=( "-0.028035070670697914545111757600619121987827831570764968417371278" "787660040466139149848126988032070568")), "100": dict(src="spectral density rho(w) by DE transport, gated 99.89 d against an independent AMFlow run", prov=dict(note_sha256_16="baa6b9903b63932f", agreement_digits="99.89", string_digits=100), certified=99, minus_rho=( "-0.007791679551459260519749415937762143296843421082607161659729934" "316087225684619465656805466749899847828")), } def _refval(key): return mp.mpf(REF[key]["val"].replace(" ", "")) def _refdigits(key): return sum(c.isdigit() for c in REF[key]["val"]) # --------------------------------------------------------------------------- # Roots of unity / classical constants (built at CURRENT precision -- never at # import time; import-time mpf runs at dps=15 and silently caps everything). # --------------------------------------------------------------------------- def roots(): Iu = mp.mpc(0, 1) r3 = mp.exp(2 * mp.pi * Iu / 3) r6 = mp.exp(2 * mp.pi * Iu / 6) cl2 = ((mp.polylog(2, r3) - mp.polylog(2, 1 / r3)) / (2 * Iu)).real # Cl2(2pi/3) return Iu, r3, r6, cl2 # --------------------------------------------------------------------------- # 1. Curve periods / geometric nome (AW eq. def_roots / def_periods) # --------------------------------------------------------------------------- def curve_nome(t): """(psi1, psi2, tau, q) on Dtilde=(t-1)^3(t-9), q=e^{i pi tau}. Euclidean t<1.""" Iu = mp.mpc(0, 1) tt = mp.mpc(t) Dt = (tt - 1) ** 3 * (tt - 9) sq = mp.sqrt(Dt) e1 = (-tt ** 2 + 6 * tt + 3 + 3 * sq) / 24 e2 = (-tt ** 2 + 6 * tt + 3 - 3 * sq) / 24 e3 = (2 * tt ** 2 - 12 * tt - 6) / 24 k2 = (e3 - e2) / (e1 - e2) kp2 = (e1 - e3) / (e1 - e2) Dq = Dt ** (mp.mpf(1) / 4) psi1 = 4 / Dq * mp.ellipk(k2) psi2 = 4 * Iu / Dq * mp.ellipk(kp2) tau = psi2 / psi1 if mp.im(tau) < 0: tau = -tau psi2 = -psi2 return psi1, psi2, tau, mp.exp(Iu * mp.pi * tau) # --------------------------------------------------------------------------- # 2. Level-6 eta-quotient hauptmodul t(q_s) (the same construction as # sunrise-evaluate.py in the sunrise bundle of this site; exact integer # coefficients, generated at runtime) + Newton path-pinned inversion = the # BRANCH GATE. Dictionary: q_AW = -q_sunrise. # --------------------------------------------------------------------------- def _euler_sparse(N): sp, k = [(0, 1)], 1 while k * (3 * k - 1) // 2 <= N: s = -1 if k % 2 else 1 for e in (k * (3 * k - 1) // 2, k * (3 * k + 1) // 2): if e <= N: sp.append((e, s)) k += 1 return sorted(sp) def _mul_sparse(a, spr, N): r = [mp.mpf(0)] * (N + 1) for e, s in spr: for i in range(0, N + 1 - e): r[i + e] += s * a[i] return r def _div_sparse(a, spr, N): r = [mp.mpf(0)] * (N + 1) for m in range(N + 1): acc = a[m] for e, s in spr[1:]: if e <= m: acc -= s * r[m - e] r[m] = acc return r _HAUPT_CACHE = {} def haupt_series(N): """t(q_s) = 9 q_s prod (1-q_s^{6n})^8 (1-q_s^n)^4 (1-q_s^{2n})^-8 (1-q_s^{3n})^-4.""" key = (mp.mp.dps, N) if key not in _HAUPT_CACHE: t = [mp.mpf(0)] * (N + 1) t[0] = mp.mpf(1) for d, r in [(6, 8), (1, 4), (2, -8), (3, -4)]: spd = [(d * e, s) for e, s in _euler_sparse(N // d)] for _ in range(abs(r)): t = _mul_sparse(t, spd, N) if r > 0 else _div_sparse(t, spd, N) _HAUPT_CACHE[key] = [mp.mpf(0)] + [9 * c for c in t[:-1]] return _HAUPT_CACHE[key] def _horner(c, q): acc = mp.mpf(0) for ck in reversed(c): acc = acc * q + ck return acc def _newton(q, tt, tser, dser, iters=80): """Newton solve of the hauptmodul t(q)=tt. RAISES on exhaustion (was a silent fall-through returning a possibly non-converged q).""" tol = mp.mpf(10) ** (-mp.mp.dps + 5) dq = None for _ in range(iters): dq = (_horner(tser, q) - tt) / _horner(dser, q) q -= dq if abs(dq) < tol: return q raise RuntimeError( f"nome Newton EXHAUSTED at hauptmodul target t={mp.nstr(mp.mpf(tt), 8)}: " f"|dq| = {mp.nstr(abs(dq), 3)} after {iters} iterations (tol " f"{mp.nstr(tol, 3)}); cusp-connected nome NOT certified -- refusing") def haupt_len(qabs, dps): """Series length for the eta quotient: coefficients grow ~10^{1.5 sqrt(n)} (partition-type), so solve n*(-log10 qa) - 1.5 sqrt(n) >= dps+15 iteratively.""" slope = -mp.log10(qabs) n = int((dps + 15) / slope) + 40 for _ in range(6): n2 = int((dps + 15 + mp.mpf(1.5) * mp.sqrt(n)) / slope) + 40 if n2 <= n: break n = n2 return max(n, 120) def nome_pathed(t, qabs_hint): """Cusp-connected sunrise nome q_s(t) by Newton on the hauptmodul, pinned to the t=0 cusp: direct seed t/9 for -3 < t < 1 (inside the convergence disk of the seed; no singular point between), real-axis geometric tracking from the t=-3 base point for t <= -3 (no singular points on t<0). Euclidean-only.""" N = haupt_len(min(qabs_hint + mp.mpf("0.05"), mp.mpf("0.95")), mp.mp.dps) tser = haupt_series(N) dser = [k * tser[k] for k in range(1, len(tser))] tval = mp.mpf(t) if tval > -3: return _newton(tval / 9, tval, tser, dser) q = _newton(mp.mpf(-3) / 9, mp.mpf(-3), tser, dser) cur = mp.mpf(-3) while cur > tval: cur = max(tval, cur * mp.mpf("1.10")) q = _newton(q, cur, tser, dser, iters=NEWTON_TRACK_ITERS) return _newton(q, tval, tser, dser) # --------------------------------------------------------------------------- # 3. ELi / Ebar q-series (AW 1607.01571 sec. 5; the same code as the original # 180-digit-validated evaluation) # --------------------------------------------------------------------------- def ELi_nm(n, m, x, y, q, N): """ELi_{n;m}(x;y;q) = sum_{j,k>=1} x^j/j^n y^k/k^m q^{jk}, truncated jk<=N.""" s = mp.mpc(0) for j in range(1, N + 1): xj = x ** j / mp.mpf(j) ** n if n != 0 else x ** j qj = q ** j kmax = N // j if kmax < 1: break inner = mp.mpc(0) for k in range(1, kmax + 1): inner += (y ** k / mp.mpf(k) ** m if m != 0 else y ** k) * qj ** k s += xj * inner return s def _sgn(n): return 1 if (n % 2 == 0) else -1 def _cn(n): return mp.mpc(1) if (n % 2 == 0) else mp.mpc(0, 1) def Ebar_1(n, m, x, y, q, N): """Depth-1 Ebar (AW eq. def_Ebar_weight_1).""" a = ELi_nm(n, m, x, y, q, N) b = ELi_nm(n, m, 1 / x, 1 / y, q, N) return (_cn(n + m) / mp.mpc(0, 1)) * (a - _sgn(n + m) * b) def _channel_table(n, m, x, y, N): """A[e] = sum_{jk=e} x^j/j^n y^k/k^m, e = 1..N-1.""" A = [mp.mpc(0)] * N for j in range(1, N): xj = (x ** j / mp.mpf(j) ** n) if n else x ** j for k in range(1, N // j + 1): e = j * k if e >= N: break A[e] += xj * ((y ** k / mp.mpf(k) ** m) if m else y ** k) return A def ELi2_o2(n1, n2, m1, m2, x1, x2, y1, y2, q, N): """ELi_{n1,n2;m1,m2;2} = sum_{e>=2} (q^e/e) sum_{e1+e2=e} A1[e1] A2[e2] (the 2o1=2 index = one dq'/q' integration = 1/(j1k1+j2k2) on the double sum).""" A1 = _channel_table(n1, m1, x1, y1, N) A2 = _channel_table(n2, m2, x2, y2, N) s = mp.mpc(0) for e in range(2, N): c = mp.mpc(0) for e1 in range(1, e): if A1[e1] != 0: c += A1[e1] * A2[e - e1] if c != 0: s += c * q ** e / mp.mpf(e) return s def Ebar_2(n1, n2, m1, m2, o1, x1, x2, y1, y2, q, N): """Depth-2 Ebar, index 2o1=2: 2^2 ELi combination (AW eqs. Ebar_* rules).""" assert o1 == 1 Iu = mp.mpc(0, 1) pre1, pre2 = _cn(n1 + m1) / Iu, _cn(n2 + m2) / Iu s = mp.mpc(0) for t1 in (0, 1): for t2 in (0, 1): X1, Y1 = (x1, y1) if t1 == 0 else (1 / x1, 1 / y1) X2, Y2 = (x2, y2) if t2 == 0 else (1 / x2, 1 / y2) coef = pre1 * ((-_sgn(n1 + m1)) ** t1) * pre2 * ((-_sgn(n2 + m2)) ** t2) s += coef * ELi2_o2(n1, n2, m1, m2, X1, X2, Y1, Y2, q, N) return s def qseries_len(qabs, dps): """ELi truncation jk<=N: channel coefficients grow polynomially (divisor-type, |A[e]| <~ e^3 here), so N(-log10 qa) - 3 log10 N >= dps+12, iterated.""" slope = -mp.log10(qabs) n = int((dps + 12) / slope) + 30 for _ in range(6): n2 = int((dps + 12 + 3 * mp.log10(n)) / slope) + 30 if n2 <= n: break n = n2 return max(n, 80) # --------------------------------------------------------------------------- # 4. Weight<=3 {0,1}-HPLs, closed forms (validated against quadrature) # --------------------------------------------------------------------------- def G011(y): return (-mp.polylog(3, 1 - y) + mp.log(1 - y) * mp.polylog(2, 1 - y) + mp.log(y) * mp.log(1 - y) ** 2 / 2 + mp.zeta(3)) def G101(y): # shuffle: G(1)G(0,1) = G(1,0,1) + 2 G(0,1,1); G(0,1) = -Li2 return mp.log(1 - y) * (-mp.polylog(2, y)) - 2 * G011(y) # --------------------------------------------------------------------------- # 4b. N/2N doubling-agreement gate for the ELi/Ebar q-series. # qseries_len() is a STARTING guess only; the accepted truncation must be # CERTIFIED by two-successive-depth agreement or the script raises. # --------------------------------------------------------------------------- def _ebar_certified(label, t, qa, compute, N0): """Certify the q-series truncation of one Ebar word by N vs 2N agreement. Accept f(N) only when |f(N) - f(2N)| < 10^-(dps+QGUARD) (dps = user digits = mp.mp.dps - GUARD). On failure N doubles by EXACT continuation of the same truncated double sum (the jk<=N terms are a strict subset of the jk<=2N terms; nothing is recomputed differently), up to NCAP_MULT*N0; at the cap raise RuntimeError with the measured residual. Returns (accepted f(N), agreement residual, accepted N).""" qtol = mp.mpf(10) ** (-(mp.mp.dps - GUARD + QGUARD)) ncap = NCAP_MULT * N0 N = N0 a = compute(N) while True: b = compute(2 * N) resid = abs(a - b) if resid < qtol: return a, resid, N if 2 * N > ncap: raise RuntimeError( f"{label} q-series NOT CONVERGED at t={mp.nstr(mp.mpf(t), 8)} " f"(|q|={float(qa):.3f}): N/2N agreement residual " f"{mp.nstr(resid, 3)} >= tol {mp.nstr(qtol, 3)} at N={N} " f"(cap {ncap} = {NCAP_MULT}*N0, N0={N0}); refusing to print an " "uncertified value") N, a = 2 * N, b # --------------------------------------------------------------------------- # 5. The equal-mass kite evaluator # --------------------------------------------------------------------------- def eval_kite(t, verbose_checks=None): """I_kite(t) at Euclidean t<1, t!=0 (current mp.dps). Returns (value, info).""" tval = mp.mpf(t) if tval >= 1: raise ValueError( "t >= 1 is the physical region: the closed form is Euclidean-derived; " "continuation needs the BSW elliptic Feynman-i0 prescription (not the " "naive per-piece t+i0) -- refused, see the header") if tval == 0: raise ValueError("t = 0 is the cusp (1/(4t) prefactor); refused") psi1, psi2, tau, q = curve_nome(tval) qa = abs(q) if qa >= mp.mpf(QMAX): raise ValueError( f"|q(t)| = {float(qa):.3f} >= QMAX = {QMAX}: too close to a cusp " "(t -> 1^- or t -> -infinity); outside this script's series domain") # --- BRANCH GATE: ellipk nome must equal the cusp-connected path-pinned # hauptmodul nome via q_AW = -q_sunrise, or we abort. q_s = nome_pathed(tval, qa) branch_err = abs(q_s - (-q)) tol = mp.mpf(10) ** (-(mp.mp.dps - GUARD + 10)) if not branch_err < tol * max(qa, mp.mpf(1)): raise ArithmeticError( f"nome BRANCH GATE FAILED at t={t}: |q_pathed - (-q_ellipk)| = " f"{mp.nstr(branch_err, 3)} (tol {mp.nstr(tol, 3)}). The ellipk nome is " "off the cusp-connected branch here; do not trust this point.") mq = -q # AW Ebar argument (-q); equals +q_sunrise N = qseries_len(qa, mp.mp.dps) # STARTING guess N0 only Iu, r3, r6, cl2 = roots() y = mp.mpc(tval) g1 = mp.log(1 - y) g011 = G011(y) g101 = G101(y) Eb1, eb1_resid, N1 = _ebar_certified( "Eb1 = Ebar_{1;-1}(z3;1;-q)", tval, qa, lambda NN: Ebar_1(1, -1, r3, mp.mpc(1), mq, NN), N) Eb2, eb2_resid, N2 = _ebar_certified( "Eb2 = Ebar_{0,2;-2,0;2}(z3,z3;1,-1;-q)", tval, qa, lambda NN: Ebar_2(0, 2, -2, 0, 1, r3, r3, mp.mpc(1), mp.mpc(-1), mq, NN), N) # the five PSLQ-clean coefficients; --mutate perturbs the last one (the # depth-2 word) by a relative 1e-30 -- a wrong closed form, to be caught # by the reference gate below (rc 3) c_eb2 = mp.mpf(-108) * (1 + mp.mpf("1e-30")) if MUTATE else mp.mpf(-108) val = (-2 * mp.pi ** 2 / 3 * g1 - 8 * g011 + 4 * g101 - 108 * cl2 * Eb1 + c_eb2 * Eb2) / (4 * mp.mpc(tval)) im_resid = abs(val.imag) / max(abs(val), mp.mpf(1)) if not im_resid < mp.mpf(10) ** (-(mp.mp.dps - GUARD + 8)): raise ArithmeticError( f"imaginary residual {mp.nstr(im_resid, 3)} too large at t={t}: " "branch inconsistency in the HPL/eMPL pieces") # value-level propagation of the two agreement residuals through the exact # prefactors (l1): |d I_kite| <= (|108 cl2| r1 + |108| r2) / |4t| err_val = (108 * cl2 * eb1_resid + 108 * eb2_resid) / (4 * abs(mp.mpf(tval))) qtol = mp.mpf(10) ** (-(mp.mp.dps - GUARD + QGUARD)) info = dict(qabs=qa, N=N, branch_err=branch_err, im_resid=im_resid, eb1_resid=eb1_resid, eb2_resid=eb2_resid, N1=N1, N2=N2, err_val=err_val, qtol=qtol) if verbose_checks is not None: verbose_checks.update(info) return val.real, info # --------------------------------------------------------------------------- # 5b. The physical region t > 1: the elliptic Feynman-i0 continuation # (PHYSICAL REGION in the header). Nothing of sections 1-5 is edited: the # served series, Newton, words and certification are CALLED with the # continued (complex) nome. Route H = Newton path-tracking of the served # hauptmodul along a polygon in the upper half t-plane; route P = Taylor # transport of the periods with the Picard-Fuchs operator # P psi'' + Q psi' + R psi = 0, P = t^3 - 10t^2 + 9t, Q = 3t^2 - 20t + 9, # R = t - 3, # verified at runtime on the served ellipk periods. # --------------------------------------------------------------------------- def _phys_polygon(t_target, h, sign=1): """Polygon t_b -> t_b + i h -> t_target + i h -> t_target (sign = -1: the lower half plane).""" tb = _parse_t(PHYS_T_BASE) tt = mp.mpf(t_target) ih = mp.mpc(0, sign * h) return [mp.mpc(tb), tb + ih, tt + ih, mp.mpc(tt)] def _phys_dist_sing(z): return min(abs(z - s) for s in PHYS_SING) def _phys_walk(path, step_frac, step_cap): """The successive points along the polygon with steps <= step_frac x (distance to {0,1,9}).""" cur = path[0] pts = [cur] for nxt in path[1:]: seg = nxt - cur for s0 in PHYS_SING: if abs(seg) > 0: lam = ((s0 - cur) * mp.conj(seg)).real / abs(seg) ** 2 if 0 <= lam <= 1 and abs(cur + lam * seg - s0) < mp.mpf("1e-12"): raise ValueError(f"path segment {mp.nstr(cur, 6)} -> {mp.nstr(nxt, 6)} passes " f"through the singular point t={s0}") while abs(nxt - cur) > 0: d = min(step_frac * _phys_dist_sing(cur), step_cap) if abs(nxt - cur) <= d: cur = nxt else: cur = cur + (nxt - cur) / abs(nxt - cur) * d pts.append(cur) return pts def _phys_route_H(path, step_frac=None, step_cap=None): """Route H: the served hauptmodul inverse tracked along the polygon (served series + served Newton).""" step_frac = mp.mpf("0.2") if step_frac is None else step_frac step_cap = mp.mpf("0.5") if step_cap is None else step_cap pts = _phys_walk(path, step_frac, step_cap) t0 = pts[0] assert t0.imag == 0 and mp.mpf(-3) < t0.real < 1 qa_sized = mp.mpf("0.5") N = haupt_len(qa_sized, mp.mp.dps) tser = haupt_series(N) dser = [k * tser[k] for k in range(1, len(tser))] q = _newton(t0.real / 9, t0.real, tser, dser) qmax = abs(q) qcap = mp.mpf(PHYS_QMAX_PATH) try: for tt in pts[1:]: q = _newton(q, tt, tser, dser, iters=NEWTON_TRACK_ITERS) qa = abs(q) qmax = max(qmax, qa) if qa >= qcap: raise ValueError( f"|q_s| on the continuation path reached {float(qa):.3f} >= {PHYS_QMAX_PATH} near " f"t = {mp.nstr(tt, 8)}: the target lies in a cusp neighbourhood (t -> 1 or t -> 9); refused") if qa + mp.mpf("0.05") > qa_sized: qa_sized = min(qa + mp.mpf("0.12"), mp.mpf("0.97")) N = haupt_len(qa_sized, mp.mp.dps) tser = haupt_series(N) dser = [k * tser[k] for k in range(1, len(tser))] q = _newton(q, tt, tser, dser, iters=NEWTON_TRACK_ITERS) q = _newton(q, pts[-1], tser, dser) except TypeError: # the served Newton's exhaustion message renders a real t; on a complex path point that # rendering itself raises -- re-raise in the served form (cannot certify) raise RuntimeError( f"nome Newton EXHAUSTED on the continuation path (|q_s| {float(qmax):.3f} so far); " "cusp-connected nome NOT certified -- refusing") resid = abs(_horner(tser, q) - pts[-1]) return q, dict(n_steps=len(pts) - 1, qmax_on_path=qmax, N_series=N, resid=resid) def _phys_pf_residual(t): """|L psi| / scale for the two served ellipk periods at a Euclidean t (mp.diff).""" f1 = lambda x: curve_nome(x)[0] f2 = lambda x: curve_nome(x)[1] out = [] for f in (f1, f2): p0, p1, p2 = f(t), mp.diff(f, t, 1), mp.diff(f, t, 2) L = t * (t - 1) * (t - 9) * p2 + (3 * t ** 2 - 20 * t + 9) * p1 + (t - 3) * p0 out.append(abs(L) / max(abs(p0), 1)) return out def _phys_taylor_coeffs(c, a0, a1, Nord): """Taylor coefficients of a Picard-Fuchs solution about c with psi(c)=a0, psi'(c)=a1.""" p = [c ** 3 - 10 * c ** 2 + 9 * c, 3 * c ** 2 - 20 * c + 9, 3 * c - 10, mp.mpf(1)] q = [3 * c ** 2 - 20 * c + 9, 6 * c - 20, mp.mpf(3)] r = [c - 3, mp.mpf(1)] a = [mp.mpc(a0), mp.mpc(a1)] + [mp.mpc(0)] * Nord for n in range(0, Nord): s = p[1] * (n + 1) * n * a[n + 1] + p[2] * n * (n - 1) * a[n] + p[3] * (n - 1) * (n - 2) * (a[n - 1] if n >= 1 else 0) s += q[0] * (n + 1) * a[n + 1] + q[1] * n * a[n] + q[2] * (n - 1) * (a[n - 1] if n >= 1 else 0) s += r[0] * a[n] + r[1] * (a[n - 1] if n >= 1 else 0) a[n + 2] = -s / (p[0] * (n + 2) * (n + 1)) return a def _phys_taylor_eval(a, u): v = mp.mpc(0) for n in range(len(a) - 1, -1, -1): v = v * u + a[n] dv = mp.mpc(0) for n in range(len(a) - 1, 0, -1): dv = dv * u + n * a[n] return v, dv def _phys_route_P(path, step_frac=None, step_cap=None): """Route P: the periods (psi1, psi2) transported along the polygon; returns -q_P (the served dictionary).""" step_frac = mp.mpf("0.3") if step_frac is None else step_frac step_cap = mp.mpf("0.6") if step_cap is None else step_cap pts = _phys_walk(path, step_frac, step_cap) t0 = pts[0].real assert t0 < 1 and t0 != 0 psi1, psi2, tau, q = curve_nome(t0) f1 = lambda x: curve_nome(x)[0] f2 = lambda x: curve_nome(x)[1] d1, d2 = mp.diff(f1, t0, 1), mp.diff(f2, t0, 1) Nord = int(2.4 * mp.mp.dps) + 20 state = [(mp.mpc(psi1), mp.mpc(d1)), (mp.mpc(psi2), mp.mpc(d2))] cur = pts[0] for tt in pts[1:]: u = tt - cur new = [] for (v, dv) in state: a = _phys_taylor_coeffs(cur, v, dv, Nord) new.append(_phys_taylor_eval(a, u)) state = new cur = tt (p1, _), (p2, _) = state tauP = p2 / p1 qP = mp.exp(mp.mpc(0, 1) * mp.pi * tauP) return -qP, dict(n_steps=len(pts) - 1, tau=tauP, psi1=p1, psi2=p2, Nord=Nord, tau0=tau) def _phys_hpl_boundary(t, sign=1): """G(1), G(0,1,1), G(1,0,1) at y = t + i sign 0, t > 1: the boundary values of the served closed forms.""" t = mp.mpf(t) Iu = mp.mpc(0, 1) Lg = mp.log(t - 1) - sign * Iu * mp.pi li2_1mt = mp.polylog(2, 1 - t) li3_1mt = mp.polylog(3, 1 - t) li2_t = mp.pi ** 2 / 3 - mp.log(t) ** 2 / 2 - mp.polylog(2, 1 / t) + sign * Iu * mp.pi * mp.log(t) g1 = Lg g011 = -li3_1mt + Lg * li2_1mt + mp.log(t) * Lg ** 2 / 2 + mp.zeta(3) g101 = Lg * (-li2_t) - 2 * g011 return g1, g011, g101 def _phys_hpl_control(t, sign=1): """The served G011 / G101 at y = t + i delta vs the boundary formulas: (digits G011, digits G101, -log10 delta).""" delta = mp.mpf(10) ** (-(mp.mp.dps - GUARD) // 2) y = mp.mpc(mp.mpf(t), sign * delta) g1, g011, g101 = _phys_hpl_boundary(t, sign) a, b = G011(y), G101(y) return _digits_c(a, g011), _digits_c(b, g101), float(-mp.log10(delta)) def _phys_prescan_qabs(t, sign=1): """|q_s(t + i0)| from a low-precision route-H pass (the cusp neighbourhoods are refused before the full run).""" saved = mp.mp.dps mp.mp.dps = PHYS_PRESCAN_DPS + GUARD try: q, _ = _phys_route_H(_phys_polygon(mp.mpf(t), mp.mpf(PHYS_HEIGHTS[0]), sign)) return abs(q) finally: mp.mp.dps = saved def eval_kite_physical(t, sign=1, naive=False): """I_kite(t + i0) at physical t > 1, t != 9 (current mp.dps): the continued nome (route H, certified against route P and a second path height), the HPL boundary values (certified against the served closed forms at t + i delta), the served words and N/2N certification with the complex nome. Returns (complex value, info). sign = -1: the --sign-plant control (t - i0). naive = True: the --naive control.""" tval = mp.mpf(t) if not tval > 1: raise ValueError("t <= 1 is not the physical region (t < 1 is the served Euclidean path; t = 1 the threshold cusp)") if tval == 9: raise ValueError("t = 9 is the pseudo-threshold cusp (|q_s| -> 1); refused") dps = mp.mp.dps - GUARD qa_pre = _phys_prescan_qabs(tval, sign) if qa_pre >= mp.mpf(QMAX): raise ValueError( f"|q_s(t + i0)| = {float(qa_pre):.3f} >= QMAX = {QMAX}: too close to a cusp " "(t -> 1^+ or t -> 9); outside this script's series domain") pf_worst = mp.mpf(0) for te in (mp.mpf(-1), mp.mpf(-3), mp.mpf(1) / 2): pf_worst = max(pf_worst, max(_phys_pf_residual(te))) pf_tol = mp.mpf(10) ** (-(dps - 8)) if not pf_worst < pf_tol: # NaN-safe raise ArithmeticError( f"Picard-Fuchs residual {mp.nstr(pf_worst, 3)} >= {mp.nstr(pf_tol, 2)} on the served " "periods at t = -1, -3, 1/2: the transport operator is not certified; refusing") h0, h1 = PHYS_HEIGHTS qH0, iH0 = _phys_route_H(_phys_polygon(tval, mp.mpf(h0), sign)) qH1, iH1 = _phys_route_H(_phys_polygon(tval, mp.mpf(h1), sign)) qP, iP = _phys_route_P(_phys_polygon(tval, mp.mpf(h0), sign)) d_HP = _digits_c(qH0, qP) d_HH = _digits_c(qH0, qH1) if not (d_HP >= dps - 6 and d_HH >= dps - 6): # NaN-safe raise ArithmeticError( f"continued nome NOT certified at t={t}: route H vs route P {d_HP:.1f} d, " f"path independence (heights {h0} / {h1}) {d_HH:.1f} d (need >= {dps - 6}); refusing") delta = mp.mpf(10) ** (-(dps // 2)) _, _, tau_naive, q_naive = curve_nome(mp.mpc(tval, sign * delta)) d_naive = _digits_c(-q_naive, qH0) q_use = -q_naive if naive else qH0 qa = abs(q_use) if qa >= mp.mpf(QMAX): raise ValueError( f"|q_s(t + i0)| = {float(qa):.3f} >= QMAX = {QMAX}: too close to a cusp " "(t -> 1^+ or t -> 9); outside this script's series domain") c011, c101, dl = _phys_hpl_control(tval, sign) if not min(c011, c101) >= dl - 3: # NaN-safe raise ArithmeticError( f"HPL boundary values NOT certified at t={t}: the served G(0,1,1) / G(1,0,1) at " f"t + i 10^-{int(dl)} vs the boundary formulas {c011:.1f} / {c101:.1f} d " f"(expect about {int(dl)}); refusing") g1, g011, g101 = _phys_hpl_boundary(tval, sign) Iu, r3, r6, cl2 = roots() N = qseries_len(qa, mp.mp.dps) # STARTING guess N0 only Eb1, eb1_resid, N1 = _ebar_certified( "Eb1 = Ebar_{1;-1}(z3;1;-q)", tval, qa, lambda NN: Ebar_1(1, -1, r3, mp.mpc(1), q_use, NN), N) Eb2, eb2_resid, N2 = _ebar_certified( "Eb2 = Ebar_{0,2;-2,0;2}(z3,z3;1,-1;-q)", tval, qa, lambda NN: Ebar_2(0, 2, -2, 0, 1, r3, r3, mp.mpc(1), mp.mpc(-1), q_use, NN), N) c_eb2 = mp.mpf(-108) * (1 + mp.mpf("1e-30")) if MUTATE else mp.mpf(-108) val = (-2 * mp.pi ** 2 / 3 * g1 - 8 * g011 + 4 * g101 - 108 * cl2 * Eb1 + c_eb2 * Eb2) / (4 * mp.mpc(tval)) err_val = (108 * cl2 * eb1_resid + 108 * eb2_resid) / (4 * abs(mp.mpf(tval))) qtol = mp.mpf(10) ** (-(mp.mp.dps - GUARD + QGUARD)) info = dict(qabs=qa, N=N, eb1_resid=eb1_resid, eb2_resid=eb2_resid, N1=N1, N2=N2, err_val=err_val, qtol=qtol, physical=True, sign=sign, naive=naive, q_s=qH0, tau_P=iP["tau"], digits_H_vs_P=d_HP, digits_path=d_HH, pf_worst=pf_worst, hpl_digits=(c011, c101, dl), naive_vs_H=d_naive, steps_H=iH0["n_steps"], steps_P=iP["n_steps"], N_haupt=iH0["N_series"], qmax_path=max(iH0["qmax_on_path"], iH1["qmax_on_path"])) return val, info def _eval_dispatch(tval): """--point dispatch: t < 1 -> the served eval_kite (untouched); t > 1 -> the continued arm; t = 1 refused by name (the threshold cusp, the branch point of the continuation).""" if tval > 1: return eval_kite_physical(tval, sign=-1 if SIGN_PLANT else 1, naive=NAIVE) if tval == 1: raise ValueError("t = 1 is the threshold cusp (|q_s| -> 1; the branch point of the " "continuation); refused") return eval_kite(tval) def _phys_line(info): """The [continued] line: the continuation gates of one physical point (all passed, else the point raised); the naive nome's distance printed for the record.""" tag = "" if info.get("sign", 1) == -1: tag = " | CONTROL t - i0" if info.get("naive"): tag += " | CONTROL naive nome" return (f" [continued] q_s(t+i0) = {mp.nstr(info['q_s'], 12)}; route H vs route P " f"{info['digits_H_vs_P']:.1f} d, path independence {info['digits_path']:.1f} d " f"({info['steps_H']} / {info['steps_P']} steps, |q_s| on the path <= {float(info['qmax_path']):.3f}); " f"Picard-Fuchs residual " f"{mp.nstr(info['pf_worst'], 2)}; HPL boundary control {info['hpl_digits'][0]:.1f} / " f"{info['hpl_digits'][1]:.1f} d (expect about {int(info['hpl_digits'][2])}); naive " f"ellipk nome vs continued {info['naive_vs_H']:.1f} d{tag}") def _phys_ref_key(p): """The reference literal a --point string coincides with, or None (matched by value).""" try: tv = _parse_t(p) except (ValueError, ZeroDivisionError): return None for key in list(REF_MINKOWSKI) + list(REF_RHO): if tv == _parse_t(key): return key return None def _phys_gate(p, v, dps): """The Minkowski reference gate at a vendored point: (lines, failures). Digits per component vs the literal, capped at min(dps, certified), bar = cap - DIGIT_BAR_SLACK; the sign of the imaginary part must agree (the i0 side). Empty lines = no literal.""" key = _phys_ref_key(p) if key is None: return [], [] lines, fails = [], [] if key in REF_MINKOWSKI: R = REF_MINKOWSKI[key] cap = min(dps, R["certified"]) bar = cap - DIGIT_BAR_SLACK ref_re = mp.mpf(R["re"]) ref_im = mp.mpf(R["im"]) d_re = min(_digits_c(v.real, ref_re), cap) d_im = min(_digits_c(v.imag, ref_im), cap) sign_ok = (v.imag < 0) == (ref_im < 0) ok = bool(d_re >= bar and d_im >= bar and sign_ok) lines.append(f" [minkowski] vs {R['src']}: Re {d_re:.1f} d, Im {d_im:.1f} d, " f"Im sign {'OK' if sign_ok else 'WRONG'}; bar min(dps {dps}, certified " f"{R['certified']}) - {DIGIT_BAR_SLACK} = {bar} -> {'PASS' if ok else 'FAIL'} " f"(goal-40 twin agrees {R['twin']['agrees_digits_re']:.2f} / " f"{R['twin']['agrees_digits_im']:.2f} d)") if not sign_ok: fails.append(f"t={p}: Im SIGN differs from the +i0 reference (the wrong i0 side)") if not d_re >= bar: # NaN-safe fails.append(f"t={p}: Re vs reference {d_re:.1f} d BELOW BAR {bar}") if not d_im >= bar: fails.append(f"t={p}: Im vs reference {d_im:.1f} d BELOW BAR {bar}") else: R = REF_RHO[key] cap = min(dps, R["certified"]) bar = cap - DIGIT_BAR_SLACK ref_im = mp.mpf(R["minus_rho"]) d_im = min(_digits_c(v.imag, ref_im), cap) sign_ok = (v.imag < 0) == (ref_im < 0) ok = bool(d_im >= bar and sign_ok) lines.append(f" [spectral] Im vs -rho({key}) ({R['src']}): {d_im:.1f} d, Im sign " f"{'OK' if sign_ok else 'WRONG'}; bar min(dps {dps}, certified {R['certified']}) - " f"{DIGIT_BAR_SLACK} = {bar} -> {'PASS' if ok else 'FAIL'}; Re: no independent " "reference at this point") if not sign_ok: fails.append(f"t={p}: Im SIGN differs from -rho({key}) (the wrong i0 side)") if not d_im >= bar: fails.append(f"t={p}: Im vs -rho({key}) {d_im:.1f} d BELOW BAR {bar}") return lines, fails def _control_banner(): if SIGN_PLANT: return " | SIGN-PLANT CONTROL (t - i0)" if NAIVE: return " | NAIVE-NOME CONTROL" return "" # --------------------------------------------------------------------------- # 6. Positive controls (runtime, at t=-1; independent of the reference literals) # --------------------------------------------------------------------------- def positive_controls(): Iu, r3, r6, cl2 = roots() t = mp.mpf(-1) psi1, psi2, tau, q = curve_nome(t) mq = -q N = qseries_len(abs(q), mp.mp.dps) out = {} # Wronskian closed form W = -12 pi i /(t(t-1)(t-9)) => kernel identity: W = -12 * mp.pi * Iu / (t * (t - 1) * (t - 9)) lhs = (1 / (Iu * mp.pi)) * psi1 ** 2 / W / t rhs = 1 - 4 * Ebar_1(0, -1, r3, mp.mpc(-1), mq, N) out["kernel_id"] = abs(lhs - rhs) # HPL <-> Ebar, depth 1: g1id = 3 * (Ebar_1(1, 0, mp.mpc(-1), mp.mpc(1), mq, N) - Ebar_1(1, 0, r6, mp.mpc(1), mq, N)) out["G1_id"] = abs(g1id - mp.log(1 - t)) # HPL <-> Ebar, depth 2 (2o=2) -- the exact code path of the Eb_02 word: one = mp.mpc(1) g11id = 9 * (Ebar_2(0, 1, -1, 0, 1, -one, -one, one, one, mq, N) - Ebar_2(0, 1, -1, 0, 1, -one, r6, one, one, mq, N) - Ebar_2(0, 1, -1, 0, 1, r6, -one, one, one, mq, N) + Ebar_2(0, 1, -1, 0, 1, r6, r6, one, one, mq, N)) out["G11_id"] = abs(g11id - mp.log(1 - t) ** 2 / 2) return out # --------------------------------------------------------------------------- def _digits(a, b): a, b = mp.mpf(a), mp.mpf(b) if a == b: return float(mp.mp.dps) return float(-mp.log10(abs((a - b) / b))) def _digits_c(a, b): """Relative agreement digits of two (complex) numbers, -log10 |a-b|/|b|; dps when equal.""" if a == b: return float(mp.mp.dps) return float(-mp.log10(abs(a - b) / abs(b))) def _parse_t(s): if "/" in s: a, b = s.split("/") return mp.mpf(a) / mp.mpf(b) return mp.mpf(s) def _cert_ok(info): """True iff the value-level l1 bound is strictly below tol (NaN-safe).""" return bool(info["err_val"] < info["qtol"]) def _cert_line(info): """Certified-bound line. The relation is COMPARED (a hardcoded '<' once printed a FALSE '1.12e-110 < tol 1.0e-110' under an escalation sabotage at t=-7/3); a '>=' here is collected by main() as a gate failure.""" rel = "<" if _cert_ok(info) else ">=" return (f" [certified] Eb1 N/2N agree {mp.nstr(info['eb1_resid'], 3)} " f"(N={info['N1']}), Eb2 N/2N agree {mp.nstr(info['eb2_resid'], 3)} " f"(N={info['N2']}); value-level l1 bound {mp.nstr(info['err_val'], 3)} " f"{rel} tol {mp.nstr(info['qtol'], 3)}") def _paper_digits_check(key, v): """Compare the computed value with the digits printed in the paper. Returns (ok, line). Agreement means |computed - printed| < one unit in the last printed decimal place, i.e. the printed string is the truncation or the rounding of the computed value.""" printed, where = PAPER_PRINTED[key] places = len(printed.split(".")[1]) unit = mp.mpf(10) ** (-places) diff = abs(mp.mpf(v) - mp.mpf(printed)) ok = bool(diff < unit) line = (f" [paper] printed {printed}... ({where}); computed " f"{mp.nstr(v, places + 4)}: |diff| = {mp.nstr(diff, 2)} " f"{'<' if ok else '>='} 10^-{places} -> {'PASS' if ok else 'FAIL'}") return ok, line def run_all(dps, extra_points): """Full pass at working precision dps+GUARD; returns results dict.""" mp.mp.dps = dps + GUARD t0 = time.time() ctl = positive_controls() t_ctl = time.time() - t0 res = {"controls": ctl, "gate": {}, "points": {}, "t_ctl": t_ctl} for key in REF: t1 = time.time() v, info = eval_kite(_parse_t(key)) res["gate"][key] = (v, info, time.time() - t1) for p in extra_points: if p in REF: continue t1 = time.time() try: v, info = _eval_dispatch(_parse_t(p)) res["points"][p] = (v, info, time.time() - t1) except (ValueError, ArithmeticError) as e: res["points"][p] = (None, str(e), time.time() - t1) res["t_total"] = time.time() - t0 return res def _fail(code, label, msgs): for msg in msgs: print(f"{label} FAIL: {msg}") print(f"RESULT: FAIL ({label}; exit {code})") return code def _point_only(args, dps): print("kite-equal-evaluate: equal-mass kite J(1,1,1,1,1) eps^0 -- POINT-ONLY mode" + _control_banner()) print(f"dps={dps} (+{GUARD} guard); control battery + reference table " "SKIPPED; nome branch gate, certified N/2N q-series agreement " "and value-level l1 bound retained per point\n") wall0 = time.time() mp.mp.dps = dps + GUARD failures = [] gated = [] for p in args.point: t1 = time.time() try: v, info = _eval_dispatch(_parse_t(p)) except ValueError as e: # named domain refusal (fail-closed) print(f"REFUSED t={p}: {e}") return EXIT_USAGE except (ArithmeticError, RuntimeError) as e: print(f"REFUSED t={p} (cannot certify): {e}") return EXIT_GATE if info.get("physical"): head = f"t={p}: I_kite(t + i0)" print(f"{head} Re = {mp.nstr(v.real, dps)}") print(f"{' ' * len(head)} Im = {mp.nstr(v.imag, dps)}") print(f" [|q_s|={float(info['qabs']):.3f} N={info['N']} " f"{time.time() - t1:.1f}s]") print(_cert_line(info)) print(_phys_line(info)) if not _cert_ok(info): failures.append( f"t={p}: certified l1 bound {mp.nstr(info['err_val'], 3)} " f">= tol {mp.nstr(info['qtol'], 3)}") lines, fails = _phys_gate(p, v, dps) for line in lines: print(line) failures.extend(fails) if lines: gated.append(p) continue print(f"t={p}: I_kite = {mp.nstr(v, dps)}") print(f" [|q|={float(info['qabs']):.3f} N={info['N']} " f"{time.time() - t1:.1f}s]") print(_cert_line(info)) if not _cert_ok(info): failures.append( f"t={p}: certified l1 bound {mp.nstr(info['err_val'], 3)} " f">= tol {mp.nstr(info['qtol'], 3)}") print(f"\nTOTAL wall time: {time.time() - wall0:.1f} s " "(measured; point-only)") if failures: return _fail(EXIT_GATE, "CERTIFIED-BOUND" if not gated else "MINKOWSKI-GATE", failures) if gated: print("RESULT: PASS (point-only: certified bounds met; Minkowski reference " f"gate at t={', '.join(gated)})") return EXIT_OK print("RESULT: PASS (point-only: certified bounds met; no reference " "comparison in this mode)") return EXIT_OK def main(): global MUTATE, SIGN_PLANT, NAIVE ap = argparse.ArgumentParser( description="Equal-mass kite top master J(1,1,1,1,1) at eps^0 -- the " "five-word modular Gamma_1(6) closed form, evaluated at " "runtime and compared with independent reference values and " "the digits printed in the paper. Domain: Euclidean t<1, " "t!=0 (|q|<0.90), and the physical region t>1, t!=9 by the " "elliptic Feynman-i0 continuation (--point t>1; complex value, " "gated at t = 4, 16 and t = 12, 50, 100). Exit codes 0/2/3/4/5, " "see the header.") ap.add_argument("--dps", type=int, default=None, help=f"decimal digits (default {DPS_DEFAULT}; --full sets " f"{DPS_FULL}, --paper {DPS_PAPER}); the run works at dps+{GUARD}") ap.add_argument("--full", action="store_true", help=f"the same suite at dps {DPS_FULL} (the held-out digit " "counts quoted in the paper are tested by --paper)") ap.add_argument("--paper", action="store_true", help=f"the same suite at dps {DPS_PAPER}: tests the held-out " f"digit counts quoted in the paper ({PAPER_HELDOUT_DIGITS['-3']} at " f"t=-3, {PAPER_HELDOUT_DIGITS['-10']} at t=-10, each capped by its " "stored reference)") ap.add_argument("--point", action="append", default=[], help="extra t (rational string, e.g. --point=-7/3); Euclidean " "t<1 on the served path, physical t>1 (t!=1, 9) on the " "continued path with the value as Re / Im, gated where a " "reference literal exists (t = 4, 16; Im at 12, 50, 100); " "repeatable") ap.add_argument("--check", action="store_true", help="two-precision rule: rerun everything at dps+60, diff") ap.add_argument("--point-only", action="store_true", help="plain point evaluation: skip the control battery " "(positive controls + reference table + paper digits); " "per-point certified q-series bounds + branch gate " "retained; requires --point; incompatible with --check") ap.add_argument("--mutate", action="store_true", help="control run: perturb the -108 coefficient of the " "depth-2 word by a relative 1e-30; the reference gate " "must then fail (exit 3) -- proves the gate is live") ap.add_argument("--sign-plant", action="store_true", help="control run (physical arm; needs a gated --point t>1): " "the i0 sign flipped -- the lower half plane and the " "conjugate boundary values; the reference gate must fail " "by name on the sign of the imaginary part (exit 3)") ap.add_argument("--naive", action="store_true", help="control run (physical arm; needs a gated --point t>1): " "the principal-branch ellipk nome at t + i delta fed to the " "same words (the naive per-piece t+i0); the reference gate " "must fail by name (exit 3)") args = ap.parse_args() if args.full and args.paper: ap.error("--full and --paper are two tiers of the same suite (dps " f"{DPS_FULL} / {DPS_PAPER}); choose one, or give --dps") dps = args.dps if args.dps is not None else ( DPS_PAPER if args.paper else DPS_FULL if args.full else DPS_DEFAULT) if dps < 10: ap.error("--dps must be at least 10") if args.point_only: # battery remains default-on for --check and bare runs (fail-closed) if args.check: ap.error("--point-only and --check are mutually exclusive: " "--check keeps the full battery by design") if not args.point: ap.error("--point-only requires at least one --point") if args.mutate: ap.error("--mutate needs the reference table (it is what catches " "the perturbation); drop --point-only") MUTATE = bool(args.mutate) if args.sign_plant and args.naive: ap.error("--sign-plant and --naive are two controls of the physical arm; " "choose one") if args.sign_plant or args.naive: if not any(_phys_ref_key(p) for p in args.point): ap.error("--sign-plant / --naive need a --point t>1 that carries a " "reference literal (t = 4, 16 or t = 12, 50, 100): the control " "is the reference gate failing by name") if args.mutate: ap.error("--mutate and the physical controls are separate controls; " "run them one at a time") SIGN_PLANT = bool(args.sign_plant) NAIVE = bool(args.naive) if args.point_only: return _point_only(args, dps) print("kite-equal-evaluate: equal-mass kite J(1,1,1,1,1) eps^0 -- five-word " "Gamma_1(6) closed form" + (" | MUTATE CONTROL" if MUTATE else "") + _control_banner()) print(f"dps={dps} (+{GUARD} guard); q-series auto-sized; nome branch-gated " "(ellipk vs path-pinned hauptmodul)") if MUTATE: print("[mutate] the -108 coefficient of the depth-2 word " "Ebar_{0,2;-2,0;2} is perturbed by a relative 1e-30 -- a wrong " "closed form; the reference gate must now fail (exit 3)") print() wall0 = time.time() try: r = run_all(dps, args.point) except ValueError as e: print(f"REFUSED: {e}") return EXIT_USAGE except (ArithmeticError, RuntimeError) as e: print(f"REFUSED (cannot certify): {e}") return EXIT_GATE mp.mp.dps = dps + GUARD failures = [] # gate failures: collected, printed, then exit 3 ctl_tol = mp.mpf(10) ** (-(dps + CONTROL_GUARD)) print("positive controls at t=-1 (1607.01571 identities, computed now; " f"each residual must be < {mp.nstr(ctl_tol, 2)}):") print(f" kernel id (psi1^2/(Wt) vs 1-4Ebar_0;-1) resid = {mp.nstr(r['controls']['kernel_id'], 3)}") print(f" G(1;y) = 3[Ebar_1;0 combo] (depth 1) resid = {mp.nstr(r['controls']['G1_id'], 3)}") print(f" G(1,1;y) = 9[Ebar_0,1;-1,0;2 combo] (d2) resid = {mp.nstr(r['controls']['G11_id'], 3)}") print(f" [{r['t_ctl']:.1f} s]\n") for name, resid in r["controls"].items(): if not resid < ctl_tol: # NaN-safe failures.append(f"positive control {name}: residual " f"{mp.nstr(resid, 3)} >= {mp.nstr(ctl_tol, 2)}") print(f"{'t':>6} | {'I_kite(t) (computed now)':>44} | reference (independent, digits MEASURED this run)") worst_held = None held = {} paper_fail = [] for key, (v, info, dt) in r["gate"].items(): ref = _refval(key) cap = min(dps, _refdigits(key)) d = min(_digits(v, ref), cap) tag = "fit-grid" if REF[key]["fitgrid"] else "HELD-OUT" if not REF[key]["fitgrid"]: held[key] = d worst_held = d if worst_held is None else min(worst_held, d) print(f"{key:>6} | {mp.nstr(v, 40):>44} | {d:6.1f} d {tag} ({REF[key]['src']}) " f"[|q|={float(info['qabs']):.3f} N={info['N']} {dt:.1f}s]") print(_cert_line(info)) if not _cert_ok(info): failures.append( f"t={key}: certified l1 bound {mp.nstr(info['err_val'], 3)} >= " f"tol {mp.nstr(info['qtol'], 3)}") bar = cap - DIGIT_BAR_SLACK if not d >= bar: # NaN-safe failures.append( f"t={key}: reference digits {d:.1f} BELOW BAR {bar:.1f} " f"(= min(dps={dps}, ref {_refdigits(key)}d) - {DIGIT_BAR_SLACK})") if key in PAPER_PRINTED: ok, line = _paper_digits_check(key, v) print(line) if not ok: paper_fail.append(f"t={key}: computed {mp.nstr(v, 12)} does not " f"reproduce the printed {PAPER_PRINTED[key][0]}...") gated = [] for key, (v, info, dt) in r["points"].items(): if v is None: print(f"{key:>6} | {'DOMAIN/BRANCH ERROR':>44} | {info}") failures.append(f"t={key}: {info}") continue if info.get("physical"): print(f"{key:>6} | {'Re ' + mp.nstr(v.real, 40):>44} | physical point t + i0 " f"[|q_s|={float(info['qabs']):.3f} N={info['N']} {dt:.1f}s]") print(f"{'':>6} | {'Im ' + mp.nstr(v.imag, 40):>44} |") print(_cert_line(info)) print(_phys_line(info)) if not _cert_ok(info): failures.append( f"t={key}: certified l1 bound {mp.nstr(info['err_val'], 3)} >= " f"tol {mp.nstr(info['qtol'], 3)}") lines, fails = _phys_gate(key, v, dps) for line in lines: print(line) failures.extend(fails) if lines: gated.append(key) continue print(f"{key:>6} | {mp.nstr(v, 40):>44} | (no stored reference; new point) " f"[|q|={float(info['qabs']):.3f} N={info['N']} {dt:.1f}s]") print(_cert_line(info)) if not _cert_ok(info): failures.append( f"t={key}: certified l1 bound {mp.nstr(info['err_val'], 3)} >= " f"tol {mp.nstr(info['qtol'], 3)}") print(f"\nmin HELD-OUT reference digits this run: {worst_held:.1f} " f"(caps: dps and stored reference lengths)") # the paper's quoted held-out digit counts (160 at t=-3, 130 at t=-10), # each capped by the stored reference: testable only when the requested # precision less the gate slack reaches the quoted count, i.e. when the # count is limited by the reference and not by --dps (--paper, or --dps >= # count + 2); the measured count is compared as printed (one decimal), the # form the quoted count was read in for key, need in PAPER_HELDOUT_DIGITS.items(): bar = dps - DIGIT_BAR_SLACK if bar < need: print(f"[paper] quoted held-out agreement at t={key}: {need} d -- not " f"testable at dps {dps} (bar {bar} < {need}; use --paper or " f"--dps >= {need + DIGIT_BAR_SLACK})") else: meas = round(held[key], 1) ok = meas >= need print(f"[paper] quoted held-out agreement at t={key}: {need} d (capped by " f"the stored {_refdigits(key)}-digit reference); this run measures " f"{meas:.1f} d -> {'PASS' if ok else 'FAIL'}") if not ok: paper_fail.append(f"t={key}: measured {meas:.1f} d < the " f"paper's quoted {need} d") print(f"TOTAL wall time: {time.time() - wall0:.1f} s (measured)") if failures: return _fail(EXIT_GATE, "REFERENCE-GATE", failures) if paper_fail: return _fail(EXIT_PAPER, "PAPER-DIGITS", paper_fail) if args.check: print(f"\n--check: rerunning everything at dps {dps + 60} ...") v1 = {k: r["gate"][k][0] for k in r["gate"]} v1.update({k: r["points"][k][0] for k in r["points"]}) t1 = time.time() try: r2 = run_all(dps + 60, args.point) except (ValueError, ArithmeticError, RuntimeError) as e: print(f"--check rerun REFUSED: {e}") return _fail(EXIT_CHECK, "TWO-PRECISION", [str(e)]) mp.mp.dps = dps + 60 + GUARD v2 = {k: r2["gate"][k][0] for k in r2["gate"]} v2.update({k: r2["points"][k][0] for k in r2["points"]}) ok = True for k in v1: if v1[k] is None or v2[k] is None: print(f" t={k}: skipped (domain error)") continue dk = _digits_c(v1[k], v2[k]) if isinstance(v1[k], mp.mpc) else _digits(v1[k], v2[k]) tgt = dps - 5 ok = ok and dk >= tgt print(f" t={k}: stable to {dk:.1f} d (target >= {tgt})") print(f"--check wall time: {time.time() - t1:.1f} s (measured)") print(f"--check verdict: {'PASS' if ok else 'FAIL'}") if not ok: return _fail(EXIT_CHECK, "TWO-PRECISION", ["dps vs dps+60 agreement below dps-5 at some point"]) print("RESULT: PASS (reference gate" + (f" with the Minkowski reference at t={', '.join(gated)}" if gated else "") + ", certified bounds, positive controls, " "paper digits" + (", two-precision check" if args.check else "") + ")") return EXIT_OK if __name__ == "__main__": sys.exit(main())