"""Simulation evidence for proof.md of 1RB1LC_1RC1LB_1LD0RB_0LE1RD_---0LF_0LA1LF.

Every numbered claim in proof.md has a section here with the same label. Run:
    ulimit -v 6000000; python3 -B check_proof.py [--quick] [--seed S]
Exit code 0 iff every check passed.  Stdlib + a C compiler (sections B2 and H compile scripts/edge_scan.c;
section B3 optionally uses tools/sim.py's C core).

Sections
  A  TM-level rules (direct TM simulation, random parameters and random surroundings)
  B  run model == TM: B1 macro step by macro step (pure Python, from step 50); B2 every left-edge
     configuration of the raw TM up to step 4*10^8 (scripts/edge_scan.c) == the model's edge returns,
     exact step counts; B3 spot checks against the tools/sim.py C core
  C  initial segment: blank tape -> step 50 = edge word [1]
  D  epoch rules in the run model (U0/U1/U2, TR, TL, TF, T, X, S), on marker tails and arbitrary tails
  E  the Theorem: E1 Inv(M,n) and the absorption law along real orbits; E2 the induction step on its own
     (adversarial: every / random / extreme words satisfying Inv(M,n), orbit words or not)
  F  top level: O_n([]) satisfies Inv, every epoch costs 2
  G  halting transition: E0 is read only in the HALT pattern
  H  twin (Eqv.v TM113): blank tape in state B -> conf([],[2,2,1]) at step 44; budget-2 edge lemma;
     every left-edge configuration up to step 4*10^8 is conf([],[2,2]+O_n([]))
"""
import os
import random
import subprocess
import sys
import tempfile
import time

TM = "1RB1LC_1RC1LB_1LD0RB_0LE1RD_---0LF_0LA1LF"
QUICK = "--quick" in sys.argv
SEED = int(sys.argv[sys.argv.index("--seed") + 1]) if "--seed" in sys.argv else 2026
FAIL = []


def ok(label, cond, info=""):
    if not cond:
        FAIL.append((label, info))
        if sum(1 for f in FAIL if f[0] == label) <= 3:
            print("  FAIL", label, info)


# ----------------------------------------------------------------------------------------------
# A tiny exact TM simulator (tape = dict, head position, state).  States A..F = 0..5.
def parse_tm(s):
    tab = {}
    for q, blk in enumerate(s.split("_")):
        for sym in (0, 1):
            t = blk[3 * sym:3 * sym + 3]
            tab[(q, sym)] = None if t[2] == "-" else (int(t[0]), 1 if t[1] == "R" else -1, ord(t[2]) - 65)
    return tab


TAB = parse_tm(TM)


class Tape:
    """cells: dict pos->1 (only 1s stored); head; state; steps."""

    def __init__(self, ones, head, state):
        self.ones = set(ones)
        self.head = head
        self.state = state
        self.steps = 0
        self.halted = False
        self.halt_read = None

    def run(self, n, stop_on_state_visit=None):
        ones, h, q = self.ones, self.head, self.state
        for _ in range(n):
            s = 1 if h in ones else 0
            t = TAB[(q, s)]
            if t is None:
                self.halted = True
                self.halt_read = (q, s)
                break
            w, d, q2 = t
            if w:
                ones.add(h)
            else:
                ones.discard(h)
            h += d
            q = q2
            self.steps += 1
        self.head, self.state = h, q


def build(blocks, head_block_end=None):
    """blocks: list of (bit, count) left to right starting at position 0.  Returns (ones, pos_of_end_of_block)."""
    ones = []
    pos = 0
    ends = []
    for bit, cnt in blocks:
        for _ in range(cnt):
            if bit:
                ones.append(pos)
            pos += 1
        ends.append(pos - 1)
    return ones, ends


def rand_word(rng, maxlen=4, maxv=6):
    """random right-hand context: runs >=1 and gaps >=1, as blocks (starting with a 1-run)."""
    out = []
    for _ in range(rng.randint(0, maxlen)):
        out += [(1, rng.randint(1, maxv)), (0, rng.randint(1, 4))]
    return out


def rand_left(rng):
    """random left context: runs >= 1 separated by single 0s, preceded by 0^inf padding."""
    out = [(0, 3)]
    for _ in range(rng.randint(0, 3)):
        out += [(1, rng.randint(1, 6)), (0, 1)]
    return out


def tm_rule_check(label, before, after, head_before_idx, head_after_idx, steps, trials_desc):
    """before/after: block lists; head_*_idx = index of the block whose LAST cell the head reads (state A).
    Blocks are aligned at position 0 and must have equal total length (padding added by caller)."""
    ones_b, ends_b = build(before)
    ones_a, ends_a = build(after)
    t = Tape(ones_b, ends_b[head_before_idx], 0)
    t.run(steps)
    good = (not t.halted and t.state == 0 and t.head == ends_a[head_after_idx] and t.ones == set(ones_a))
    ok(label, good, trials_desc)
    return good


