"""check_proof.py: simulation checks for every claim in proof.md (machine #10).

    1RB0RF_1RC1RF_1LD1RA_---0LE_1LA0LE_0RA0LC

usage:  python3 -B check_proof.py [--quick] [--seed S]
        (from anywhere; imports scripts/deepen_model.py next to this file)

Each check prints one line "[PASS|FAIL|INFO] <id> <what> : <numbers>" and a JSON summary at the end.
INFO lines report data and never fail. Simulation is evidence, not proof; the ids match the section numbers
of proof.md (revision 2, after the skeptic review of 2026-10-10).
"""
import os, sys, json, random, re, time
HERE = os.path.dirname(os.path.abspath(__file__))
sys.path.insert(0, os.path.join(HERE, 'scripts'))
from deepen_model import (trim, carry, inc, eat, epoch, bounce_phase, word as dword,
                          Tstr, sym, regroup_tail, normalized, parse)

QUICK = '--quick' in sys.argv
SEED = int(sys.argv[sys.argv.index('--seed') + 1]) if '--seed' in sys.argv else 1
random.seed(SEED)
TM = '1RB0RF_1RC1RF_1LD1RA_---0LE_1LA0LE_0RA0LC'
TR = {}
for i, g in enumerate(TM.split('_')):
    for s in (0, 1):
        t = g[3 * s:3 * s + 3]
        TR[(i, s)] = None if t == '---' else (int(t[0]), 1 if t[1] == 'R' else -1, ord(t[2]) - 65)
TAB = [TR[(q, s)] for q in range(6) for s in (0, 1)]
RESULTS = {}


def report(cid, ok, what, info):
    RESULTS[cid] = {'pass': bool(ok), 'what': what, 'info': info}
    print('[%s] %s %s : %s' % ('PASS' if ok else 'FAIL', cid, what, info), flush=True)


# ---------------------------------------------------------------- TM simulator
def run(left, state, right, stop, maxs):
    """Tape = 0^inf left right 0^inf, head on right[0] in `state` (letter).
    stop(kind, pos, q) is called on 'C>' (C entered by a right move) and 'LREC' (new leftmost cell).
    Returns (status, steps, tape, head, q, base)."""
    pad = 4 * (len(left) + len(right)) + 64
    t = bytearray(pad) + bytearray(int(c) for c in left + right) + bytearray(pad)
    h = pad + len(left); base = h; q = ord(state) - 65; lo = pad; step = 0
    while step < maxs:
        if h <= 2 or h >= len(t) - 3:
            ext = bytearray(len(t)); t = ext + t + ext; h += len(ext); base += len(ext); lo += len(ext)
        tr = TAB[2 * q + t[h]]
        if tr is None:
            return 'HALT', step, t, h, q, base
        t[h] = tr[0]; h += tr[1]; q = tr[2]; step += 1
        if h < lo:
            lo = h
            if stop('LREC', h - base, q):
                return 'LREC', step, t, h, q, base
        elif q == 2 and tr[1] == 1 and stop('C>', h - base, q):
            return 'C>', step, t, h, q, base
    return 'MAX', step, t, h, q, base


def cells(t, a, b):
    return ''.join('01'[x] for x in t[a:b])


def Khub(a, b):
    return '110' * a + '0' + '011' * b


# ------------------------------------------------- 2.1 low-level shift rules
def block_run(q, cl, h, limit=10 ** 5):
    t = list(cl); s = 0; q = ord(q) - 65
    while 0 <= h < len(t) and s < limit:
        tr = TR[(q, t[h])]
        if tr is None:
            return None
        t[h] = tr[0]; h += tr[1]; q = tr[2]; s += 1
    return 'ABCDEF'[q], ('R' if h >= len(t) else 'L'), ''.join(map(str, t)), s


SHIFT = [  # (state, dir, block, block', state', steps per block)
    ('A', 'R', '001', '111', 'A', 3), ('A', 'R', '010', '110', 'A', 3),
    ('B', 'R', '010', '111', 'B', 3), ('F', 'R', '001', '011', 'F', 3),
    ('E', 'L', '110', '010', 'E', 5), ('E', 'L', '1', '0', 'E', 1)]
bad = 0; n_inst = 0
for (q, d, w, w2, q2, spk) in SHIFT:
    for k in range(0, 41):
        cl = [int(x) for x in w * k]
        if k == 0:
            continue
        r = block_run(q, cl, 0 if d == 'R' else len(cl) - 1)
        n_inst += 1
        if r is None or r != (q2, d, w2 * k, spk * k):
            bad += 1
report('2.1a', bad == 0, 'shift rules A>001^k, A>010^k, B>010^k, F>001^k, 110^k<E, 1^k<E (k=1..40)',
       '%d/%d instances' % (n_inst - bad, n_inst))
# turn rules used by the carry: A> 000 -> <E 001 (5 steps); B> 000 -> <E 010 (3); 111 <E -> <E 000 (3)
TURN = [('A', 0, '000', ('E', 'L', '001', 5)), ('B', 0, '000', ('E', 'L', '010', 3)),
        ('E', 2, '111', ('E', 'L', '000', 3)), ('A', 0, '100', ('B', 'R', '001', 3)),
        ('B', 0, '100', ('B', 'R', '101', 3))]
okt = all(block_run(q, [int(x) for x in w], h) == exp for (q, h, w, exp) in TURN)
report('2.1b', okt, 'group turns A>000-><E 001, B>000-><E 010, 111<E-><E 000, A>100->001 B>, B>100->101 B>',
       '%d/%d' % (sum(block_run(q, [int(x) for x in w], h) == exp for (q, h, w, exp) in TURN), len(TURN)))
