#!/usr/bin/env python3
# SPDX-License-Identifier: MIT
# Copyright (c) 2026 Anthropic, PBC. Part of BootLoops 1.0; see LICENSE and NOTICE beside this file.

"""kite-boundary-general.py -- the (1,1,x)-LINE FORM of kite-boundary.py (the (1,1,2) boundary script beside it)
generalized to the (1,1,x) line: the same engine (graded 42-state Frobenius bounded branch at s = 0 + adaptive
exact-Taylor stepping, flint-arb backend with midpoint trim) with these differences ONLY:
  (1) CONN is loaded from --conn (kite-connection-A14-x3.json, or kite-connection-A14.json for the x = 2 positive control)
      instead of the embedded (1,1,2) literal; the DE singular points (the transporter's pole list) are READ from the connection's denominators
      (kite-boundary.py hardcodes {0, 1, 2, 6 -/+ 4 sqrt2}); (2) closed_form_seed is the x-general table of kite-seed-general.py
      (M10 = (2+x) T1^2 + x^2 V(1,1,x) -- derived there; --neg-m10 = the alternative 2x coefficient as a NEGATIVE CONTROL through the kernel gate);
  (3) the Frobenius exit s_exit is a list (--s-exits; kite-boundary.py's -1/10 plus a second exit inside the nearest-DE-pole radius: for x = 3 the DE pole
      7 - 4 sqrt3 = 0.0718 lies INSIDE kite-boundary.py's |s_exit| = 0.1, so the bounded-branch radius must be MEASURED, not assumed: coefficient growth
      rate + exit-independence of the transported values are reported); (4) transport stops at every --points value (the two AMFlow points and
      the N_kite point -2), walls per leg; (5) the independent-reference gate is --ref (kite-boundary-derived.json for the x = 2 control) instead of the sibling file;
  (6) a JSON receipt (--receipt) with every number, the sha256 of the inputs and UTC time stamps. mp.dps is set inside main() as in kite-boundary.py.
  (7) the connection strings are parsed by the sibling kite_exact_parse.py (sha256-pinned; no sympy at runtime); the eps-graded
      coefficient lists come from <conn stem>-graded.json beside the connection when it exists (the shipped ones are pinned and checked entry by
      entry against the strings at every start -- the lists a sympy parse produced, kept in that exact form) and otherwise from the parser's own
      reduced form; the DE singular points come from the parser's exact factorization of the s-only denominator letters (python-flint
      fmpz_poly.factor; the census strings, floats and mpf values are in the form a sympy census prints).

WHAT IT COMPUTES.  The 42-component boundary vector (eps^-2..eps^0 of the 14 masters) of the unequal-mass kite
on the (1,1,x) line -- propagator masses (1, 0, 1, 0, sqrt x): D1 = k1^2 - 1, D2 = k2^2, D3 = (k1-k2)^2 - 1,
D4 = (k1-p)^2, D5 = (k2-p)^2 - x, p^2 = s (Euclidean s < 0), measure int d^dk/(i pi^{d/2}) per loop, d = 4 - 2 eps,
mu = 1 -- at every --points value, from classical constants alone: the unique solution of the exact rational IBP
connection dM/ds = A(4 - 2 eps, s) M that stays bounded at the regular singular point s = 0, with the s = 0 value
given in closed form by the two-loop vacuum ring (tadpoles T(M) and sunsets V(M1, M2, M3); the table of the 14
masters is in closed_form_seed below and derived in kite-seed-general.py), summed by the Frobenius series at every
--s-exits point and carried on by adaptive Taylor steps of the same exact connection (SAFETY 0.25, order
1.7*dps + 20; flint-arb ball radii trimmed to midpoints, correctness coming from the gates, not the balls).  At
x = 2 with --conn kite-connection-A14.json and --ref kite-boundary-derived.json it reproduces kite-boundary.py's
s = -2 vector (the positive control: the gate line prints the worst agreement of the 42 components and of
N_kite = -8 M13[eps^0](-2)).  At x = 3 it is the engine that derives the shipped seed kite-boundary-derived-x3.json
(exit -1/40, transported to s = -5/3), the seed kite-evaluate.py --masses 1,1,3 transports from.  The printed
checks: the xi identity of the line (xi(1,1,2) + 8 Catalan = 0, or xi(1,1,3) + (8/sqrt3) Cl2(pi/3) = 0), V(0,0,M)
against its exact Gamma form, the eps^-1 cancellation of the top seed, |B0 y0| (the seed in the kernel of the
residue matrix), the measured coefficient growth of the Frobenius series at every exit and the agreement of the
same point reached from two exits.  Belongs to the unequal-mass kite of the BootLoops results paper (M. D.
Schwartz, "Polylogarithmic, elliptic, and Calabi-Yau Feynman integrals from a hybrid bootstrap", the table of
thirty integrals, row 14) and the site page bootloops.ai/diagrams/kite.

INPUTS: --conn (the connection JSON); <conn stem>-graded.json beside it when present (the shipped ones are
sha256-pinned and checked against the strings); the pinned sibling kite_exact_parse.py; the optional --ref gate
file (kite-boundary-derived.json schema, at s = -2).  OUTPUT: the JSON receipt named by --receipt (the transported
vector at every exit and point, the checks, the DE singular points, the walls, the input hashes and a UTC stamp)
and the printed summary.  mp.dps is set inside main() after argparse (module-level mpf at the import precision
would silently truncate every constant).

Requires: python >= 3.9, mpmath >= 1.3, python-flint >= 0.6 (the singular-point census needs it; arb is also
used for the transport unless --no-flint).  No AMFlow, Kira, Mathematica or network.

CLI (write --s-exits=-1/10 and --points=-2 with '=': a bare leading-dash value confuses argparse):
  python3 kite-boundary-general.py --conn kite-connection-A14.json --x 2 --dps 30 --s-exits=-1/10 --points=-2 --ref kite-boundary-derived.json --receipt receipt_x2.json
      # the x = 2 positive control (seconds)
  python3 kite-boundary-general.py --conn kite-connection-A14-x3.json --x 3 --dps 190 --s-exits=-1/40 --points=-5/3 --receipt receipt_x3.json
      # the x = 3 seed as shipped (190 digits)

Exit codes: 0 (the gate and exit-independence figures are printed, not raised); 2 mpmath not installed or the
--conn file not found; 3 a pinned sibling altered or missing, or a connection entry outside the parser's grammar;
4 python-flint missing.
"""