# ----------------------------------------------------------------------------------------------
def section_A(rng):
    print("A  TM-level rules")
    N = 120 if QUICK else 500
    cnt = {"R1": 0, "R3": 0, "R2": 0, "R2e": 0}
    for _ in range(N):
        n, m, x = rng.randint(0, 25), rng.randint(0, 25), rng.randint(0, 25)
        L = rand_left(rng)
        r = rand_word(rng)
        # R1:  l 1^(n+2) <A 0 1^(m+1) 0 0 1^x 0 r  -->  l 1^(n+1) 0 1^(m+1) <A 0 1^(x+1) 0 0 r   (3m+2x+15)
        before = L + [(1, n + 2), (0, 1), (1, m + 1), (0, 2), (1, x), (0, 1)] + r + [(0, 3)]
        after = L + [(1, n + 1), (0, 1), (1, m + 1), (0, 1), (1, x + 1), (0, 2)] + r + [(0, 3)]
        hb = len(L)
        cnt["R1"] += tm_rule_check("A1 R1", before, after, hb, hb + 2, 3 * m + 2 * x + 15, (n, m, x))
        # R3:  0^inf <A 0 1^(k+1) 0 0 1^m 0 r  -->  0^inf 1^(k+2) <A 0 1^(m+1) 0 0 r   (3k+2m+11)
        k = n
        before = [(0, 4), (0, 1), (1, k + 1), (0, 2), (1, m), (0, 1)] + r + [(0, 3)]
        after = [(0, 3), (1, k + 2), (0, 1), (1, m + 1), (0, 2)] + r + [(0, 3)]
        ones_b, ends_b = build(before)
        ones_a, ends_a = build(after)
        t = Tape(ones_b, ends_b[0], 0)
        t.run(3 * k + 2 * m + 11)
        g = not t.halted and t.state == 0 and t.head == ends_a[1] and t.ones == set(ones_a)
        ok("A2 R3", g, (k, m))
        cnt["R3"] += g
        # R2:  l 1^(j+1) 0 1^(k+1) 0 1 <A 0 r  -->  l 1^(j+1) <A 0 1^(k+2) 0 0 r   (k+10)
        j = rng.randint(0, 25)
        before = L + [(1, j + 1), (0, 1), (1, k + 1), (0, 1), (1, 1), (0, 1)] + r + [(0, 3)]
        after = L + [(1, j + 1), (0, 1), (1, k + 2), (0, 2)] + r + [(0, 3)]
        ones_b, ends_b = build(before)
        ones_a, ends_a = build(after)
        t = Tape(ones_b, ends_b[len(L) + 4], 0)
        t.run(k + 10)
        g = not t.halted and t.state == 0 and t.head == ends_a[len(L)] and t.ones == set(ones_a)
        ok("A3 R2", g, k)
        cnt["R2"] += g
        # R2e: 0^inf 1^(k+1) 0 1 <A 0 r  -->  0^inf <A 0 1^(k+2) 0 0 r   (k+10)
        before = [(0, 5), (1, k + 1), (0, 1), (1, 1), (0, 1)] + r + [(0, 3)]
        after = [(0, 5), (1, k + 2), (0, 2)] + r + [(0, 3)]
        ones_b, ends_b = build(before)
        ones_a, ends_a = build(after)
        t = Tape(ones_b, ends_b[3], 0)
        t.run(k + 10)
        g = not t.halted and t.state == 0 and t.head == ends_a[0] - 1 and t.ones == set(ones_a)
        ok("A4 R2e", g, k)
        cnt["R2e"] += g
    # HALT: 0^inf 1 <A 0 r  halts reading E0 after 3 transitions
    hc = 0
    for _ in range(40):
        r = rand_word(rng)
        ones_b, ends_b = build([(0, 5), (1, 1), (0, 1)] + r + [(0, 3)])
        t = Tape(ones_b, ends_b[1], 0)
        t.run(10)
        g = t.halted and t.steps == 3 and t.halt_read == (4, 0)
        ok("A5 HALT", g)
        hc += g
    print(f"   R1 {cnt['R1']}/{N}  R3 {cnt['R3']}/{N}  R2 {cnt['R2']}/{N}  R2e {cnt['R2e']}/{N}  HALT {hc}/40"
          "  (n,m,x,k in 0..25, random left/right context)")


# ----------------------------------------------------------------------------------------------
# Run model.  L = list of runs left of A (top = L[-1]); R = [c1, g1, c2, g2, ..., cr].
def mstep(L, R):
    """one macro step, in place.  returns (rule, tm_steps)."""
    if not L or L[-1] >= 2:
        c1 = R[0]
        if len(R) >= 3 and R[1] == 2:
            x = R[2]
            newR = [x + 1] + ([R[3] + 1] + R[4:] if len(R) >= 4 else [])
        else:
            x = 0
            newR = [1] + ([R[1] - 1] + R[2:] if len(R) >= 2 else [])
        if not L:
            L.append(c1 + 1)
            rule, st = "R3", 3 * (c1 - 1) + 2 * x + 11
        else:
            L[-1] -= 1
            L.append(c1)
            rule, st = "R1", 3 * (c1 - 1) + 2 * x + 15
        R[:] = newR
        return rule, st
    if len(L) == 1:
        return "HALT", 3
    a = L[-2]
    del L[-2:]
    R[:] = [a + 1, 2] + R
    return ("R2" if L else "R2e"), a + 9


def model_ok(L, R):
    return (all(v >= 1 for v in L) and len(R) % 2 == 1 and all(v >= 1 for v in R)
            and (len(R) < 2 or R[1] >= 2))


def to_model(ones, head):
    """TM tape (state A, head reading the last cell of the left part) -> (L, R) or None."""
    lo, hi = min(ones | {head}), max(ones | {head})
    left = "".join("1" if p in ones else "0" for p in range(lo, head + 1)).lstrip("0")
    right = "".join("1" if p in ones else "0" for p in range(head + 1, hi + 2))
    L = [len(s) for s in left.split("0")] if left else []
    if any(v == 0 for v in L) or not right.startswith("0"):
        return None
    right = right[1:].rstrip("0")
    R = []
    i = 0
    while i < len(right):
        j = i
        while j < len(right) and right[j] == "1":
            j += 1
        R.append(j - i)
        if j >= len(right):
            break
        k2 = j
        while right[k2] == "0":
            k2 += 1
        R.append(k2 - j)
        i = k2
    return L, R


def section_C():
    print("C  initial segment")
    t = Tape([], 0, 0)
    t.run(50)
    lr = to_model(t.ones, t.head) if t.state == 0 else None
    ok("C1 step50", lr == ([], [3, 2, 1]), (t.state, lr))
    print(f"   blank --50 steps--> state {'ABCDEF'[t.state]}, (L,R) = {lr}  i.e. 0^inf <A 0 1^3 0^2 1 0^inf = edge word [1]")
    return t