# left pass over the bouncer blocks: 0 (011)^k <D 1  ->  <A 0 1 (010)^(k-1) 001, 5k-2 steps
bad = 0
for k in range(1, 41):
    cl = [int(x) for x in '011' * k + '1']
    r = block_run('D', cl, 3 * k - 1)
    if r != ('A', 'L', '1' + '010' * (k - 1) + '001', 5 * k - 2):
        bad += 1
report('2.1c', bad == 0, 'left pass (011)^k <D 1 -> <A 1 (010)^(k-1) 001 in 5k-2 steps (k=1..40)', '%d/40' % (40 - bad))

# ------------------------------------------------- 2.2 high-level rules R, B, E vs the TM
NR = 150 if QUICK else 1000
G3 = ['000', '001', '010']; GT = ['100', '010', '000']
cnt = {'R': [0, 0], 'B': [0, 0], 'E': [0, 0]}; fails = []
for _ in range(NR):   # R restart
    n = random.choice([1, 2, 3, random.randint(1, 40), random.randint(40, 1500)])
    Q = trim([random.choice(GT) for _ in range(random.randint(0, 30))])
    st, step, t, h, q, base = run('00000', 'E', '0' + '010' * n + Tstr(Q),
                                  lambda k, p, q: k == 'C>' and p == 3 * n + 3, 10 ** 8)
    ok = (st == 'C>' and cells(t, 0, h).lstrip('0') == Khub(1, n) and
          cells(t, h, len(t)).rstrip('0') == dword(1, Q) and step == 3 * n + 5)
    cnt['R'][0] += 1; cnt['R'][1] += ok
    if not ok: fails.append(('R', n, sym(Q)))
for _ in range(NR):   # B bounce
    k1 = random.choice([0, 1, 2, random.randint(0, 25)]); k2 = random.choice([0, 1, 2, random.randint(0, 25)])
    P = [random.choice(G3) for _ in range(k1)]
    if random.random() < 0.5:
        P = ['001' if g == '000' else g for g in P]
    G = trim(P + [random.choice(GT) for _ in range(k2)]); x0 = random.choice([0, 1, 1]); X = dword(x0, G)
    a = random.choice([1, 2, random.randint(1, 60), random.randint(60, 2000)])
    b = random.choice([1, 2, random.randint(1, 60), random.randint(60, 2000)])
    st, step, t, h, q, base = run(Khub(a, b + 1), 'C', X, lambda k, p, q: k == 'C>' and p == 0, 10 ** 9)
    y0, H = inc(x0, G)
    ok = (st == 'C>' and cells(t, 0, h).lstrip('0') == Khub(a + 1, b) and
          (cells(t, h, len(t)).rstrip('0') or '0') == dword(y0, H))
    cnt['B'][0] += 1; cnt['B'][1] += ok
    if not ok: fails.append(('B', a, b, X))
for _ in range(NR):   # E eat (normalized words)
    G = trim([random.choice(G3) for _ in range(random.choice([0, 1, 2, random.randint(0, 30)]))])
    if random.random() < 0.5:
        G = trim(['001' if (g == '000' and random.random() < 0.7) else g for g in G])
    x0 = random.choice([0, 1]); X = dword(x0, G)
    n = random.choice([1, 2, random.randint(1, 60), random.randint(60, 2000)])
    st, step, t, h, q, base = run(Khub(n, 1), 'C', X, lambda k, p, q: k == 'LREC' and q == 4, 10 ** 9)
    n2, Q2, _ = eat(n, x0, G)
    ok = st == 'LREC' and cells(t, h, len(t)).rstrip('0') == ('0' + '010' * n2 + Tstr(Q2)).rstrip('0')
    cnt['E'][0] += 1; cnt['E'][1] += ok
    if not ok: fails.append(('E', n, X))
report('2.2', not fails, 'rules R (restart, 3n+5 steps), B (bounce = inc), E (eat) on random instances vs TM',
       {k: '%d/%d' % (v[1], v[0]) for k, v in cnt.items()})

# exhaustive inc check on short words: every X in x0{000,001,010}*{100,010,000}* with <= 5 groups, a,b small
from itertools import product
tot = okc = 0
for L1 in range(0, 4):
    for L2 in range(0, 6 - L1):
        if L1 + L2 > 5: continue
        for P in product(G3, repeat=L1):
            for T in product(GT, repeat=L2):
                G = trim(list(P) + list(T))
                for x0 in (0, 1):
                    X = dword(x0, G)
                    st, step, t, h, q, base = run(Khub(2, 3), 'C', X, lambda k, p, q: k == 'C>' and p == 0, 10 ** 7)
                    y0, H = inc(x0, G)
                    ok = (st == 'C>' and cells(t, 0, h).lstrip('0') == Khub(3, 2) and
                          (cells(t, h, len(t)).rstrip('0') or '0') == dword(y0, H))
                    tot += 1; okc += ok
report('2.2x', okc == tot, 'bounce rule B exhaustive over all words with <=5 groups (a=2,b=3)', '%d/%d' % (okc, tot))

# ------------------------------------------------- 2.3x inc on all words of L_X with <= 7 groups, value lemma, carry reading
def carry_wrong(G):
    """WRONG reading of the S case (carry resumes after the toggled segment); used only as a mutation test."""
    G = list(G); i = 0
    while True:
        if i == len(G): G.append('001'); return G
        g = G[i]
        if g == '010': i += 1
        elif g == '001': G[i] = '000'; i += 1
        elif g == '000': G[i] = '001'; return G
        else:   # '100'
            G[i] = '010'; j = i + 1
            while True:
                if j == len(G): G.append('010'); break
                if G[j] == '010': G[j] = '000'
                elif G[j] == '000': G[j] = '010'; break
                j += 1
            i = j + 1