import argparse
import json
import os
import time

import hashlib
import importlib.util
import sys
from fractions import Fraction

try:
    import mpmath as mp
except ImportError:
    print("ERROR: this script needs mpmath (pip install mpmath)")
    sys.exit(2)

try:
    from flint import arb, arb_mat, ctx as _fctx, fmpz
    HAVE_FLINT = True
except Exception:
    HAVE_FLINT = False

# ---------------------------------------------------------------------------
# The exact rational 14x14 IBP connection A(d, s) of the selected line: loaded
# in main() from --conn (kite-connection-A14-x3.json, or kite-connection-A14.json
# for the x = 2 positive control).
# ---------------------------------------------------------------------------
CONN = None   # loaded in main() from --conn

# ---------------------------------------------------------------------------
# Independent-reference gate (GATE ONLY -- never used in the construction): the
# kite-boundary-derived.json-schema file named by --ref, compared at s = -2 (the
# x = 2 control reads kite-boundary-derived.json).  Without --ref the derivation
# still runs and no gate line is printed.
# ---------------------------------------------------------------------------
HERE = os.path.dirname(os.path.abspath(__file__))
GATE_JSON = None   # set in main() from --ref
SCRIPT = "kite-boundary-general"
# sha256 of the pinned siblings: the parser module and the
# two shipped graded files; a graded file of another name (a connection of another line) is used unpinned
PINS = {
    "kite_exact_parse.py": "1724d2b18d66801ba4dca2a7335fa6c244e3978dc84b8fa63eca96f87393a877",
    "kite-connection-A14-graded.json": "ca2bb12bde134dc0b90498fad3220279ad705fc0bfa7c9c981837dc7c28a9493",
    "kite-connection-A14-x3-graded.json": "19405d4ac20b2449ec02092b085ca1ce7a8f5e94bd948b15525ed5ba8e1cd41c",
}


def _pinned_path(name_or_path):
    """Path of a sibling file; a PINNED file (the exact parser module, a shipped graded file) is refused on
    any byte change (exit 3, the recorded and recomputed sha256 named)."""
    path = name_or_path if os.path.isabs(name_or_path) else os.path.join(HERE, name_or_path)
    want = PINS.get(os.path.basename(path))
    if want is not None:
        if not os.path.exists(path):
            sys.stderr.write(f"{SCRIPT} REFUSED: pinned file {os.path.basename(path)} 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: {os.path.basename(path)} integrity pin mismatch (recorded {want}, "
                             f"recomputed {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_exact_parse.py")   # the exact parser of the connection strings (no sympy)
GRADED_PATH = None   # set in main() to <conn stem>-graded.json beside --conn when it exists


def load_gate():
    """Return ({(master, eps_order): mpf_string}, dps) or (None, 0) if absent.
    Only eps^-2..eps^0 entries are gate targets (higher orders in the stored
    file are outside the eps^-2..eps^0 grading derived here)."""
    if GATE_JSON is None or not os.path.exists(GATE_JSON):
        return None, 0
    bnd = json.load(open(GATE_JSON))
    assert int(bnd["s"]) == -2
    ref = {}
    for i_str, co in bnd["boundary"].items():
        for o_str, v in co.items():
            if KMIN <= int(o_str) <= KMAX:
                ref[(int(i_str), int(o_str))] = v
    return ref, int(bnd["dps"])

KMIN, KMAX = -2, 0
NEPS = KMAX - KMIN