def section_B(t50):
    print("B  run model == TM")
    nmac = 3000 if QUICK else 40000
    t = Tape(set(t50.ones), t50.head, t50.state)
    L, R = [], [3, 2, 1]
    kinds = {}
    good = True
    for i in range(nmac):
        rule, st = mstep(L, R)
        ok("B1 no-halt", rule != "HALT")
        kinds[rule] = kinds.get(rule, 0) + 1
        t.run(st)
        lr = to_model(t.ones, t.head) if (t.state == 0 and not t.halted) else None
        if lr != (L, R) or not model_ok(L, R):
            ok("B1 model==TM", False, (i, rule, L, R[:8], lr))
            good = False
            break
    print(f"   B1 {nmac if good else i} macro steps from step 50 identical to the TM (TM step {50 + t.steps:,});"
          f" rules used {kinds}")
    # B2: EVERY left-edge configuration of the raw TM (state A, all cells at or left of the head blank)
    # up to TM step `limit`, found by a direct C scan, against the run model's edge returns.
    limit = 3 * 10 ** 7 if QUICK else 4 * 10 ** 8
    hits = edge_scan("A", limit)
    pred = [(0, "blank")] + model_edges([], [2], 10, limit)      # conf([],[2]) at step 10 (proof.md section 5)
    ok("B2 edge scan == model", hits == pred,
       next(((i, h, p) for i, (h, p) in enumerate(zip(hits, pred)) if h != p), (len(hits), len(pred))))
    edges = [(s, R) for s, R in pred if isinstance(R, list) and s >= 50]
    words = [R[2:] for s, R in edges]
    ok("B2 edge form", all(R[:2] == [3, 2] for s, R in edges))
    top = top_orbit(min(len(words) - 1, 30 if QUICK else 38))
    ok("B2 words == O_n([])", words[:len(top)] == top)
    print(f"   B2 raw TM to step {limit:,} (scripts/edge_scan.c): {len(hits)} left-edge configurations = blank (step 0),"
          f" conf([],[2]) (step 10) and Edge(O_n([])) for n = 0..{len(edges) - 1}, each at the model's exact TM step"
          f" (last {edges[-1][0]:,}); nothing else; words == O_n([]) for n <= {len(top) - 1}")
    # B3: spot checks against an independent engine (tools/sim.py C core), optional
    try:
        sys.path.insert(0, "[local working file]")
        import sim
        m = sim.parse(TM)
    except Exception as e:  # pragma: no cover
        print("   B3 skipped (tools/sim.py unavailable):", e)
        return
    checked = 0
    for n, (step, R) in enumerate(edges):
        if n < 6 or n % 8 == 0 or n == len(edges) - 1:
            res = sim.run(m, max_steps=step)
            c = res.final
            tape = c.tape if isinstance(c.tape, (bytes, bytearray)) else bytes(c.tape)
            ones = {c.lo + i for i, b in enumerate(tape) if b}
            lr = to_model(ones, c.pos) if (c.state == "A" and not res.halted) else None
            ok("B3 edge==sim.py", lr == ([], R), (n, step))
            checked += 1
    print(f"   B3 {checked} of those edge returns re-checked with the tools/sim.py C core (independent engine)")


def edge_scan(start, limit):
    """compile and run scripts/edge_scan.c; returns [(step, R or 'blank' or raw fields)]."""
    here = os.path.dirname(os.path.abspath(__file__))
    with tempfile.TemporaryDirectory() as td:
        exe = os.path.join(td, "edge_scan")
        subprocess.run(["cc", "-O2", "-o", exe, os.path.join(here, "scripts", "edge_scan.c")], check=True)
        out = subprocess.run([exe, TM, start, str(limit)], check=True, capture_output=True, text=True).stdout
    hits = []
    for line in out.splitlines():
        f = line.split()
        if f[0] == "HALT":
            hits.append(("HALT", int(f[1])))
            continue
        s, c1, c2, runs = int(f[0]), int(f[1]), int(f[2]), [int(v) for v in f[3:]]
        if c1 == 0 and c2 == 0 and not runs:
            hits.append((s, "blank"))
        elif c1 == 0 and c2 == 1:
            hits.append((s, runs))
        else:
            hits.append((s, ("odd", c1, c2, runs)))
    return hits


def model_edges(L, R, step, limit):
    """run the model from conf(L,R) at TM step `step`; list (step, R) at every L = [] up to TM step limit."""
    L, R = list(L), list(R)
    out = [(step, list(R))] if not L else []
    while True:
        rule, st = mstep(L, R)
        if rule == "HALT" or step + st > limit:
            return out
        step += st
        if not L:
            out.append((step, list(R)))


TOP = [[1]]


def top_orbit(nmax):
    """[O_0([]), ..., O_nmax([])], cached."""
    while len(TOP) <= nmax:
        TOP.append(Ep(TOP[-1])[1])
    return TOP[:nmax + 1]


# ----------------------------------------------------------------------------------------------
# Epochs
def Ep(W, cap=5 * 10 ** 7, count=False):
    """edge epoch: (L=[], R=[3,2]+W) -> (L=[], R=[r,2]+W'); returns (cost=5-r, W') [, macro steps]; None on HALT."""
    L, R = [], [3, 2] + list(W)
    k = 0
    while True:
        rule, _ = mstep(L, R)
        k += 1
        if rule == "HALT":
            return None
        if not L:
            assert len(R) >= 2 and R[1] == 2
            return ((5 - R[0], R[2:], k) if count else (5 - R[0], R[2:]))
        if k > cap:
            raise RuntimeError("cap")


def run_until(L, R, nsteps):
    for _ in range(nsteps):
        rule, _ = mstep(L, R)
        if rule == "HALT":
            return False
    return True