def value(x0, G):
    v = x0; k = 0
    for g in G:
        if g in ('000', '001'):
            k += 1
            if g == '001': v += 1 << k
    return v


def in_LX(G):
    s = ''.join({'000': '0', '001': '1', '010': 'i', '100': 'S'}[g] for g in G)
    return re.fullmatch(r'[01i]*[Si0]*', s) is not None


MAXG = 6 if QUICK else 7
seen = set(); tot = okc = wrong_dis = val_ok = clos_ok = 0
for L1 in range(0, MAXG + 1):
    for L2 in range(0, MAXG + 1 - L1):
        for P in product(G3, repeat=L1):
            for T in product(GT, repeat=L2):
                G = trim(list(P) + list(T))
                for x0 in (0, 1):
                    X = dword(x0, G)
                    if X in seen: continue
                    seen.add(X)
                    st, step, t, h, q, base = run(Khub(2, 3), 'C', X, lambda k, p, q: k == 'C>' and p == 0, 10 ** 7)
                    y0, H = inc(x0, G)
                    tmX = (cells(t, h, len(t)).rstrip('0') or '0') if st == 'C>' else None
                    ok = st == 'C>' and cells(t, 0, h).lstrip('0') == Khub(3, 2) and tmX == dword(y0, H)
                    tot += 1; okc += ok
                    if x0 == 1 and tmX != dword(0, trim(carry_wrong(G))): wrong_dis += 1
                    val_ok += value(y0, H) == value(x0, G) + 1
                    clos_ok += in_LX(H)
report('2.3x', okc == tot and val_ok == tot and clos_ok == tot and wrong_dis > 0,
       'inc vs TM on every word of L_X with <=%d groups; value(inc X)=value(X)+1; inc keeps L_X; the other reading of the 100 bullet is refuted' % MAXG,
       'TM %d/%d, value %d/%d, closure %d/%d; "resume after the toggled segment" disagrees with the TM on %d words' % (okc, tot, val_ok, tot, clos_ok, tot, wrong_dis))


# ------------------------------------------------- 2.5 L1 lemmas (proof.md 2.5), each checked as an exact window rule
def win(cl, pos, q, until=None, limit=10 ** 6):
    """Run the TM on the finite window cl (list of 0/1) from head pos in state q.
    Stops when the head leaves the window ('L'/'R') or until(h, q) is true ('STOP'). Returns (how, state, cells, steps, h)."""
    t = list(cl); h = pos; q = ord(q) - 65; s = 0
    while s < limit:
        if h < 0: return 'L', 'ABCDEF'[q], ''.join(map(str, t)), s, h
        if h >= len(t): return 'R', 'ABCDEF'[q], ''.join(map(str, t)), s, h
        if until and s > 0 and until(h, 'ABCDEF'[q]): return 'STOP', 'ABCDEF'[q], ''.join(map(str, t)), s, h
        tr = TR[(q, t[h])]
        if tr is None: return 'HALT', 'ABCDEF'[q], ''.join(map(str, t)), s, h
        t[h] = tr[0]; h += tr[1]; q = tr[2]; s += 1
    return 'MAX', 'ABCDEF'[q], ''.join(map(str, t)), s, h


def bits(w): return [int(x) for x in w]


def tog(tau): return [{'010': '000', '100': '100'}[g] for g in tau]


L1 = {}
# L0 additions used by the L1 lemmas: (window, head index, state) -> (exit side, state, window', steps)
L0X = [('1101', 3, 'E', ('L', 'E', '0100', 6)),     # S-restore: 1 101 <E -> <E 0 100 (E reads the 1 left of 101)
       ('1100', 2, 'E', ('L', 'E', '0100', 5)),     # S-restore tail used by TG: 1 1 0 <E 0 -> <E 0 1 0 0
       ('000', 1, 'E', ('R', 'A', '110', 4)),        # placeholder restart: 0 0 <E 0 -> 110 A>
       ('000', 2, 'E', ('R', 'F', '011', 3)),        # 000 <E -> 011 F>   (end of JW-010)
       ('010', 0, 'F', ('R', 'A', '000', 3)),        # F> 010 -> 000 A>   (start of JW-010)
       ('000', 0, 'F', ('R', 'C', '011', 3)),        # F> 000 -> 011 C>   (JW-000)
       ('001', 0, 'F', ('R', 'F', '011', 3)),        # F> 001 -> 011 F>   (JW-001)
       ('00110', 4, 'C', ('R', 'F', '11011', 9)),    # even entry: 0 011 C> 0 -> 110 11 F>
       ('00111', 4, 'E', ('R', 'C', '11011', 9))]    # odd exit:  0 011 1 <E (E on x0=1) -> 110 11 C>
ok0 = sum(win(bits(w), p, q)[:3] + (win(bits(w), p, q)[3],) == exp for (w, p, q, exp) in L0X)
L1['L0x'] = '%d/%d' % (ok0, len(L0X))
# TG toggle lemma: 1 B> tau 000  -->  <E 0 tog(tau) 010, tau in {010,100}^k (exhaustive k<=8); a sentinel 1 after the window must not be read
tg_bad = tg_tot = 0
for k in range(0, 9):
    for tau in product(['010', '100'], repeat=k):
        r = win(bits('1' + ''.join(tau) + '000' + '1'), 1, 'B')
        exp = '0' + ''.join(tog(tau)) + '010' + '1'
        tg_tot += 1; tg_bad += not (r[0] == 'L' and r[1] == 'E' and r[2] == exp)
L1['TG'] = '%d/%d' % (tg_tot - tg_bad, tg_tot)
# CA carry lemma: A> w 0^inf --> <E carry(w) 0^inf for w in {000,001,010}*{100,010,000}* (exhaustive <=6 groups + random <=30)
ca_bad = ca_tot = 0
def ca_check(w):
    r = win(bits(''.join(w) + '000' * (len(w) + 3)), 0, 'A')
    return r[0] == 'L' and r[1] == 'E' and r[2].rstrip('0') == ''.join(trim(carry(w)[0])).rstrip('0')