# ---------------------------------------------------------------------------
# 1. Closed-form vacuum seed (the derived boundary condition at s = 0)
# ---------------------------------------------------------------------------
def T_tadpole(M, order=2):
    """[eps^-1 .. eps^order] of T(M) = Gamma(1+eps) M^(1-eps)/(eps(1-eps)).
    Exact exponential form: T = (M/eps) exp((1-g-lnM) eps + sum_{k>=2}(1+(-1)^k zeta(k)) eps^k/k)."""
    n = order + 1
    with mp.extradps(15):
        Mv = mp.mpf(M)
        cs = [mp.mpf(0), 1 - mp.euler - mp.log(Mv)]
        for k in range(2, n + 1):
            cs.append((1 + (-1) ** k * mp.zeta(k)) / k)
        p = [mp.mpf(1)]
        for j in range(1, n + 1):
            p.append(sum(k * cs[k] * p[j - k] for k in range(1, j + 1)) / j)
        return [Mv * pj for pj in p]


def _cl2(theta):
    return mp.polylog(2, mp.exp(mp.mpc(0, 1) * theta)).imag


def _xi(x, y, z):
    """Martin (2.20)/(2.21) at Q^2=1 for x,y <= z, all > 0; lambda<0 via the
    Clausen (FJJ 4.22) form. xi(1,1,2) = -8 Catalan."""
    lam = x * x + y * y + z * z - 2 * (x * y + y * z + z * x)
    if lam == 0:
        return mp.mpf(0)
    if lam > 0:
        R = mp.sqrt(lam)
        u = (z + x - y - R) / (2 * z)
        v = (z + y - x - R) / (2 * z)
        return R * (2 * mp.log(u) * mp.log(v) - mp.log(x / z) * mp.log(y / z)
                    - 2 * mp.polylog(2, u) - 2 * mp.polylog(2, v) + mp.pi ** 2 / 3)
    r = mp.sqrt(-lam)
    ax = mp.acos((y + z - x) / (2 * mp.sqrt(y * z)))
    ay = mp.acos((z + x - y) / (2 * mp.sqrt(z * x)))
    az = mp.acos((x + y - z) / (2 * mp.sqrt(x * y)))
    return -2 * r * (_cl2(2 * ax) + _cl2(2 * ay) + _cl2(2 * az))


def _I_finite(x, y, z):
    """Martin hep-ph/0111209 eqs. (2.19), (2.26)-(2.28) at Q^2=1; x<=y<=z."""
    nz = sum(1 for m in (x, y, z) if m == 0)
    if nz == 3:
        return mp.mpf(0)
    if nz == 2:
        L = mp.log(z)
        return z * (-L * L / 2 + 2 * L - mp.mpf(5) / 2 - mp.pi ** 2 / 6)
    if nz == 1:
        a, b = z, y
        La, Lb = mp.log(a), mp.log(b)
        if a == b:
            return -a * La * La + 4 * a * La - 5 * a
        return ((a - b) * (mp.polylog(2, b / a) - mp.log(a / b) * mp.log(a - b)
                           + La * La / 2 - mp.pi ** 2 / 6)
                - mp.mpf(5) / 2 * (a + b) + 2 * a * La + 2 * b * Lb - a * La * Lb)
    Lx, Ly, Lz = mp.log(x), mp.log(y), mp.log(z)
    return ((x - y - z) / 2 * Ly * Lz + (y - x - z) / 2 * Lx * Lz
            + (z - x - y) / 2 * Lx * Ly + 2 * (x * Lx + y * Ly + z * Lz)
            - mp.mpf(5) / 2 * (x + y + z) - _xi(x, y, z) / 2)


def V_vacuum(M1, M2, M3):
    """{-2,-1,0}: Laurent of the 2-loop vacuum sunset V(M1,M2,M3) in this script's
    normalization (int d^dk/(i pi^{d/2}) per loop, mu = 1, no e^{gamma_E eps}).
    V = -e^{-2 gamma eps} I_bold|_{Q^2=1}; I_bold from TSIL eq. (2.34)."""
    with mp.extradps(20):
        x, y, z = sorted(mp.mpf(m) for m in (M1, M2, M3))
        c = x + y + z
        A = [m * (mp.log(m) - 1) if m > 0 else mp.mpf(0) for m in (x, y, z)]
        Aeps = [m * (-1 - mp.pi ** 2 / 12 + mp.log(m) - mp.log(m) ** 2 / 2) if m > 0
                else mp.mpf(0) for m in (x, y, z)]
        # I_bold Laurent: [-2] = -c/2, [-1] = sum A - c/2, [0] = I + sum Aeps
        Ib2 = -c / 2
        Ib1 = sum(A) - c / 2
        Ib0 = _I_finite(x, y, z) + sum(Aeps)
        g = mp.euler
        # V = -(1 - 2g eps + 2g^2 eps^2 + O(eps^3)) * (Ib2/eps^2 + Ib1/eps + Ib0 + ...)
        out = {-2: -Ib2,
               -1: -(Ib1 - 2 * g * Ib2),
               0: -(Ib0 - 2 * g * Ib1 + 2 * g * g * Ib2)}
        return {k: +v for k, v in out.items()}