def Ivec(N):
    out, s = [], 0
    for i in range(0, len(N), 2):
        s += N[i] - 1
        out.append(s)
    return out


def is_marker_word(N):
    return len(N) % 2 == 0 and all(N[i] >= 2 and N[i + 1] == 1 for i in range(0, len(N), 2))


def rand_marker(rng, mmax=4, gaps=(2, 2, 4, 4, 6, 8, 10, 12, 14)):
    M = []
    for _ in range(rng.randint(0, mmax)):
        M += [rng.choice(gaps), 1]
    return M


def rand_tail(rng, kmax=3):
    """a random *tail* (proof.md section 2): [g_1, c_2, g_2, c_3, ...], even length, entries >= 1, g_1 >= 2."""
    N = []
    for i in range(rng.randint(0, kmax)):
        N += [rng.randint(2 if i == 0 else 1, 7), rng.randint(1, 6)]
    return N


def orbit(M, n):
    W = [1] + list(M)
    for _ in range(n):
        W = Ep(W)[1]
    return W


def section_D(rng):
    print("D  epoch rules (run model)")
    N = 150 if QUICK else 600
    # U0, U1, U2 (9 macro steps each)
    r = Ep([1], count=True)
    ok("D1 U0", r == (2, [2, 2, 1], 9), r)
    u1 = u2 = 0
    for _ in range(N):
        g = rng.randint(4, 30)
        W = [1, g, rng.randint(1, 5)]
        for _ in range(rng.randint(0, 3)):
            W += [rng.randint(1, 5), rng.randint(1, 5)]
        r = Ep(W, count=True)
        u1 += ok_ret("D1 U1", r == (2, [2, 2, 1, g - 2] + W[2:], 9), (W, r))
        b = rng.randint(1, 30)
        M = rand_marker(rng)
        r = Ep([1, 2, b] + M, count=True)
        u2 += ok_ret("D1 U2", r == (2, [b + 2, 2, 1] + M, 9), (b, M, r))
    print(f"   U0 ok; U1 {u1}/{N} (g in 4..30, arbitrary tail); U2 {u2}/{N} (b in 1..30, marker tail)")
    # TR: transient, 8 macro steps; TF: final, 1 macro step
    tr = 0
    for it in range(N):
        a, b = rng.randint(2, 40), rng.randint(1, 40)
        N_ = rand_tail(rng) if it % 2 else rand_marker(rng)
        L, R = [], [3, 2, a, 2, b] + N_
        run_until(L, R, 8)
        tr += ok_ret("D2 TR", (L, R) == ([2, a], [b + 2, 2, 1] + N_), (a, b, N_, L, R))
        z = rng.randint(2, 40)
        L, R = [2, 1], [z, 2] + N_
        mstep(L, R)
        ok("D4 TF", (L, R) == ([], [3, 2, z, 2] + N_))
    print(f"   TR {tr}/{N} (a,b in 2..40 / 1..40; half marker tails, half arbitrary tails); TF ok")
    # TL: translation lemma
    tl = 0
    pool = []
    for _ in range(N):
        M = rand_marker(rng, 3)
        pool.append(orbit(M, rng.randint(0, 9)))
    for _ in range(N // 3):  # arbitrary (non-orbit) words too
        W = [rng.randint(1, 6)]
        for _ in range(rng.randint(0, 4)):
            W += [rng.randint(1, 5), rng.randint(1, 6)]
        pool.append(W)
    for W in pool:
        e = Ep(W, count=True)
        if e is None:
            continue
        c, W2, k = e
        L0 = [rng.randint(1, 9) for _ in range(rng.randint(0, 3))]
        y = rng.randint(2, 9)
        z = c + 1 + rng.randint(0, 6)
        L, R = L0 + [y], [z, 2] + W
        depth_ok = True
        for i in range(k):
            mstep(L, R)
            if i < k - 1 and len(L) <= len(L0) + 1:
                depth_ok = False
        tl += ok_ret("D3 TL", depth_ok and (L, R) == (L0 + [y - 1], [z - c + 1, 2] + W2), (W, c, z, L0, y))
    print(f"   TL {tl}/{len(pool) - sum(1 for W in pool if Ep(W) is None)} (random stacks L0, y>=2, z>=c+1; orbit words and arbitrary words)")
    # X: exhaustion rule, 10 macro steps
    xc = 0
    for _ in range(N):
        a2, b2 = rng.randint(2, 40), rng.randint(1, 40)
        rest = []
        for _ in range(rng.randint(0, 3)):
            rest += [rng.randint(1, 8), rng.randint(1, 4)]
        L, R = [2, 2], [2, 2, a2, 2, b2] + rest
        run_until(L, R, 10)
        restp = ([rest[0] + 2] + rest[1:]) if rest else []
        xc += ok_ret("D5 X", (L, R) == ([], [2, 2, a2 + 2, 2, b2 + 2] + restp), (a2, b2, rest, L, R))
    print(f"   X {xc}/{N} (a''>=2, b''>=1, arbitrary tail; tail's first gap +2)")
    # T and X at epoch level (marker tails and arbitrary tails, W_i(N) = Ep^i([1]+N)), S rule
    tc = xcnt = sc = 0
    tot = totarb = 0
    for it in range(2 * N):
        a, b = rng.randint(2, 14), rng.randint(1, 16)
        arb = it % 2 == 1
        N_ = rand_tail(rng) if arb else rand_marker(rng, 3)
        W = [a, 2, b] + N_
        e = Ep(W)
        # replay prediction
        zeta, Wi, pred = b + 2, [1] + N_, None
        for i in range(a - 1):
            ei = Ep(Wi)
            if ei is None:
                pred = "inner-halt"
                break
            ci, Wn = ei
            if zeta < ci + 1:
                if i == a - 2 and zeta == 2 and Wi[1:2] == [2] and Wi[0] >= 2 and len(Wi) >= 3:
                    a2, b2, rest = Wi[0], Wi[2], Wi[3:]
                    pred = (3, [a2 + 2, 2, b2 + 2] + (([rest[0] + 2] + rest[1:]) if rest else []))
                    xcnt += 1
                else:
                    pred = "other-exhaustion"
                break
            zeta, Wi = zeta - ci + 1, Wn
        else:
            pred = (2, [zeta, 2] + Wi)
        if isinstance(pred, tuple):
            tot += 1
            totarb += arb
            tc += ok_ret("D6 T/X", e == pred, (W, e, pred))
        bS = rng.randint(1, 20)
        aS = rng.randint(1, 20)
        M = rand_tail(rng) if arb else rand_marker(rng)
        e = Ep([2, 2, aS, 2, bS] + M, count=True)
        sc += ok_ret("D7 S", e == (2, [aS + 1, 2, bS + 2, 2, 1] + M, 18), (aS, bS, M, e))
    print(f"   T/X epoch level {tc}/{tot} random T-words ({totarb} with arbitrary tails; {xcnt} by rule X);"
          f" S {sc}/{2 * N} (a,b in 1..20, half marker tails, half arbitrary tails; 18 macro steps)")


def ok_ret(label, cond, info=""):
    ok(label, cond, info)
    return 1 if cond else 0


# ----------------------------------------------------------------------------------------------
# The invariant Inv(M, n) of proof.md section 3
def inv_fail(M, I, n, W):
    m = len(I)
    alpha = sum(1 for v in I if v < n)
    Ia = I[alpha - 1] if alpha else 0
    js = alpha + 1 if alpha < m else None
    T = (n - Ia) % 2 == 1
    if T:
        if not (len(W) >= 3 and W[1] == 2 and W[0] >= 2):
            return "T-shape"
        a, b, N = W[0], W[2], W[3:]
        if b != n + alpha:
            return "b"
        off = 3
        absnow = js is not None and I[js - 1] == n
    else:
        if not (len(W) >= 5 and W[:2] == [2, 2] and W[3] == 2 and W[2] >= 2):
            return "S-shape"
        a, b, N = W[2], W[4], W[5:]
        if b != n - 1 + alpha:
            return "bS"
        off = 4
        absnow = js is not None and I[js - 1] == n + 1
    if not is_marker_word(N) or any(N[i] % 2 for i in range(0, len(N), 2)):
        return "marker-word/even"
    J = Ivec(N)
    aff = b + 1 - a
    if aff < -1 or (aff == -1 and not absnow):
        return "aff"
    if len(J) != aff + (m - alpha):
        return "count"
    if aff >= 1 and J[aff - 1] > a - 2:
        return "INT"
    if js is not None:
        for j in range(js, m + 1):
            idx = aff + 1 + (j - js)
            if j == js and absnow:
                if idx >= 1 and J[idx - 1] > a - 3:
                    return "REACH"
            elif idx < 1 or J[idx - 1] != a - off + I[j - 1] - n:
                return "EXT"
    return None


def check_orbit(M, nmax, stats):
    I = Ivec(M)
    W = [1] + list(M)
    for n in range(nmax + 1):
        if n > 0:
            f = inv_fail(M, I, n, W)
            stats["words"] += 1
            if f:
                return (n, f, W)
        e = Ep(W)
        if e is None:
            return (n, "HALT")
        c, W = e
        stats["epochs"] += 1
        if (c == 3) != (n in I):
            return (n, "cost", c)
        stats["c3"] += (c == 3)
    return None


def section_E(rng):
    print("E1 invariant Inv(M,n) + absorption law along real orbits (random even-gap marker words)")
    trials = 300 if QUICK else 2500
    nmax_cap = 22 if QUICK else 28
    stats = {"words": 0, "epochs": 0, "c3": 0, "T": 0}
    bad = 0
    t0 = time.time()
    for _ in range(trials):
        M = rand_marker(rng, 6)
        I = Ivec(M)
        nmax = min(nmax_cap, (I[-1] if I else 0) + rng.randint(2, 8))
        r = check_orbit(M, nmax, stats)
        if r:
            bad += 1
            ok("E Inv/AL", False, (M, r))
    # exhaustive small cases: every even-gap M with gaps in {2,4,6,8} and m <= 3
    ex = 0
    gs = [2, 4, 6, 8]
    words = [[]]
    for mm in range(1, 4):
        new = []
        for w in words:
            if len(w) // 2 == mm - 1:
                for g in gs:
                    new.append(w + [g, 1])
        words += new
    for M in words:
        I = Ivec(M)
        r = check_orbit(M, (I[-1] if I else 0) + 6, stats)
        ex += 1
        if r:
            bad += 1
            ok("E Inv/AL exhaustive", False, (M, r))
    print(f"   {trials} random + {ex} exhaustive (gaps 2..8, m<=3) orbits; {stats['words']} words checked against Inv,"
          f" {stats['epochs']} epochs (cost-3: {stats['c3']}); failures {bad}  [{time.time() - t0:.0f}s]")


def section_F():
    print("F  top level O_n([])")
    nmax = 30 if QUICK else 37
    top = top_orbit(nmax)
    hist = []
    for n, W in enumerate(top):
        if n > 0:
            f = inv_fail([], [], n, W)
            ok("F Inv(0,n)", f is None, (n, f, W))
            if n % 2 == 0 and n >= 6:
                hist.append(tuple(W[2:5]) + (tuple(W[5:]),))
        e = Ep(W)
        ok("F cost2", e is not None and e[0] == 2, (n, e))
    print(f"   n = 0..{nmax}: every epoch costs 2, Inv([],n) holds; even-n words (a,2,b|markers): {hist[:5]} ...")
    # the observer's Fs history: Fs(k,P,Q) = S(2k, 2k+1) + [P,1,Q,1]
    fs = [(w[0] // 2, w[3][0], w[3][2]) for w in hist if len(w[3]) == 4 and w[0] + 1 == w[2]]
    print(f"   Fs(k,P,Q) seen: {fs[:9]}  (observe history from step 18,826: (4,2,6) (5,6,4) (6,4,8) (7,4,10) (8,10,6))")
    ok("F history", fs[1:6] == [(4, 2, 6), (5, 6, 4), (6, 4, 8), (7, 4, 10), (8, 10, 6)], fs[:6])


# ----------------------------------------------------------------------------------------------
# E2: the induction step of the Theorem on its own.  The step uses only "W satisfies Inv(M,n)" plus the
# induction hypothesis on the inner orbits of W's marker tail N (real orbits, smaller measure).  So build
# words that satisfy Inv(M,n) directly from the clauses (they need not lie on any orbit) and check that
# Ep(W) costs 2 + [n in I(M)] and that the result satisfies Inv(M, n+1).
def inv_setup(M, n):
    I = Ivec(M)
    m = len(I)
    alpha = sum(1 for v in I if v < n)
    Ia = I[alpha - 1] if alpha else 0
    js = alpha + 1 if alpha < m else None
    T = (n - Ia) % 2 == 1
    if T:
        b, off = n + alpha, 3
        absnow = js is not None and I[js - 1] == n
    else:
        b, off = n - 1 + alpha, 4
        absnow = js is not None and I[js - 1] == n + 1
    return I, m, alpha, js, T, b, off, absnow


def inv_cons(M, n, a):
    """constraints on the marker positions I_1(N) < I_2(N) < ... of a word in Inv(M,n) with this a, as a list of
    (kind, value): 'int' (<= a-2, INT), 'reach' (<= a-3, the absorption exception), 'eq' (exact, EXT).
    None when the clauses allow no word with this a (aff < -1, or aff = -1 outside the allowed step)."""
    I, m, alpha, js, T, b, off, absnow = inv_setup(M, n)
    aff = b + 1 - a
    if aff < -1 or (aff == -1 and not absnow):
        return None
    cons = [("int", a - 2)] * max(aff, 0)
    if js is not None:
        for j in range(js, m + 1):
            t = aff + 1 + (j - js)
            if t >= 1:
                cons.append(("reach", a - 3) if (j == js and absnow) else ("eq", a - off + I[j - 1] - n))
    return cons


def inv_positions(cons, pick, rng=None):
    """yield increasing position lists p_1 < ... < p_K with p_t = t (mod 2) (even gaps) meeting cons.
    pick: 'all' (every list), 'rand' (one random list), or a pair (int_pick, reach_pick) of 'lo'/'hi'
    (extreme placements: everything packed at the bottom or at the top of its allowed range)."""
    K = len(cons)
    ub = [0] * K
    for t in range(K, 0, -1):
        kind, v = cons[t - 1]
        u = v if t == K else min(v, ub[t] - 1)
        if (u - t) % 2:
            u -= 1
        if kind == "eq" and u != v:
            return  # parity or order impossible
        ub[t - 1] = u
    if K and ub[0] < 1:
        return

    def rec(t, prev, acc):
        if t > K:
            yield list(acc)
            return
        kind, v = cons[t - 1]
        if kind == "eq":
            ch = [v]
        else:
            ch = list(range(prev + 1, ub[t - 1] + 1, 2))
            if pick == "rand":
                ch = [rng.choice(ch)]
            elif pick != "all":
                side = pick[0] if kind == "int" else pick[1]
                ch = [ch[0] if side == "lo" else ch[-1]]
        for p in ch:
            acc.append(p)
            yield from rec(t + 1, p, acc)
            acc.pop()

    yield from rec(1, 0, [])


def inv_word(M, n, a, P):
    I, m, alpha, js, T, b, off, absnow = inv_setup(M, n)
    N, prev = [], 0
    for p in P:
        N += [p - prev + 1, 1]
        prev = p
    return ([a, 2, b] if T else [2, 2, a, 2, b]) + N


EPC = {}


def ep_long(W):
    """Ep for long epochs: the Python model up to 2*10^6 macro steps, then the C port of the same model
    (scripts/epc.c, mode 'ep'), compiled once.  Section E2 cross-checks the two on short epochs."""
    try:
        return Ep(W, cap=2 * 10 ** 6)
    except RuntimeError:
        pass
    return ep_c(W)


def ep_c(W):
    if "exe" not in EPC:
        EPC["dir"] = tempfile.TemporaryDirectory()
        EPC["exe"] = os.path.join(EPC["dir"].name, "epc")
        here = os.path.dirname(os.path.abspath(__file__))
        subprocess.run(["cc", "-O2", "-o", EPC["exe"], os.path.join(here, "scripts", "epc.c")], check=True)
    EPC["calls"] = EPC.get("calls", 0) + 1
    out = subprocess.run([EPC["exe"], "ep"] + [str(v) for v in W], check=True, capture_output=True,
                         text=True).stdout.split()
    assert out[0] == "r" and out[2] == "W'"
    r = int(out[1])
    return None if r < 0 else (5 - r, [int(v) for v in out[3:]])


def step_check(M, n, W, st):
    """the induction step on one word W in Inv(M,n); returns a failure description or None."""
    I = Ivec(M)
    f = inv_fail(M, I, n, W)
    if f:
        return ("generator", f)
    e = ep_long(W)
    if e is None:
        return ("HALT",)
    c, W2 = e
    st["epochs"] += 1
    st["c3"] += c == 3
    if c != 2 + (n in I):
        return ("cost", c)
    f = inv_fail(M, I, n + 1, W2)
    if f:
        return ("Inv(n+1)", f, W2)
    _, m, alpha, js, T, b, off, absnow = inv_setup(M, n)
    aff = b + 1 - (W[0] if T else W[2])
    st["T" if T else "S"] += 1
    st["exception"] += absnow and aff >= 0
    st["merged"] += aff == -1
    return None


def section_E2(rng):
    print("E2 induction step alone: words built from the Inv(M,n) clauses (orbit words or not)")
    t0 = time.time()
    tot_bad = 0
    # (i) exhaustive: every Inv(M,n)-word with a <= amax, every even-gap M with gaps in {2,4,6,8} and m <= 3
    amax = 10 if QUICK else 14
    Ms = [[]]
    for mm in range(1, 4):
        Ms += [w + [g, 1] for w in Ms if len(w) // 2 == mm - 1 for g in (2, 4, 6, 8)]
    st = dict(epochs=0, c3=0, T=0, S=0, exception=0, merged=0)
    bad = 0
    for M in Ms:
        for n in range(1, 2 * amax + 1):
            for a in range(2, amax + 1):
                cons = inv_cons(M, n, a)
                if cons is None:
                    continue
                for P in inv_positions(cons, "all"):
                    W = inv_word(M, n, a, P)
                    r = step_check(M, n, W, st)
                    if r:
                        bad += 1
                        ok("E2 exhaustive", False, (M, n, W, r))
    print(f"   exhaustive (a <= {amax}, {len(Ms)} words M with gaps in 2..8, m <= 3): {st['epochs']} words "
          f"({st['T']} T, {st['S']} S; {st['exception']} absorption exceptions, {st['merged']} merged aff=-1), "
          f"{st['c3']} of cost 3; failures {bad}  [{time.time() - t0:.0f}s]")
    tot_bad += bad
    # (ii) random: random M (m <= 6, gaps 2..14), random n, random a <= 26, one random placement
    t1 = time.time()
    st = dict(epochs=0, c3=0, T=0, S=0, exception=0, merged=0)
    bad = 0
    target = 300 if QUICK else 2000
    amax = 18 if QUICK else 26
    while st["epochs"] < target:
        M = rand_marker(rng, 6)
        I = Ivec(M)
        n = rng.randint(1, (I[-1] if I else 0) + 10)
        a = rng.randint(2, amax)
        cons = inv_cons(M, n, a)
        if cons is None:
            continue
        for P in inv_positions(cons, "rand", rng):
            r = step_check(M, n, inv_word(M, n, a, P), st)
            if r:
                bad += 1
                ok("E2 random", False, (M, n, a, P, r))
    print(f"   random (a <= {amax}, m <= 6, gaps 2..14): {st['epochs']} words ({st['T']} T, {st['S']} S;"
          f" {st['exception']} exceptions, {st['merged']} merged), {st['c3']} of cost 3; failures {bad}"
          f"  [{time.time() - t1:.0f}s]")
    tot_bad += bad
    # (iii) extreme placements: internal markers packed at the bottom / top of their range, the within-reach
    # marker at its lowest / highest allowed position
    t1 = time.time()
    st = dict(epochs=0, c3=0, T=0, S=0, exception=0, merged=0)
    bad = 0
    target = 400 if QUICK else 3000
    amax = 16 if QUICK else 22
    while st["epochs"] < target:
        M = rand_marker(rng, 7)
        I = Ivec(M)
        n = rng.randint(1, (I[-1] if I else 0) + 6)
        seen = set()
        for a in range(2, amax + 1):
            cons = inv_cons(M, n, a)
            if cons is None:
                continue
            for pick in (("lo", "lo"), ("lo", "hi"), ("hi", "lo"), ("hi", "hi")):
                for P in inv_positions(cons, pick):
                    W = inv_word(M, n, a, P)
                    if tuple(W) in seen:
                        continue
                    seen.add(tuple(W))
                    r = step_check(M, n, W, st)
                    if r:
                        bad += 1
                        ok("E2 extreme", False, (M, n, W, r))
    print(f"   extreme placements (a <= {amax}, m <= 7): {st['epochs']} words ({st['T']} T, {st['S']} S;"
          f" {st['exception']} exceptions, {st['merged']} merged), {st['c3']} of cost 3; failures {bad}"
          f"  [{time.time() - t1:.0f}s]")
    tot_bad += bad
    print(f"   E2 total failures {tot_bad}: every word cost 2 + [n in I(M)] and landed in Inv(M, n+1)"
          f"  ({EPC.get('calls', 0)} epochs longer than 2*10^6 macro steps ran on the C port scripts/epc.c)")
    # the C port agrees with the Python model on short epochs (random Inv words and arbitrary words)
    agree = tot = 0
    for _ in range(60 if QUICK else 300):
        M = rand_marker(rng, 5)
        I = Ivec(M)
        n = rng.randint(1, (I[-1] if I else 0) + 6)
        a = rng.randint(2, 12)
        cons = inv_cons(M, n, a)
        W = None
        if cons is not None:
            for P in inv_positions(cons, "rand", rng):
                W = inv_word(M, n, a, P)
        if W is None:
            W = [rng.randint(1, 6)] + rand_tail(rng, 4)
        e = Ep(W)
        tot += 1
        agree += ok_ret("E2 epc==python", ep_c(W) == e, (W, e))
    print(f"   C port vs Python model on short epochs: {agree}/{tot} identical (cost and result word)")
    # (iv) sensitivity: move one clause's bound by 2 (INT, the reach exception, or an exact EXT position);
    # the words this adds lie outside Inv(M,n), and the step must visibly fail on them.  (A shift by 1 adds
    # nothing: positions have fixed parity.)
    t1 = time.time()
    res = []
    for kind, sh in (("int", 2), ("reach", 2), ("eq", 2), ("eq", -2)):
        outside = fails = 0
        for M in ([], [2, 1], [4, 1], [6, 1], [2, 1, 4, 1], [4, 1, 2, 1], [8, 1, 2, 1]):
            I = Ivec(M)
            for n in range(1, 16):
                for a in range(2, 11):
                    cons = inv_cons(M, n, a)
                    if cons is None:
                        continue
                    cons = [(k, v + (sh if k == kind else 0)) for k, v in cons]
                    for P in inv_positions(cons, "all"):
                        W = inv_word(M, n, a, P)
                        if inv_fail(M, I, n, W) is None:
                            continue
                        outside += 1
                        e = Ep(W)
                        fails += e is None or e[0] != 2 + (n in I) or inv_fail(M, I, n + 1, e[1]) is not None
        ok("E2 sensitivity", outside > 0 and fails > 0, (kind, sh, outside, fails))
        res.append(f"{kind}{sh:+d}: {fails}/{outside}")
    print(f"   sensitivity (clause bound moved by 2; words outside Inv whose step fails): {', '.join(res)}"
          f"  [{time.time() - t1:.0f}s]")


# ----------------------------------------------------------------------------------------------
# H: the twin.  Eqv.v TM113.tm = 1RB1LA_1LC0RA_0LD1RC_---0LE_0LF1LE_1RA1LB is this machine with states renamed
# (their A..F = our B,C,D,E,F,A), so its blank start is our blank tape in state B.
TWIN113 = "1RB1LA_1LC0RA_0LD1RC_---0LE_0LF1LE_1RA1LB"


def section_H(rng):
    print("H  twin: Eqv.v TM113.tm (= this TM started in state B; Eqv.v gives row 833)")
    # relabeling check: TM113.tm's state X is our state REN[X]
    ren = dict(zip("ABCDEF", "BCDEFA"))
    ours = TM.split("_")
    rel = []
    for X in "ABCDEF":
        blk = ours["ABCDEF".index(ren[X])]
        out = ""
        for k in (0, 3):
            t = blk[k:k + 3]
            out += t if t[2] == "-" else t[:2] + {v: u for u, v in ren.items()}[t[2]]
        rel.append(out)
    ok("H relabel", "_".join(rel) == TWIN113, "_".join(rel))
    # initial segment: blank, state B, 44 steps -> conf([], [2,2,1])
    t = Tape([], 0, 1)
    t.run(44)
    lr = to_model(t.ones, t.head) if t.state == 0 and not t.halted else None
    ok("H init44", lr == ([], [2, 2, 1]), lr)
    # budget-2 edge lemma (model): Ep(W) = (c, W') with c <= 2  =>  conf([],[2,2]+W) -> conf([],[4-c,2]+W'),
    # same k macro steps, L != [] in between.  (c = 3 would leave value 1 on the edge frame: HALT.)
    n2 = 0
    pool = top_orbit(24 if QUICK else 30) + [orbit(rand_marker(rng, 3), rng.randint(0, 9)) for _ in range(300)]
    for _ in range(200):
        W = [rng.randint(1, 6)]
        for _ in range(rng.randint(0, 4)):
            W += [rng.randint(1, 5), rng.randint(1, 6)]
        pool.append(W)
    for W in pool:
        e = Ep(W, count=True)
        if e is None or e[0] > 2:
            continue
        c, W2, k = e
        L, R = [], [2, 2] + W
        inner_ok = True
        for i in range(k):
            rule, _ = mstep(L, R)
            if rule == "HALT" or (i < k - 1 and not L):
                inner_ok = False
                break
        n2 += ok_ret("H edge2", inner_ok and (L, R) == ([], [4 - c, 2] + W2), (W, c))
    print(f"   relabeling ok; blank in state B --44 steps--> conf([],[2,2,1]); budget-2 edge lemma {n2} cases"
          f" (top orbit, marker orbits, arbitrary words with cost <= 2)")
    # every left-edge configuration of the raw TM from blank in state B up to `limit`
    limit = 3 * 10 ** 7 if QUICK else 4 * 10 ** 8
    hits = edge_scan("B", limit)
    pred = model_edges([], [1], 8, limit)    # conf([],[1]) at step 8 (R3 R1 R2e -> step 44)
    ok("H scan == model", hits == pred,
       next(((i, h, p) for i, (h, p) in enumerate(zip(hits, pred)) if h != p), (len(hits), len(pred))))
    edges = [(s, R) for s, R in pred if s >= 44]
    words = [R[2:] for s, R in edges]
    ok("H edge form", all(R[:2] == [2, 2] for s, R in edges))
    top = top_orbit(len(words) - 1)
    ok("H words == O_n([])", words == top)
    print(f"   raw TM from blank in state B to step {limit:,}: {len(hits)} left-edge configurations = conf([],[1])"
          f" (step 8) and conf([],[2,2]+O_n([])) for n = 0..{len(edges) - 1} (last at step {edges[-1][0]:,}),"
          f" exact model step counts; nothing else")


def section_G():
    print("G  E0 (the undefined transition) is entered only from D0; scan the TM run for E-reads-0")
    t = Tape([], 0, 0)
    n = 2 * 10 ** 5 if QUICK else 2 * 10 ** 6
    ones, h, q = t.ones, 0, 0
    e_reads = 0
    for _ in range(n):
        s = 1 if h in ones else 0
        tr = TAB[(q, s)]
        if q == 4:
            e_reads += 1
            ok("G E reads 1", s == 1)
        if tr is None:
            ok("G no halt", False)
            break
        w, d, q2 = tr
        (ones.add(h) if w else ones.discard(h))
        h += d
        q = q2
    print(f"   {n:,} TM steps from blank: E entered {e_reads:,} times, always reading 1 (E1 = 0LF); no halt")


if __name__ == "__main__":
    rng = random.Random(SEED)
    t0 = time.time()
    section_A(rng)
    t50 = section_C()
    section_B(t50)
    section_D(rng)
    section_F()
    section_G()
    section_H(rng)
    section_E(rng)
    section_E2(rng)
    print(f"\n{'ALL CHECKS PASSED' if not FAIL else 'FAILURES: ' + str(len(FAIL))}  (seed {SEED}, {'quick' if QUICK else 'full'}, {time.time() - t0:.0f}s)")
    sys.exit(1 if FAIL else 0)