for L1n in range(0, 7):
    for L2n in range(0, 7 - L1n):
        for P in product(G3, repeat=L1n):
            for T in product(GT, repeat=L2n):
                ca_tot += 1; ca_bad += not ca_check(list(P) + list(T))
for _ in range(300 if QUICK else 3000):
    w = [random.choice(G3) for _ in range(random.randint(0, 15))] + [random.choice(GT) for _ in range(random.randint(0, 15))]
    ca_tot += 1; ca_bad += not ca_check(w)
L1['CA'] = '%d/%d' % (ca_tot - ca_bad, ca_tot)
# BL bounce-left: (110)^a 0 (011)^b 0 [E]1 0 Y --> (110)^(a+1) 0 (011)^b C> Y in 8b+6 steps (Y never read); split BLa (left pass) + BLb (right sweep)
bl_bad = 0
for b in range(1, 41):
    ra = win(bits('0' + '011' * b + '010'), 3 * b + 2, 'E', until=lambda h, q: h == 0 and q == 'A')
    rb = win(bits('01' + '010' * b + '00'), 0, 'A')
    bl_bad += not (ra[0] == 'STOP' and ra[2] == '01' + '010' * b + '00' and ra[3] == 5 * b + 2)
    bl_bad += not (rb[0] == 'R' and rb[1] == 'C' and rb[2] == '1100' + '011' * b and rb[3] == 3 * b + 4)
    for a in (1, 2, 7, 40):
        rf = win(bits('110' * a + '0' + '011' * b + '010'), 3 * a + 3 * b + 2, 'E')
        bl_bad += not (rf[0] == 'R' and rf[1] == 'C' and rf[2] == '110' * (a + 1) + '0' + '011' * b and rf[3] == 8 * b + 6)
L1['BL'] = '%d/%d' % (40 * 6 - bl_bad, 40 * 6)
# FT final turn: (110)^k 1 1 C> 0 --> <E (010)^k 0 0 1  in 5k+3 steps
ft_bad = 0
for k in range(1, 41):
    r = win(bits('110' * k + '11' + '0'), 3 * k + 2, 'C')
    ft_bad += not (r[0] == 'L' and r[1] == 'E' and r[2] == '010' * k + '001' and r[3] == 5 * k + 3)
L1['FT'] = '%d/40' % (40 - ft_bad)
# JW-010: F> 010 w 0^inf --> 011 F> carry(w) 0^inf (w normalized; head arrives at the start of carry(w))
jw_bad = jw_tot = 0
for _ in range(300 if QUICK else 3000):
    w = trim([random.choice(G3) for _ in range(random.randint(0, 25))])
    r = win(bits('010' + ''.join(w) + '000' * (len(w) + 3)), 0, 'F', until=lambda h, q: h == 3 and q == 'F')
    jw_tot += 1; jw_bad += not (r[0] == 'STOP' and r[2].rstrip('0') == ('011' + ''.join(trim(carry(w)[0]))).rstrip('0'))
L1['JW010'] = '%d/%d' % (jw_tot - jw_bad, jw_tot)
# E-odd / E-even intermediate configuration (110)^(a+1) 1 1 {C>,F>} on the full TM: K(a,1,1G) -> (110)^(a+1) 11 C> carry(G); K(a,1,0G) -> (110)^(a+1) 11 F> G
def first_at(left, right, pos_off, qwant, maxs=10 ** 7):
    pad = 8 * (len(left) + len(right)) + 64
    t = bytearray(pad) + bytearray(int(c) for c in left + right) + bytearray(pad)
    h = pad + len(left); q = 2; s = 0; target = pad + len(left) + pos_off
    while s < maxs:
        tr = TAB[2 * q + t[h]]
        if tr is None: return None
        t[h] = tr[0]; h += tr[1]; q = tr[2]; s += 1
        if h == target and q == qwant:
            return cells(t, 0, h).lstrip('0'), cells(t, h, len(t)).rstrip('0')
    return None
ee_bad = ee_tot = 0
for _ in range(200 if QUICK else 1000):
    a = random.choice([1, 2, random.randint(1, 50), random.randint(50, 500)])
    G = trim([random.choice(G3) for _ in range(random.randint(0, 25))])
    r1 = first_at(Khub(a, 1), dword(1, G) if G else '1', 1, 2)
    r0 = first_at(Khub(a, 1), '0' + ''.join(G), 1, 5)
    ee_tot += 2
    ee_bad += not (r1 and r1[0] == '110' * (a + 1) + '11' and r1[1] == ''.join(trim(carry(G)[0])).rstrip('0'))
    ee_bad += not (r0 and r0[0] == '110' * (a + 1) + '11' and r0[1] == ''.join(G).rstrip('0'))
L1['Eentry'] = '%d/%d' % (ee_tot - ee_bad, ee_tot)
report('2.5', ok0 == len(L0X) and tg_bad == 0 and ca_bad == 0 and bl_bad == 0 and ft_bad == 0 and jw_bad == 0 and ee_bad == 0,
       'L1 lemmas TG (toggle), CA (carry), BL (bounce left part), FT (final turn), JW-001/010/000, eat entries, as exact window runs', L1)