X_STR = "3"          # the line's x, set from --x in main()
NEG_M10 = False      # --neg-m10: the alternative coefficient 2x for M10 (a negative control)


def closed_form_seed(mutate=None):
    """The 42-component boundary vector y(0) on the (1,1,x) line (derivation: kite-seed-general.py docstring):
    M0=M2=T(1)T(x); M1=M3=T(1)T(x)/x; M4=M7=T(1)^2; M5=-M6=V(0,1,0); M8=V(1,1,x); M9=T(1)^2+xV(1,1,x); M10=(2+x)T(1)^2+x^2V(1,1,x);
    M11=V(0,1,x); M12=[V(0,1,x)-V(0,1,0)]/x; M13=[V(1,1,x)-V(1,1,0)-V(0,1,x)+V(0,1,0)]/x (eps^-2, eps^-1 = 0)."""
    X = mp.mpf(X_STR)
    T1 = T_tadpole(1, 2)
    Tx = T_tadpole(X, 2)

    def prod2(a, b):
        out = {k: mp.mpf(0) for k in (-2, -1, 0)}
        for i, xx in enumerate(a[:3]):
            for j, yy in enumerate(b[:3]):
                k = (i - 1) + (j - 1)
                if k <= 0:
                    out[k] += xx * yy
        return out

    V010 = V_vacuum(0, 1, 0)
    V110 = V_vacuum(1, 1, 0)
    V01x = V_vacuum(0, 1, X)
    V11x = V_vacuum(1, 1, X)
    T1Tx = prod2(T1, Tx)
    T1T1 = prod2(T1, T1)
    h = 1 / X
    c10 = (2 * X) if NEG_M10 else (2 + X)
    seed = {0: T1Tx, 1: {k: v * h for k, v in T1Tx.items()}, 2: dict(T1Tx),
            3: {k: v * h for k, v in T1Tx.items()}, 4: T1T1, 5: V010,
            6: {k: -v for k, v in V010.items()}, 7: dict(T1T1), 8: V11x,
            9: {k: T1T1[k] + X * V11x[k] for k in (-2, -1, 0)},
            10: {k: c10 * T1T1[k] + X * X * V11x[k] for k in (-2, -1, 0)},
            11: V01x, 12: {k: (V01x[k] - V010[k]) * h for k in (-2, -1, 0)},
            13: {-2: mp.mpf(0), -1: mp.mpf(0),
                 0: (V11x[0] - V01x[0] - V110[0] + V010[0]) * h}}
    if mutate is not None:
        i, K, delta = mutate
        seed[i][K] += mp.mpf(delta)
    return seed


# ---------------------------------------------------------------------------
# 2. Graded system, exact Taylor data (exact rationals: the parsed connection)
# ---------------------------------------------------------------------------
def build_graded():
    """The eps-graded connection: Aco[(i, j, k)] = (pc, qc), the eps^k coefficient of
    A_ij(4 - 2 eps, s) as the exact rational function P(s)/Q(s) (ascending Fraction lists):
    READ from the graded file beside --conn when it exists (the shipped ones pinned; the lists
    a sympy parse produced, kept in that exact form so that the numerics consume exactly
    these integers) and CHECKED against the connection strings by the exact parser at
    every call (a mismatch is REFUSED by name, exit 3); else the parser's own reduced form."""
    N = len(CONN["masters"])
    if GRADED_PATH is None:
        Aco = kx.grade_canonical(CONN["A"], NEPS)
    else:
        try:
            Aco = kx.load_graded(CONN["A"], GRADED_PATH, NEPS)
        except kx.ParseError as ex:
            sys.stderr.write(f"{SCRIPT} REFUSED: {os.path.basename(GRADED_PATH)} vs the connection: {ex}\n")
            raise SystemExit(3)
    STATE = [(i, K) for i in range(N) for K in range(KMIN, KMAX + 1)]
    NS = len(STATE)
    ENT = []
    for n, (i, K) in enumerate(STATE):
        for m, (j, Kp) in enumerate(STATE):
            k = K - Kp
            if 0 <= k <= NEPS and (i, j, k) in Aco:
                pc, qc = Aco[(i, j, k)]
                ENT.append((n, m, pc, qc))
    IDX = {sK: n for n, sK in enumerate(STATE)}
    return STATE, IDX, NS, ENT


def b_taylor_exact(pc, qc, nmax):
    """Exact Taylor coefficients (Fractions), orders 0..nmax, of s*P/Q at s=0."""
    v = 0
    while qc[v] == 0:
        v += 1
    qt = qc[v:]
    c = [Fraction(0)] * (nmax + 2)
    for n in range(nmax + 2):
        t = pc[n] if n < len(pc) else Fraction(0)
        for k in range(1, min(n, len(qt) - 1) + 1):
            t -= qt[k] * c[n - k]
        c[n] = t / qt[0]
    shift = 1 - v
    out = [Fraction(0)] * (nmax + 1)
    for n in range(nmax + 1 - shift):
        out[n + shift] = c[n]
    return out