# ------------------------------------------------- 4. initial segment
st, step, t, h, q, base = run('', 'A', '0', lambda k, p, q: False, 27)
right = cells(t, h, len(t)).rstrip('0'); left = cells(t, 0, h).lstrip('0')
ok = st == 'MAX' and q == 4 and left == '' and t[h] == 0 and right == '0' + '010' * 2 + '001'
report('4', ok, 'blank tape ->(27 steps) L(2, empty) = 0^inf <E (010)^2 001 0^inf', {'state': 'ABCDEF'[q], 'right': right})
st2, step2, *_ = run('', 'A', '0', lambda k, p, q: False, 27)
report('4b', st2 == 'MAX', 'no halt in the first 27 steps', step2)


# ------------------------------------------------- 3. arithmetic (n, c) model
def Cw(c):
    return sum(x << j for j, x in enumerate(c))


def layout_of(G):
    c = [0]
    for g in G:
        if g == '010': c[-1] += 1
        else: c.append(0)
    while len(c) > 1 and c[-1] == 0: c.pop()
    return c


def step_nc(n, c):
    """(n, c) at an eat -> (n', c', D, M).  D = normalization depth of the next residue."""
    c = list(c)
    if n & 1:
        H = c; n2 = n + 1; M = (n + 1) // 2 + Cw(c)
    else:
        V = n // 2 + Cw(c); t = 0
        while (V >> t) & 1: t += 1
        H = c[t + 1:]; n2 = n + 2 + t + sum(c[:t + 1]); M = V >> (t + 1)
    Js = [j for j, x in enumerate(H) if x > 0]
    D = 1 + sum(1 for j in range(max(Js)) if not (M >> j) & 1) if Js else 0
    c2 = [0]
    for j in range(max(len(H), M.bit_length()) + 1):
        c2[-1] += H[j] if j < len(H) else 0
        if (M >> j) & 1: c2[-1] += 1
        else: c2.append(0)
    while len(c2) > 1 and c2[-1] == 0: c2.pop()
    return n2, c2, D, M


def B(x):
    return x.bit_length()


# 3a: (n,c) model == deepen epoch map (n, Q) incl. D, along the true orbit
E1 = 3000 if QUICK else 20000
n, Q = 2, []
x0, G, _ = bounce_phase(n, Q); c = layout_of(G)
mism = 0
for tt in range(E1):
    n2, Q2, info = epoch(n, Q)
    x0b, Gb, na = bounce_phase(n2, Q2)
    Dm = (na + 1).bit_length() - 1 if na is not None else None
    nn, cc, D, M = step_nc(n, c)
    if not (nn == n2 and cc == layout_of(Gb) and D == Dm): mism += 1
    n, Q, c = n2, Q2, cc
report('3a', mism == 0, '(n,c) closed-form map == deepen epoch map (n, layout, D) along the orbit',
       '%d/%d epochs, final n=%d' % (E1 - mism, E1, n))

# 3b: closed form == epoch map on random reachable-shape states (Q = image of a random normalized counter)
NS = 2000 if QUICK else 20000
mism = tot = 0
for _ in range(NS):
    L = random.randint(0, 30)
    H0 = [random.choice(G3) for _ in range(L)] + ['001']
    Q = regroup_tail(H0)
    Bn = random.randint(2, 30); n = random.randint(max(2, 1 << (Bn - 1)), (1 << Bn) - 1)
    x0, G, na = bounce_phase(n, Q)
    if na is None: continue     # not normalized at this eat: outside the domain
    c = layout_of(G)
    n2, Q2, info = epoch(n, Q)
    x0b, Gb, nb = bounce_phase(n2, Q2)
    Dm = (nb + 1).bit_length() - 1 if nb is not None else None
    nn, cc, D, M = step_nc(n, c)
    tot += 1
    # D from the closed form is the depth even when it exceeds B(n')-1 (then bounce_phase reports None)
    ok = nn == n2 and (cc == layout_of(Gb) if nb is not None else True) and (D == Dm if Dm is not None else D > B(n2) - 1)
    mism += not ok
report('3b', mism == 0, 'closed form (M, D, n\', c\') == epoch map on random normalized states', '%d/%d' % (tot - mism, tot))

# 3c: TM left records (1e9 steps) == deepen (n, T) map
recs = []
for line in open(os.path.join(HERE, 'data', 'left_records_1e9.txt')):
    p = line.split()
    if len(p) >= 4 and int(p[1]) >= 2: recs.append((int(p[0]), int(p[1]), p[3]))
okr = 0
for (s1, n1, T1), (s2, n2, T2) in zip(recs, recs[1:]):
    Q = parse('1' + T1[3:])[1] if len(T1) > 3 else []
    nn, QQ, info = epoch(n1, Q)
    okr += (nn == n2 and QQ is not None and Tstr(QQ) == T2)
report('3c', okr == len(recs) - 1, 'epoch map vs TM left records to 1e9 steps (file data/left_records_1e9.txt)',
       '%d/%d epochs, n %d -> %d' % (okr, len(recs) - 1, recs[0][1], recs[-1][1]))

# ------------------------------------------------- 5. what is proved about the map (proof.md 5.2-5.5), sanity-checked
bad1 = bad2 = n1c = n2c = 0
for _ in range(NS):
    Bn = random.randint(3, 40); n = random.randint(1 << (Bn - 1), (1 << Bn) - 1)
    L = random.randint(1, Bn + 2)
    c = [random.choice([0, 0, 1, 1, 2, 3, random.randint(0, 40)]) for _ in range(L)]
    n2, c2, D, M = step_nc(n, c)
    J = max([j for j, x in enumerate(c) if x > 0], default=-1)
    if n % 2 == 0 and J <= B(n) - 1:
        n1c += 1; bad1 += D > B(n2) - 1
    if n % 2 == 1 and Cw(c) < (1 << (B(n + 1) - 1)):
        n2c += 1; bad2 += D > B(n2) - 1
report('5.2', bad1 == 0 and bad2 == 0, 'Lemma S-even (J(c)<=B(n)-1) and Lemma S-odd (C<2^(B(n+1)-1)) imply D<=B(n\')-1',
       'even %d/%d, odd %d/%d random states' % (n1c - bad1, n1c, n2c - bad2, n2c))
# 5.3 the simple bounds are NOT inductive: explicit witness
Bw = 12; nw = (1 << (Bw - 1)) + 1; cw = [(1 << (Bw - 2)) - 1]
s1 = step_nc(nw, cw); s2 = step_nc(s1[0], s1[1])
wit = {'n': nw, 'C': Cw(cw), 'step1': (s1[0], Cw(s1[1]), 'J', len(s1[1]) - 1), 'step2': (s2[0], Cw(s2[1])),
       'C2/2^(B-1)': Cw(s2[1]) / 2 ** (B(s2[0]) - 1), 'D_each_step': (s1[2], s2[2])}
report('5.3', Cw(s2[1]) >= (1 << (B(s2[0]) - 1)) and s2[0] % 2 == 1,
       'witness: Inv_S (odd: C<2^(B-1); even: J<=B-1) is not closed under the map', wit)


def Pot(n, c):
    return (n + 1) // 2 + Cw(c)


def Zb(M):
    return M.bit_length() - bin(M).count('1') if M > 0 else -1


def Jc(c):
    return max([j for j, x in enumerate(c) if x > 0], default=-1)


def rand_c(Bn):
    L = random.randint(1, max(1, Bn - 2)); k = random.random()
    if k < 0.3: return [random.choice([0, 1]) for _ in range(L)]
    if k < 0.6: return [random.choice([0, 0, 1, 2, 3]) for _ in range(L)]
    return [random.randint(0, 6) if j < 3 else random.choice([0, 1]) for j in range(L)]


def cascade_start(Bn, c):
    """odd n with M=(n+1)/2+C(c)=2^(B-1), B=B(n+1); shrinks c until C<=2^(B-2)-1. Returns (n, c) or None."""
    c = list(c)
    while Cw(c) > (1 << (Bn - 2)) - 1 and any(c):
        k = max(j for j, x in enumerate(c) if x > 0); c[k] -= 1
    n = (1 << Bn) - 2 * Cw(c) - 1
    return (n, c) if n > 0 and n % 2 == 1 and B(n + 1) == Bn else None


# 5.4 cascade lemma: (i) V'=2^B, M'=2^(B-1), D'=B-1, slack' = 0 if n''<2^B else 1;
#     (ii) c0 odd and n''+1<2^B: n'' odd, C''=3*2^(B-2)+(C-c0)/2 >= 2^(B(n''+1)-1) (Inv_S fails), 2M''=5*2^(B-1)+3-C, D''<=B-1
#     (iii) any c0: D''<=B-1, so the eat two epochs later is normalized
def cascade_check(n, c):
    Bn = B(n + 1); C = Cw(c)
    n1, c1, D0, M0 = step_nc(n, c)
    n2, c2, D1, M1 = step_nc(n1, c1)
    n3, c3, D2, M2 = step_nc(n2, c2)
    ok = (M0 == 1 << (Bn - 1) and n1 == n + 1 and n1 // 2 + Cw(c1) == 1 << Bn and M1 == 1 << (Bn - 1) and D1 == Bn - 1
          and B(n2) - 1 - D1 == (0 if n2 < (1 << Bn) else 1) and D2 <= Bn - 1 and B(n3) - 1 - D2 >= 0)
    odd = c[0] % 2 == 1 and n2 + 1 < (1 << Bn)
    if odd:
        ok = ok and n2 % 2 == 1 and Cw(c2) == 3 * (1 << (Bn - 2)) + (C - c[0]) // 2 \
             and Cw(c2) >= 1 << (B(n2 + 1) - 1) and 2 * M2 == 5 * (1 << (Bn - 1)) + 3 - C
    return ok, odd


casc = bad = nodd = sl1 = 0
for _ in range(NS):
    Bn = random.randint(4, 30)
    st = cascade_start(Bn, rand_c(Bn))
    if st is None: continue
    ok, odd = cascade_check(*st)
    casc += 1; bad += not ok; nodd += odd
report('5.4', bad == 0, 'cascade lemma (i)-(iii): slack exactly 0 iff n\'\'<2^B; with c0 odd Inv_S fails two epochs later; that eat is still normalized',
       '%d/%d constructed cascades (%d with c0 odd and n\'\'+1<2^B)' % (casc - bad, casc, nodd))

# 5.5 Lemma Z (J(c') = max(D-1, Z(M))) and Lemma R (slack_t>=0 and P_t<2^B(n_{t+1}) => slack_{t+1}>=0) on arbitrary states
zb = zt = rb = rt = 0
for _ in range(NS * 2):
    Bn = random.randint(3, 36); n = random.randint(1 << (Bn - 1), (1 << Bn) - 1)
    c = [random.choice([0, 0, 1, 1, 2, 3, random.randint(0, 9), random.randint(0, 200)]) for _ in range(random.randint(1, Bn + 2))]
    if random.random() < 0.3:     # bias towards sparse M (near cascades)
        st = cascade_start(Bn, c) if random.random() < 0.5 else None
        if st: n, c = st
    n1, c1, D0, M0 = step_nc(n, c)
    Hinert = D0 > 0
    zt += 1; zb += Jc(c1) != max(D0 - 1 if Hinert else -1, Zb(M0))
    if B(n1) - 1 - D0 >= 0 and Pot(n, c) < (1 << B(n1)):
        n2, c2, D1, M1 = step_nc(n1, c1)
        rt += 1; rb += B(n2) - 1 - D1 < 0