# ---------------------------------------------------------------------------
# 3. Numeric backends (flint arb preferred; pure mpmath fallback)
# ---------------------------------------------------------------------------
class Backend:
    def __init__(self, use_flint):
        self.flint = use_flint and HAVE_FLINT

    def setup(self):
        if self.flint:
            _fctx.prec = int(mp.mp.dps * 3.3322) + 80

    def from_rat(self, r):
        if self.flint:
            return arb(fmpz(r.numerator)) / arb(fmpz(r.denominator))
        return mp.mpf(r.numerator) / mp.mpf(r.denominator)

    def from_mpf(self, x):
        x = mp.mpf(x)
        if not self.flint:
            return x
        sgn, man, e, _ = x._mpf_
        if man == 0:
            return arb(0)
        a = arb(fmpz(int(man))) * (arb(2) ** int(e))
        return -a if sgn else a

    def to_mpf(self, a):
        if not self.flint:
            return a
        return mp.mpf(a.str(mp.mp.dps + 15, radius=False))

    def zero(self):
        return arb(0) if self.flint else mp.mpf(0)

    def trim(self, x):
        """Discard the arb error ball (keep midpoint). Avoids the ball-radius
        blow-up of deep recurrences (dependency problem); correctness is
        guaranteed by the two-precision rule + independent-reference gate, not by balls."""
        return x.mid() if self.flint else x

    def trim_vec(self, v):
        if not self.flint:
            return v
        return [x.mid() for x in v]

    def solve_shifted(self, B0rows, n, rhs, NS):
        """solve (n I - B0) x = rhs. B0rows: dense list of lists (backend numbers)."""
        if self.flint:
            M = arb_mat(NS, NS)
            for i in range(NS):
                for j in range(NS):
                    M[i, j] = -B0rows[i][j]
                M[i, i] += arb(n)
            R = arb_mat([[v] for v in rhs])
            S = M.solve(R)
            return [S[i, 0] for i in range(NS)]
        M = mp.matrix(NS, NS)
        for i in range(NS):
            for j in range(NS):
                M[i, j] = -B0rows[i][j]
            M[i, i] += n
        sol = mp.lu_solve(M, mp.matrix(rhs))
        return [sol[i] for i in range(NS)]


def frobenius_eval(be, STATE, NS, ENT, seed_vec, s_exit, nmax):
    """Sum the Frobenius series of the bounded solution at s_exit (|s_exit| < 1)."""
    serx = [(r, c, b_taylor_exact(pc, qc, nmax)) for (r, c, pc, qc) in ENT]
    B0rows = [[be.zero() for _ in range(NS)] for _ in range(NS)]
    for (r, c, t) in serx:
        if t[0] != 0:
            B0rows[r][c] += be.from_rat(t[0])
    ser = [(r, c, [be.from_rat(x) for x in t[1:]]) for (r, c, t) in serx]
    Y = [list(seed_vec)]
    for n in range(1, nmax + 1):
        rhs = [be.zero()] * NS
        for (r, c, t) in ser:
            kmax = min(n, len(t))
            acc = be.zero()
            for k in range(1, kmax + 1):
                acc += t[k - 1] * Y[n - k][c]
            rhs[r] += acc
        Y.append(be.trim_vec(be.solve_shifted(B0rows, n, rhs, NS)))
    sa = be.from_mpf(s_exit)
    acc = [be.zero()] * NS
    p = None
    for n in range(len(Y)):
        if p is None:
            p = sa * 0 + 1  # backend one
        for i in range(NS):
            acc[i] += Y[n][i] * p
        p = be.trim(p * sa)
    # coefficient max-norm profile (log10) for the measured-radius report (read at two dps, never one)
    prof = []
    for n in range(len(Y)):
        mx = max(abs(be.to_mpf(v)) for v in Y[n])
        prof.append(float(mp.log10(mx)) if mx != 0 else None)
    return be.trim_vec(acc), B0rows, prof


def _shift_poly(be, coeffs, x0):
    c = list(coeffs)[::-1]
    n = len(c)
    out = []
    for j in range(n):
        for i in range(1, n - j):
            c[i] += c[i - 1] * x0
        out.append(c[n - 1 - j])
        c = c[:n - 1 - j]
    return out


def _series_rat(be, pc, qc, x0, order):
    ps, qs = _shift_poly(be, pc, x0), _shift_poly(be, qc, x0)
    c = [be.zero()] * (order + 1)
    q0 = qs[0]
    for n in range(order + 1):
        t = ps[n] if n < len(ps) else be.zero()
        for k in range(1, min(n, len(qs) - 1) + 1):
            t -= qs[k] * c[n - k]
        c[n] = be.trim(t / q0)
    return c


POLES = []