report('5.5', zb == 0 and rb == 0, 'Lemma Z: J(c\')=max(D-1, Z(M)); Lemma R: slack_t>=0 and P_t<2^B(n_{t+1}) imply slack_{t+1}>=0',
       'Z %d/%d, R %d/%d arbitrary states with the hypotheses' % (zt - zb, zt, rt - rb, rt))

# 5.6 evidence (not a proof): cascade and near-cascade aftermaths, 40 epochs each
NC = 3000 if QUICK else 30000
viol = tot = 0; exlen = {}
for i in range(NC):
    Bn = random.randint(6, 24)
    if i % 2 == 0:
        st = cascade_start(Bn, rand_c(Bn))
    else:   # near cascade: odd n, M = 2^(B-1)+2^k, inert groups only at slots <= k
        k = random.randint(0, Bn - 3); c = [random.choice([0, 0, 1, 2]) for _ in range(k + 1)]
        while Cw(c) > (1 << (Bn - 2)) and any(c):
            kk = max(j for j, x in enumerate(c) if x > 0); c[kk] -= 1
        n = 2 * ((1 << (Bn - 1)) + (1 << k) - Cw(c)) - 1
        st = (n, c) if n > 0 and B(n + 1) == Bn else None
    if st is None: continue
    n, c = st; tot += 1; run_ = best = 0
    for t in range(40):
        n2, c2, D, M = step_nc(n, c)
        if B(n2) - 1 - D < 0: viol += 1; break
        if Pot(n, c) >= (1 << B(n2)): run_ += 1; best = max(best, run_)
        else: run_ = 0
        n, c = n2, c2
    exlen[best] = exlen.get(best, 0) + 1
report('5.6', viol == 0, 'evidence: constructed cascade / near-cascade starts (B=6..24, small c), 40 epochs each: no non-normalized eat',
       {'starts': tot, 'violations': viol, 'longest excursion per start (epochs with P>=2^B(n\')) -> count': dict(sorted(exlen.items()))})

# ------------------------------------------------- 6. the orbit
E2 = 60000 if QUICK else 400000
n, Q = 2, []
x0, G, _ = bounce_phase(n, Q); c = layout_of(G)
minslack = 99; tight = []; invS_fail = []; exc = []; maxCodd = 0.0; t_maxCodd = None; notR = []
prev_exc = True
t0 = time.time()
for tt in range(E2):
    n2, c2, D, M = step_nc(n, c)
    sl = B(n2) - 1 - D
    minslack = min(minslack, sl)
    if sl <= 1 and tt >= 16: tight.append((tt, n, sl))
    if tt >= 20:
        J = Jc(c)
        if n & 1:
            r = Cw(c) / 2 ** (B(n) - 1)
            if r > maxCodd: maxCodd, t_maxCodd = r, tt
            if not Cw(c) < (1 << (B(n + 1) - 1)): invS_fail.append(tt)
        elif not J <= B(n) - 1: invS_fail.append(tt)
    if tt >= 1 and prev_exc: notR.append((tt, sl))       # epochs not covered by Lemma R (previous epoch was an excursion)
    ex = Pot(n, c) >= (1 << B(n2))
    if ex and tt >= 20: exc.append(tt)
    prev_exc = ex
    n, c = n2, c2
late_notR = [x for x in notR if x[0] >= 20]
report('6a', minslack >= 0, 'orbit: counter normalized before every eat (slack B(n_{t+1})-1-D_t >= 0)',
       'epochs %d, final n=%d, min slack %d; t>=16 with slack<=1: %s; %.1fs' % (E2, n, minslack, tight[:6], time.time() - t0))
report('6b', invS_fail == ([] if E2 <= 350274 else [350274]), 'orbit t>=20: Inv_S fails exactly at t=350,274 (so it is NOT an invariant of the orbit); odd-n peak of C/2^(B-1)',
       'Inv_S fails at t=%s; max C/2^(B-1) at odd n = %.4f (t=%s)' % (invS_fail[:5], maxCodd, t_maxCodd))
report('6c', all(s >= 0 for _, s in notR), 'orbit: excursions (P_t>=2^B(n_{t+1})) and the epochs Lemma R does not cover, checked directly',
       'excursion epochs t>=20: %s; epochs after an excursion (t, slack): %s' % (exc[:6], late_notR[:6]))
RESULTS['6a']['tight_epochs(t,n,slack)'] = tight[:20]

# ------------------------------------------------- 7. other orbits of the map (small random starts)
RUNS = 200 if QUICK else 2000; EP = 1500 if QUICK else 3000
viol = sl0 = sl0_casc = sl1 = steps = invf = 0; maxCr = maxP = 0.0; exn = exent_c = exent_o = exmax = 0
for r in range(RUNS):
    n = random.randint(2, 1 << random.randint(2, 16))
    c = [random.choice([0, 0, 0, 1, 2]) for _ in range(random.randint(0, max(1, B(n) - 3)))] or [0]
    while Cw(c) >= (1 << (B(n) - 1)) and len(c) > 1: c.pop()
    pcasc = False; pex = False; run_ = 0
    for tt in range(EP):
        n2, c2, D, M = step_nc(n, c)
        sl = B(n2) - 1 - D; steps += 1
        if sl < 0: viol += 1; break
        casc_now = (n & 1) == 1 and M == 1 << (B(n2) - 1)
        ex = Pot(n, c) >= (1 << B(n2))
        if tt >= 50:
            if sl == 0: sl0 += 1; sl0_casc += pcasc
            if sl == 1: sl1 += 1
            ok_inv = (Cw(c) < (1 << (B(n + 1) - 1))) if n & 1 else (Jc(c) <= B(n) - 1)
            invf += not ok_inv
            maxCr = max(maxCr, Cw(c) / 2 ** (B(n) - 1)); maxP = max(maxP, Pot(n, c) / 2 ** B(n))
            if ex:
                exn += 1; run_ += 1; exmax = max(exmax, run_)
                if not pex:
                    if pcasc: exent_c += 1
                    else: exent_o += 1
            else: run_ = 0
        pcasc, pex = casc_now, ex
        n, c = n2, c2