def de_singular_points():
    """Real roots in s of every irreducible s-factor of the denominators of A(d,s) (the DE's singular points; s = 0 always included),
    exact (a rational, or a + b*sqrt(n)) from the parser's factorization of the connection strings (kite_exact_parse.singular_points;
    python-flint fmpz_poly.factor).  Factors that mix d and s are reported (none expected)."""
    try:
        return kx.singular_points(CONN["A"])
    except ImportError:
        sys.stderr.write(f"{SCRIPT}: python-flint is required for the singular-point census of the connection "
                         "(pip install python-flint); not serving\n")
        raise SystemExit(4)
    except kx.ParseError as ex:
        sys.stderr.write(f"{SCRIPT} REFUSED: the singular-point census of the connection: {ex}\n")
        raise SystemExit(3)


def make_transporter(be, ENT, NS, order):
    RATb = [(r, c, [be.from_rat(x) for x in pc], [be.from_rat(x) for x in qc])
            for (r, c, pc, qc) in ENT]
    poles = list(POLES)   # the DE singular points read from the connection denominators (de_singular_points)

    def step(y0, x0m, hm):
        x0 = be.from_mpf(x0m)
        h = be.from_mpf(hm)
        tc = [(r, cc, _series_rat(be, pc, qc, x0, order)) for (r, cc, pc, qc) in RATb]
        a = [list(y0)]
        for n in range(order + 1):
            ssum = [be.zero()] * NS
            for (r, cc, t) in tc:
                acc = be.zero()
                for k in range(n + 1):
                    acc += t[k] * a[n - k][cc]
                ssum[r] += acc
            a.append(be.trim_vec([v / (n + 1) for v in ssum]))
        y = [be.zero()] * NS
        hp = h * 0 + 1
        for n in range(order + 2):
            for i in range(NS):
                y[i] += a[n][i] * hp
            hp = be.trim(hp * h)
        return be.trim_vec(y)

    def transport(y0, x0, x1, safety=mp.mpf("0.25")):
        y, t, x1 = list(y0), mp.mpf(x0), mp.mpf(x1)
        dirn = 1 if x1 > t else -1
        while abs(x1 - t) > mp.mpf(10) ** (-mp.mp.dps + 8):
            r = min(abs(t - p) for p in poles)
            h = dirn * min(safety * r, abs(x1 - t))
            y = step(y, t, h)
            t += h
        return y

    return transport


# ---------------------------------------------------------------------------
# 4. Full pipeline at current mp.dps: derive the 42-component seed at s = -2
# ---------------------------------------------------------------------------
def run_pipeline(mutate=None, use_flint=True, s_exits=("-1/10",), points=("-2",)):
    """Returns dict with the transported 42-vectors at every point from every exit, the checks, the walls."""
    t0 = time.time()
    X = mp.mpf(X_STR)
    be = Backend(use_flint)
    be.setup()
    STATE, IDX, NS, ENT = build_graded()
    seedd = closed_form_seed(mutate=mutate)
    seq = [be.from_mpf(seedd[i][K]) for (i, K) in STATE]
    checks = {}
    if X == 2:
        checks["xi_check"] = ("xi(1,1,2)+8G", float(-mp.log10(abs(_xi(mp.mpf(1), mp.mpf(1), mp.mpf(2)) + 8 * mp.catalan) + mp.mpf(10) ** (-mp.mp.dps))))
    elif X == 3:
        checks["xi_check"] = ("xi(1,1,3)+(8/sqrt3)Cl2(pi/3)", float(-mp.log10(abs(_xi(mp.mpf(1), mp.mpf(1), mp.mpf(3)) + 8 / mp.sqrt(3) * _cl2(mp.pi / 3)) + mp.mpf(10) ** (-mp.mp.dps))))
    V22 = V_vacuum(0, 0, 5)
    with mp.extradps(20):
        ex = mp.taylor(lambda e: -mp.power(5, 1 - 2 * e) * mp.gamma(1 + e) * mp.gamma(1 - e) ** 2
                       * mp.gamma(1 + 2 * e) / (2 * (2 * e - 1) * mp.gamma(2 - e)), 0, 2)
    checks["V00M_vs_Gamma"] = float(min(-mp.log10(abs(V22[k] - ex[k + 2]) + mp.mpf(10) ** (-mp.mp.dps)) for k in (-2, -1, 0)))
    V010 = V_vacuum(0, 1, 0); V110 = V_vacuum(1, 1, 0); V01x = V_vacuum(0, 1, X); V11x = V_vacuum(1, 1, X)
    checks["top_seed_epsm1_cancel"] = float(-mp.log10(abs(V11x[-1] - V01x[-1] - V110[-1] + V010[-1]) + mp.mpf(10) ** (-mp.mp.dps)))
    t_build = time.time() - t0
    nmax = mp.mp.dps + 30
    order = max(60, int(1.7 * mp.mp.dps) + 20)
    transport = make_transporter(be, ENT, NS, order)
    out = {"exits": {}, "checks": checks, "timings": {"build": t_build}, "backend": "flint-arb" if be.flint else "mpmath", "nmax": nmax, "order": order}
    for sx in s_exits:
        s_exit = kx.mpf_rat(sx)
        t1 = time.time()
        yexit_b, B0rows, prof = frobenius_eval(be, STATE, NS, ENT, seq, s_exit, nmax)
        if "B0y0_max" not in checks:
            r = [sum((B0rows[i][j] * seq[j] for j in range(NS)), be.zero()) for i in range(NS)]
            checks["B0y0_max"] = float(max(abs(be.to_mpf(x)) for x in r))
        t_frob = time.time() - t1
        # measured growth: slope of log10|Y_n| over the last 20 orders -> rate per order; R_est = 10^-slope
        tail = [(n, v) for n, v in enumerate(prof) if v is not None][-21:]
        slope = (tail[-1][1] - tail[0][1]) / (tail[-1][0] - tail[0][0]) if len(tail) > 1 else None
        leg = {"s_exit": sx, "frobenius_s": t_frob, "coef_log10_profile_every10": {str(n): prof[n] for n in range(0, len(prof), 10)},
               "growth_slope_last20": slope, "R_est_from_slope": (10 ** (-slope) if slope is not None else None), "points": {}}
        y = list(yexit_b); t = s_exit
        for pt in points:
            sv = kx.mpf_rat(pt); t2 = time.time()
            y = transport(y, t, sv); t = sv
            vec = {(i, K): be.to_mpf(y[IDX[(i, K)]]) for (i, K) in STATE}
            leg["points"][pt] = {"vector": vec, "wall_s": time.time() - t2}
        out["exits"][sx] = leg
    out["timings"]["total"] = time.time() - t0
    return out