report('7', viol == 0, 'other orbits of the map (small random starts): no non-normalized eat',
       {'orbits x epochs': '%d x %d (%d eats)' % (RUNS, EP, steps), 'slack-0 eats after burn-in': sl0,
        'of which right after a cascade': sl0_casc, 'slack-1 eats': sl1, 'Inv_S failures': invf,
        'max C/2^(B-1)': round(maxCr, 4), 'max P/2^B(n)': round(maxP, 4), 'excursion epochs': exn,
        'excursion entries (after cascade / other)': (exent_c, exent_o), 'longest excursion': exmax})


# ------------------------------------------------- 8. arbitrary (non-reachable) starts: the map alone does not stay safe
def report_info(cid, what, info):
    RESULTS[cid] = {'pass': True, 'what': what, 'info': info, 'informational': True}
    print('[INFO] %s %s : %s' % (cid, what, info), flush=True)


RUNS8 = 2000 if QUICK else 20000; EP8 = 60 if QUICK else 200
when = {}
for r in range(RUNS8):
    Bn = random.randint(4, 20); n = random.randint(1 << (Bn - 1), (1 << Bn) - 1)
    c = [random.choice([0, 0, 1, 1, 2, 3, random.randint(0, 50)]) for _ in range(random.randint(1, Bn + 1))]
    for tt in range(EP8):
        n2, c2, D, M = step_nc(n, c)
        if B(n2) - 1 - D < 0:
            when[tt] = when.get(tt, 0) + 1; break
        n, c = n2, c2
report_info('8', 'arbitrary starts (any c, entries up to 50): first non-normalized eat by epoch (larger searches find epochs 5-7 too)',
            {'runs': RUNS8, 'epochs': EP8, 'violations_by_epoch': dict(sorted(when.items()))})


def first_violation(n, c, K):
    for tt in range(K):
        n2, c2, D, M = step_nc(n, c)
        if B(n2) - 1 - D < 0: return tt, n
        n, c = n2, c2
    return None


ex1 = first_violation(20, [0, 0, 0, 22], 20); ex2 = first_violation(16384, [0] * 13 + [1295], 5)
# super-cascade: odd n, C >= 2^(B-1), M = 2^B exactly, n'' < 2^B  ->  slack -1 one epoch later
Bs = 10; ns = 1019; cs = [0, 1] + [0] * 6 + [2]      # (n+1)/2 = 510, C = 2 + 2*256 = 514, M = 1024; slack 0 at epoch 0
ex3 = first_violation(ns, cs, 5)
report('8b', ex1 == (5, 56) and ex2 == (1, 16386) and ex3 is not None and ex3[0] == 1,
       'explicit violating starts of the map: (20,[0,0,0,22]) at epoch 5; (16384,[0]*13+[1295]) at epoch 1; super-cascade (1019,[0,1,0,0,0,0,0,0,2]) at epoch 1',
       {'ex1': ex1, 'ex2': ex2, 'super-cascade': ex3, 'super-cascade C/2^(B-1)': Cw(cs) / 2 ** (Bs - 1)})


# 8c the TM itself reaches a non-normalized eat from K(20, 1, X(20, [0,0,0,22])) at eat 5
def eat_word(n, c):
    m = n >> 1; G = []
    for j in range(max(len(c), B(n) - 1)):
        G += ['010'] * (c[j] if j < len(c) else 0)
        G.append('001' if (m >> j) & 1 else '000')
    return n & 1, trim(G)


x0s, Gs = eat_word(20, [0, 0, 0, 22])
pad = 4000
tp = bytearray(pad) + bytearray(int(ch) for ch in Khub(20, 1) + dword(x0s, Gs)) + bytearray(pad)
h = pad + len(Khub(20, 1)); lo = pad; q = 2; s = 0; eats = []; want = None; HUB = bytes([1, 1, 0, 0, 0, 1, 1])
while s < 2 * 10 ** 6 and len(eats) < 8:
    tr = TAB[2 * q + tp[h]]
    if tr is None: eats.append('HALT'); break
    tp[h] = tr[0]; h += tr[1]; q = tr[2]; s += 1
    if h < lo:
        lo = h
        if q == 4:          # left record L(n, Q): the eat of this epoch is the first K(n, 1, X) after it
            Rs = cells(tp, h + 1, len(tp)); want = 0
            while Rs[3 * want:3 * want + 3] == '010': want += 1
    elif q == 2 and tr[1] == 1 and want and bytes(tp[h - 7:h]) == HUB and cells(tp, 0, h).lstrip('0') == Khub(want, 1):
        Xs = cells(tp, h, len(tp)).rstrip('0'); x0c, Gc = parse(Xs) if len(Xs) > 1 else (int(Xs or '0'), [])
        eats.append((s, want, normalized(Gc))); want = None
        if not normalized(Gc): break
report('8c', len(eats) == ex1[0] + 1 and all(e[2] for e in eats[:-1]) and eats[-1][2] is False,
       'TM from K(20,1,X(20,[0,0,0,22])): the first non-normalized eat is the one the map predicts (slack<0 at map epoch 5 = TM eat 6)',
       {'eats after each left record (step, n, normalized)': eats})

allok = all(v['pass'] for v in RESULTS.values() if not v.get('informational'))
print(json.dumps({'all_pass': allok, 'quick': QUICK, 'seed': SEED,
                  'results': {k: {'pass': v['pass'], 'info': v['info']} for k, v in RESULTS.items()}}, default=str))
sys.exit(0 if allok else 1)