def digits_agree(a, ref):
    """Agreement digits vs a reference; ref == 0 (or absent) -> absolute digits
    (the component is exactly zero analytically; report -log10|computed|)."""
    a = mp.mpf(a)
    ref = mp.mpf(ref) if ref is not None else mp.mpf(0)
    if ref == 0:
        if a == 0:
            return float(mp.mp.dps)
        return float(-mp.log10(abs(a)))
    if a == ref:
        return float(mp.mp.dps)
    return float(-mp.log10(abs((a - ref) / ref)))


def main():
    global CONN, GATE_JSON, X_STR, NEG_M10, POLES, GRADED_PATH
    ap = argparse.ArgumentParser(
        description="Unequal-mass kite on the (1,1,x) line: the s = 0 vacuum seed, the Frobenius bounded branch "
                    "and the certified transport to every --points value (the (1,1,x)-line form of kite-boundary.py)")
    ap.add_argument("--conn", required=True,
                    help="the exact rational connection JSON (kite-connection-A14-x3.json, or "
                         "kite-connection-A14.json for the x = 2 control)")
    ap.add_argument("--x", default="3", help="the third squared mass x of the (1,1,x) line (default 3)")
    ap.add_argument("--dps", type=int, default=30,
                    help="target precision in decimal digits (default 30; the run works at dps+25)")
    ap.add_argument("--s-exits", default="-1/10,-1/40",
                    help="comma-separated Frobenius exit points inside the convergence disk at s = 0 "
                         "(default -1/10,-1/40; write --s-exits=-1/10)")
    ap.add_argument("--points", default="-5/3,-2,-10/3",
                    help="comma-separated transport targets, visited in order from each exit "
                         "(default -5/3,-2,-10/3; write --points=-2)")
    ap.add_argument("--ref", default=None,
                    help="kite-boundary-derived.json-schema file at s = -2 (the independent-reference gate of "
                         "the x = 2 control)")
    ap.add_argument("--receipt", required=True,
                    help="path of the JSON receipt to write (every number, the input hashes, the walls)")
    ap.add_argument("--no-flint", action="store_true",
                    help="use the pure-mpmath transport backend (the singular-point census still needs python-flint)")
    ap.add_argument("--neg-m10", action="store_true",
                    help="negative control: the alternative M10 coefficient 2x (the kernel check must flag it)")
    ap.add_argument("--mutate", action="store_true",
                    help="control: perturb the top-seed component by 1e-30 (the gate must collapse)")
    args = ap.parse_args()
    import subprocess
    sha = lambda q: hashlib.sha256(open(q, "rb").read()).hexdigest()
    mp.mp.dps = args.dps + 25
    X_STR = args.x; NEG_M10 = args.neg_m10; GATE_JSON = args.ref
    if not os.path.isfile(args.conn):
        sys.stderr.write(f"{SCRIPT}: --conn file {args.conn} not found; not serving\n")
        raise SystemExit(2)
    CONN = json.load(open(args.conn))
    _graded = (args.conn[:-5] if args.conn.endswith(".json") else args.conn) + "-graded.json"
    GRADED_PATH = _pinned_path(_graded) if os.path.exists(_graded) else None
    pts, mixed = de_singular_points()
    POLES = [r.mpf(mp.mp.dps + 10) for r in pts]   # each root evaluated to dps + 10 digits, then at the working precision
    mutate = (13, 0, "1e-30") if args.mutate else None
    print(f"x={args.x} dps={args.dps} (+25 guard) conn={args.conn} backend={'mpmath' if (args.no_flint or not HAVE_FLINT) else 'flint-arb'}{'  [NEG-M10]' if NEG_M10 else ''}{'  [MUTATED]' if mutate else ''}")
    print("DE singular points (from the connection denominators):", [str(r) for r in pts], "mixed d,s factors:", mixed)
    r = run_pipeline(mutate=mutate, use_flint=not args.no_flint, s_exits=args.s_exits.split(","), points=args.points.split(","))
    print("checks:", r["checks"])
    ref, ref_dps = load_gate()
    rec = {"receipt": f"TRANSPORT (1,1,{args.x}) at dps {args.dps} (working {mp.mp.dps}): closed-form s=0 seed -> Frobenius bounded branch -> exact-DE Taylor transport",
           "stamp_utc": subprocess.run(["date", "-u", "+%Y-%m-%dT%H:%M:%SZ"], capture_output=True, text=True).stdout.strip(),
           "conn": args.conn, "conn_sha256": sha(args.conn), "x": args.x, "dps": args.dps, "mp_dps_working": mp.mp.dps, "backend": r["backend"],
           "nmax_frobenius": r["nmax"], "taylor_order": r["order"], "neg_m10_control": NEG_M10, "mutated": bool(mutate),
           "de_singular_points": [str(q) for q in pts], "de_singular_points_float": [float(q) for q in pts], "mixed_d_s_denominator_factors": mixed,
           "checks": r["checks"], "timings": r["timings"], "exits": {}}
    exits = list(r["exits"].keys())
    for sx, leg in r["exits"].items():
        E = {k: v for k, v in leg.items() if k != "points"}; E["points"] = {}
        for pt, P in leg["points"].items():
            vec = P["vector"]
            E["points"][pt] = {"wall_s": P["wall_s"], "M13_eps0": mp.nstr(vec[(13, 0)], args.dps), "N_kite_like_-8_M13_eps0": mp.nstr(-8 * vec[(13, 0)], args.dps),
                               "vector": {f"{i},{K}": mp.nstr(vec[(i, K)], args.dps) for (i, K) in sorted(vec)}}
            if pt == "-2" and ref is not None:
                worst = None
                for (i, K), v in vec.items():
                    d = digits_agree(v, ref.get((i, K))); worst = d if worst is None else min(worst, d)
                E["points"][pt]["gate_vs_ref_worst_of_42_digits"] = worst
                E["points"][pt]["gate_vs_ref_N_kite_digits"] = digits_agree(-8 * vec[(13, 0)], -8 * mp.mpf(ref[(13, 0)]))
                print(f"exit {sx}: GATE at s=-2 vs {os.path.basename(args.ref)} ({ref_dps} d): worst of 42 = {worst:.1f} d; N_kite {E['points'][pt]['gate_vs_ref_N_kite_digits']:.1f} d")
        rec["exits"][sx] = E
    # exit-independence witness (the measured-radius question): the same point reached from two exits
    if len(exits) >= 2:
        rec["exit_independence_digits"] = {}
        for pt in r["exits"][exits[0]]["points"]:
            v0 = r["exits"][exits[0]]["points"][pt]["vector"]; v1 = r["exits"][exits[1]]["points"][pt]["vector"]
            worst = min(digits_agree(v0[k], v1[k]) if abs(v1[k]) > mp.mpf(10) ** (-(args.dps - 15)) else digits_agree(v0[k], None) for k in v0)
            rec["exit_independence_digits"][pt] = worst
            print(f"exit-independence at s={pt}: exits {exits[0]} vs {exits[1]} agree to {worst:.1f} d (worst of 42)")
    for sx, leg in r["exits"].items():
        print(f"exit {sx}: frobenius {leg['frobenius_s']:.1f} s, growth slope (last 20 orders) {leg['growth_slope_last20']}, R_est {leg['R_est_from_slope']}; "
              + "; ".join(f"s={pt}: {P['wall_s']:.1f} s, M13[eps0]={mp.nstr(P['vector'][(13,0)], 25)}" for pt, P in leg["points"].items()))
    print(f"timings: build {r['timings']['build']:.1f} s, TOTAL {r['timings']['total']:.1f} s ({r['backend']})")
    rec["PRODUCER"] = {"script": os.path.abspath(__file__), "sha256": sha(__file__), "copy_of": "kite-boundary.py", "stamp_source": "date -u",
                       "parser": os.path.basename(kx.__file__), "parser_sha256": sha(kx.__file__),
                       "graded_lists": (os.path.basename(GRADED_PATH) + " " + sha(GRADED_PATH) + " (checked against the connection strings)")
                       if GRADED_PATH else "the parser's reduced form (no graded file beside --conn)"}
    json.dump(rec, open(args.receipt, "w"), indent=1)
    print("receipt ->", args.receipt)


if __name__ == "__main__":
    main()