"""dualquad.py -- where do the two ends of a quadratic base part?

QUESTION. beta a root of x^2 + px + q, -1 <= p <= q, q >= 2, digits 0,
..., n - 1, n = q = |N(beta)|, which writes every element of O = Z[beta]
once (Katai and Kovacs for irreducible f; the run checks T_J distinct
and T_j complete mod beta^j O). In the coordinates (x, y) = x + y beta,
beta acts by (x, y) -> (-qy, x - py) and beta-bar = -p - beta. The
span-j leading tile is a translate of T_j = {sum_(k<j) d_k beta^k}, the
depth-j trailing tile a coset of beta^j O. Depth j reads an ideal I iff
beta^j is in I; span j hides it iff T_j equidistributes mod I. At base
2i they part one way, a digit early (dual2i_split.py). Where do they
part at the other quadratic bases, and which way?

THE ARGUMENT (written before the engine).
  (1) THE PRODUCT. For a character chi of O/I the sum of e(chi) over
      T_j is prod_(k<j) sum_(d<n) e(d chi(beta^k)), and factor k is 0
      iff n chi(beta^k) is in Z and chi(beta^k) is not: span j hides I
      iff every nonzero chi has such a k < j.
  (2) T_j is a complete system mod beta^j O, so read implies hidden;
      hidden needs [O : I] | n^j.
  (3) THE LAW L, beta of prime norm in a maximal order: span j >= 1
      hides I iff beta^j or beta^(j-1) beta-bar lies in I. Ramified,
      beta-bar O = beta O and L is one condition. Split, the ideals of
      p-power index are P^a Pbar^b, O/I = Z/p^a x Z/p^b with beta ->
      (r, u), v_p(r) = 1, u a unit: (0, 1/p^b) never dies at b >= 2;
      at b = 1 every character dies by k = a, and at j = a the one
      character (1, c2) cancelling at k = a - 1 survives.
  (4) TRANSPLANTS, printed: L at the composite norm of -3 + i, and at
      the non-maximal Z[2i], where beta-bar = -beta makes L the read.

PREDICTIONS (frozen before the engine).
  C1 T_J has n^J distinct points and T_j is a complete system mod
     beta^j O, every base, j <= J, n^J <= 10^4.
  C2 2i: 23 ideals of 2-power index <= 256, gaps {0: 13, 1: 10},
     (4, 2, 1) span 1 depth 2, as dual2i.py printed.
  P0 the character criterion agrees with the box count at every base,
     ideal and j; plant (beta-bar^k for beta^k) > 0 off.
  P1 read at j implies hidden at j, 0 off.
  P2 L agrees with the count at -1 + i, -2 + i, -4 + i, -1 + sqrt(-2)
     and (-1 + sqrt(-7))/2, every ideal of index <= 300 or dividing
     n^J, every j <= J; plant (beta-bar squared) > 0 off.
  P3 at each split base (beta-bar) is hidden from span 1, read at no
     j <= J.
  P4 2i: L is the read at every ideal, so it misses exactly the split
     ideals.
  P5 -3 + i: L's disagreements printed.

KILL (frozen). L dies if P2 prints one disagreement; the criterion (1)
dies if P0 does. C1 and C2 are read first; a red control voids the run.

THE SWEEP (frozen after the first run, before its run). L held at
-3 + i as well. The if half of L holds in every Z[beta]: beta = -p mod
beta-bar, so O/(beta-bar) = Z/n is generated by 1, D is complete there,
and T_j is complete mod beta^(j-1) beta-bar. Span 1 hides I iff
O/I = Z/m, m | q, generated by 1, I = (m, beta - r) with
r(r + p) = 0 mod m, and L there is r = 0 or -p: a third root needs a
prime l | p with l^2 | q, Dedekind's test for Z[beta] non-maximal at l.
  P6 L agrees with the count at every base x^2 + px + q, q = 2..12,
     -1 <= p <= q, with no prime l | p and l^2 | q; the others' offs
     printed, 2i among them. C1 at every base.
KILL (frozen). P6 dies if one such base prints L off > 0.

L' (frozen after the sweep, before its run). The sweep's misses mix
the roots 0 and -p prime by prime of q. J_r = (q, beta - r) for
r(r + p) = 0 mod q has O/J_r = Z/q generated by 1, J_0 = (beta),
J_(-p) = (beta-bar), and T_j is complete mod beta^(j-1) J_r: the if
half of span j >= 1 hides I iff I contains some beta^(j-1) J_r; span
1's only if is the root count.
  P7 L' agrees with the count at every base of the sweep, 2i and the
     non-maximal ones included.
KILL (frozen). L' dies if one base prints L' off > 0.

THE COLON LAW (frozen after L' fired, before its run). L' missed at
q = 8, 4 | p: beta = -2 mod 4 roots r(r + p) mod 4 and lifts to no root
mod 8. Span 1 hides J iff O/J = Z/m, m | q, generated by 1, an
upward-closed family; so span j hides I iff I contains beta^(j-1) J
for such a J, iff span 1 hides (I : beta^(j-1)), the if half proved.
  P8 the colon law agrees with the count at every base of the sweep
     and the seven named.
KILL (frozen). The colon law dies if one base prints off > 0.

THE ONLY IF (proved by hand after P8, its two checks frozen before
their run). The step: at an ideal L of l-power index, v = v_l(q), a
character psi' = beta psi != 0 outside B (B the characters with
psi(1) != 0 and q psi(1) = 0) has a lift outside B; otherwise
beta-bar psi = (x, 0) with l^v x != 0 puts l^v characters (y, 0) in the
kernel of beta, so that kernel is (1/l^v)Z x 0 and holds a lift with
psi(1) = 0. Lifting from (I_l : beta^(j-1)) down to each l-primary part
I_l, and the CRT, give the only if; the gap follows: with s the least
span, O/(I : beta^(s-1)) = Z/m, beta -> r, the least depth is s - 1 + e,
e the least e with r^e = 0 mod m, and never if no power of r is 0.
  S the step holds at every ideal of one prime's power dividing n^J,
     0 dead ends; at mixed indices it dead-ends, only at q of two
     primes, so the proof needs the CRT.
  G the least depth to J is s - 1 + e at every ideal with s >= 1, 0
     off; plant (s + e) > 0 off.
KILL (frozen). The proof is void if S prints one dead end at one
prime; the gap law dies if G prints one off.

FINDINGS (entered after the run, from its printed output).
  C  green: T_J distinct and T_j complete mod beta^j O at all 106
     runs, 100 distinct bases; 2i's 23 ideals, gaps {0: 13, 1: 10},
     (4, 2, 1) span 1 depth 2.
  P0 the character criterion is the count at every base, 0 off; the
     plant (beta-bar^k) 142 off. P1 read implies hidden, 0 off.
  P2 L 0 off at the five prime-norm bases; plant 156 off. P3 (beta-bar)
     is (2, 0, 1) at -1 + i, ramified, span 1 depth 1, and at the four
     split bases span 1, read at no depth: (5, 4, 1), (17, 8, 1),
     (3, 2, 1), (2, 1, 1). -2 + i: 6 ideals one condition, 5 hidden and
     unread to J = 5.
  P4 2i: L is the read; it misses 16 (ideal, j) pairs, one per split
     ideal of the family. P5 -3 + i: L 0 off (2 | p leaves one prime
     of q with a choice).
  P6 FIRES: L misses at 13 of the 80 sweep bases with no l | p,
     l^2 | q (and at 17 of the 19 others), at most 30 off (x^2 - x + 6,
     disc -23, maximal): the CRT mixes the roots 0 and -p prime by prime
     of q.
  P7 FIRES: L' misses at x^2 + 8, x^2 + 4x + 8, x^2 + 8x + 8 only, 4
     off each, (4, 2, 1) at span 1 and (16, 8, 2) at span 2 first:
     beta = -2 mod 4 roots r(r + p) mod 4 and lifts to no root mod 8.
  P8 the colon law HOLDS: span j hides I iff span 1 hides
     (I : beta^(j-1)), at 100 distinct bases, 17144 (ideal, j) pairs,
     1 <= j <= J, n^J <= 10^4, whose index divides n^j, 9085 of them
     hidden, 0 off; the plant (beta^(j-2)) 3505 off. Both halves are
     proved, the only if by the step S checks.
  S 0 dead ends at one prime's power; 1574 at mixed indices, at
     q = 6, 10, 12 alone.
  G the gap law 0 off, plant 1844 off; gaps printed {0: 1615, 1: 213,
     2: 16} over the 100 bases.
  K  L dies; L' dies; the colon law stands, proved. P6 and P7 now
     assert the misses, which the colon law's roots mod every m | q
     repair. Two frozen sentences are false: the sweep's "a third root
     needs a prime l | p with l^2 | q" (x^2 - x + 6 has four) and L''s
     "span 1's only if is the root count" (roots mod q alone).

RUN RECORD. First run, the named bases: 15/15, 0.9 s, 12.1 MB. With
the sweep: 115/116, P6 FAIL, 14.9 s. With L': 115/117, P6 and P7 FAIL.
With the colon law: 116/118. P6 and P7 asserting the misses: 118/118,
15.1 s, peak commit 22.8 MB under a memory guard. With P8's plant in
the script and the distinct-base count: 119/119, 15.6 s, 23.2 MB. With
S and G: 123/123, 19.6 s, 26.5 MB. A review found S's second clause,
"only at q of two primes", true by the index alone (a mixed index
dividing n^J needs two primes of q); cut, the check now asks for a dead
end at some mixed ideal: 123/123, 19.6 s, 27.1 MB. A scratch copy with m | q^2 in
place of m | q passes P8, since an ideal containing beta^(j-1) m and
beta^(j-1)(beta - r) contains beta^(j-1) (gcd(m, q), beta - r).
"""

import sys

CHECKS = []
NMAX = 300
POINTS = 10 ** 4

# name, p, q
BASES = [
    ("-1 + i", 2, 2),
    ("-2 + i", 4, 5),
    ("-4 + i", 8, 17),
    ("-1 + sqrt(-2)", 2, 3),
    ("(-1 + sqrt(-7))/2", 1, 2),
    ("2i", 0, 4),
    ("-3 + i", 6, 10),
]
PRIME_NORM = {"-1 + i", "-2 + i", "-4 + i", "-1 + sqrt(-2)",
              "(-1 + sqrt(-7))/2"}


def check(name, ok, detail=""):
    CHECKS.append(bool(ok))
    print(f"  {'PASS' if ok else 'FAIL'}  {name}" +
          (f"  ({detail})" if detail else ""))


def section(title):
    print()
    print(title)


def mul(u, v, p, q):
    """(x1 + y1 beta)(x2 + y2 beta), beta^2 = -p beta - q."""
    x1, y1 = u
    x2, y2 = v
    return x1 * x2 - q * y1 * y2, x1 * y2 + x2 * y1 - p * y1 * y2


def powers(g, p, q, kmax):
    out = [(1, 0)]
    for _ in range(kmax):
        out.append(mul(out[-1], g, p, q))
    return out


class Lat:
    """Z(a, 0) + Z(c, d), 0 <= c < a: every full-rank sublattice once."""

    def __init__(self, a, c, d):
        self.a, self.c, self.d = a, c, d
        self.n = a * d

    def __contains__(self, v):
        x, y = v
        if y % self.d:
            return False
        return (x - (y // self.d) * self.c) % self.a == 0

    def key(self, v):
        x, y = v
        qq, r = divmod(y, self.d)
        return (x - qq * self.c) % self.a, r

    def characters(self):
        """chi = (A, G)/N, chi(x, y) = (Ax + Gy)/N mod 1: the N
        characters trivial on the lattice."""
        a, c, d, n = self.a, self.c, self.d, self.n
        return [(d * k1, (-c * k1 + a * k2) % n)
                for k1 in range(a) for k2 in range(d)]


def hnf(v1, v2):
    """The lattice spanned by v1, v2 as Lat(a, c, d)."""
    (x1, y1), (x2, y2) = v1, v2
    g, s, t = egcd(y1, y2)
    cx = s * x1 + t * x2
    ax = (y2 // g) * x1 - (y1 // g) * x2
    a, d = abs(ax), abs(g)
    return Lat(a, cx % a, d)


def egcd(u, v):
    if v == 0:
        return (u, 1, 0) if u >= 0 else (-u, -1, 0)
    g, s, t = egcd(v, u % v)
    return g, t, s - (u // v) * t


def divisors(m):
    return [k for k in range(1, m + 1) if m % k == 0]


def ideals(p, q, indices):
    out = []
    for nn in sorted(indices):
        for a in divisors(nn):
            d = nn // a
            for c in range(a):
                L = Lat(a, c, d)
                if (0, a) in L and (-q * d, c - p * d) in L:
                    out.append(L)
    return out


def tiles(p, q, n, jmax):
    pw = powers((0, 1), p, q, jmax)
    T = [[(0, 0)]]
    for j in range(1, jmax + 1):
        bx, by = pw[j - 1]
        T.append([(x + dd * bx, y + dd * by)
                  for (x, y) in T[-1] for dd in range(n)])
    return T


def hidden_by_count(L, Tj):
    if len(Tj) % L.n:
        return False
    cnt = {}
    for v in Tj:
        k = L.key(v)
        cnt[k] = cnt.get(k, 0) + 1
    return len(cnt) == L.n and len(set(cnt.values())) == 1


def least_hidden_by_characters(L, pw, n, jmax):
    """The least j <= jmax at which every nonzero character has a
    vanishing factor, or None."""
    N = L.n
    worst = 0
    for A, G in L.characters():
        if A == 0 and G == 0:
            continue
        first = None
        for k in range(jmax):
            x, y = pw[k]
            v = (A * x + G * y) % N
            if v and (n * v) % N == 0:
                first = k
                break
        if first is None:
            return None
        worst = max(worst, first + 1)
    return worst


def run_base(name, p, q):
    n = q
    J = 0
    while n ** (J + 1) <= POINTS:
        J += 1
    pw = powers((0, 1), p, q, J + 1)
    bar = (-p, -1)
    pwbar = powers(bar, p, q, J + 1)
    inds = set(range(1, NMAX + 1)) | set(divisors(n ** J))
    Is = ideals(p, q, inds)
    T = tiles(p, q, n, J)
    section(f"base {name}: x^2 + {p}x + {q}, n = {n}, J = {J}, "
            f"{len(Is)} ideals")

    # C1
    pts = set(T[J])
    okc = len(pts) == n ** J
    for j in range(J + 1):
        B = hnf(pw[j], mul(pw[j], (0, 1), p, q))
        okc &= B.n == n ** j and len({B.key(v) for v in T[j]}) == n ** j
    check(f"C1 {name}: T_J distinct, T_j complete mod beta^j O", okc,
          f"{len(pts)} points")

    rows = []
    off0 = off1 = offL = plantL = offR = offC = plantC = live = hid = 0
    mroots = [(m, r) for m in divisors(q) for r in range(m)
              if (r * (r + p)) % m == 0]
    roots = [r for r in range(q) if (r * (r + p)) % q == 0]
    badR = []
    for L in Is:
        cnt = [hidden_by_count(L, T[j]) for j in range(J + 1)]
        read = [pw[j] in L for j in range(J + 1)]
        s = next((j for j in range(J + 1) if cnt[j]), None)
        t = next((j for j in range(J + 1) if read[j]), None)
        if least_hidden_by_characters(L, pw, n, J) != s:
            off0 += 1
        off1 += sum(1 for j in range(J + 1) if read[j] and not cnt[j])
        bad = []
        for j in range(J + 1):
            if j == 0:
                law = L.n == 1
                lawp = law
            else:
                law = read[j] or mul(pw[j - 1], bar, p, q) in L
                lawp = read[j] or mul(pw[j - 1], mul(bar, bar, p, q),
                                      p, q) in L
            if j == 0:
                lawr = L.n == 1
            else:
                lawr = any(mul(pw[j - 1], (q, 0), p, q) in L and
                           mul(pw[j - 1], (-r, 1), p, q) in L
                           for r in roots)
            offR += lawr != cnt[j]
            if lawr != cnt[j]:
                badR.append(((L.a, L.c, L.d), j))
            if j == 0:
                lawc = L.n == 1
            else:
                lawc = any(mul(pw[j - 1], (m, 0), p, q) in L and
                           mul(pw[j - 1], (-r, 1), p, q) in L
                           for m, r in mroots)
            offC += lawc != cnt[j]
            if j >= 1:
                lawq = any(mul(pw[max(j - 2, 0)], (m, 0), p, q) in L and
                           mul(pw[max(j - 2, 0)], (-r, 1), p, q) in L
                           for m, r in mroots)
                plantC += lawq != cnt[j]
                if (n ** j) % L.n == 0:
                    live += 1
                    hid += cnt[j]
            offL += law != cnt[j]
            plantL += lawp != cnt[j]
            if law != cnt[j]:
                bad.append((j, cnt[j], law))
        rows.append((L, s, t, bad))
    # the plant on the character criterion
    plant0 = sum(1 for L, s, _, _ in rows
                 if least_hidden_by_characters(L, pwbar, n, J) != s)
    return dict(name=name, n=n, J=J, Is=Is, rows=rows, off0=off0,
                offR=offR, roots=roots, offC=offC, mroots=mroots,
                plantC=plantC, live=live, hid=hid, badR=badR,
                off1=off1, offL=offL, plantL=plantL, plant0=plant0,
                pw=pw, bar=bar, p=p, q=q)


def nonmaximal_at_q(p, q):
    return any(q % (l * l) == 0 and p % l == 0
               for l in range(2, q + 1)
               if all(l % k for k in range(2, l)))


def classify(R):
    tally = {}
    for L, s, t, _ in R["rows"]:
        if s is None and t is None:
            k = "neither"
        elif t is None:
            k = "hidden, unread to J"
        else:
            k = f"gap {t - s}"
        tally[k] = tally.get(k, 0) + 1
    return tally


def prime_of(m):
    """The prime l when m > 1 is a power of l, else None."""
    l = next(k for k in range(2, m + 1) if m % k == 0) if m > 1 else None
    while m > 1 and l and m % l == 0:
        m //= l
    return l if m == 1 else None


def step_dead_ends(R):
    """(primary, mixed): the live psi' = beta psi with no lift outside B,
    psi' != 0 outside B, over the ideals of index dividing n^J of one
    prime, then of two or more."""
    p, q, n, J = R["p"], R["q"], R["n"], R["J"]
    dead = [0, 0]
    for L in R["Is"]:
        N = L.n
        if N == 1 or (n ** J) % N:
            continue
        inB = lambda c: c[0] != 0 and (q * c[0]) % N == 0
        pre = {}
        for A, G in L.characters():
            c = (A % N, G % N)
            pre.setdefault((c[1], (-q * c[0] - p * c[1]) % N), []).append(c)
        bad = sum(1 for t, cs in pre.items() if t != (0, 0) and not inB(t)
                  and all(inB(c) for c in cs))
        dead[prime_of(N) is None] += bad
    return dead


def gap_offs(R, shift):
    """Ideals whose least depth to J is not s - 1 + e + shift, e the
    least e >= 1 with r^e = 0 mod m, O/(I : beta^(s-1)) = Z/m with
    beta -> r; None when r is nilpotent mod no power."""
    p, q, J, pw = R["p"], R["q"], R["J"], R["pw"]
    off = 0
    for L, s, t, _ in R["rows"]:
        if not s:
            continue
        m, r = min(((m, r) for m, r in R["mroots"]
                    if mul(pw[s - 1], (m, 0), p, q) in L and
                    mul(pw[s - 1], (-r, 1), p, q) in L))
        e = next((e for e in range(1, m + 1) if pow(r, e, m) == 0), None)
        want = None if e is None else s - 1 + e + shift
        off += (want if want is not None and want <= J else None) != t
    return off


def main():
    results = [run_base(*b) for b in BASES]

    section("C2  2i against dual2i.py")
    R = next(r for r in results if r["name"] == "2i")
    two = [(L, s, t) for L, s, t, _ in R["rows"]
           if L.n <= 256 and L.n & (L.n - 1) == 0]
    gaps = {}
    for L, s, t in two:
        gaps[t - s] = gaps.get(t - s, 0) + 1
    at = {(L.a, L.c, L.d): (s, t) for L, s, t in two}
    print(f"  {len(two)} ideals, gaps {dict(sorted(gaps.items()))}; "
          f"(4, 2, 1) {at.get((4, 2, 1))}")
    check("C2 the 2i table is dual2i.py's",
          len(two) == 23 and gaps == {0: 13, 1: 10}
          and at.get((4, 2, 1)) == (1, 2))

    section("P0, P1  the character criterion; read implies hidden")
    for R in results:
        print(f"  {R['name']:20s} criterion off {R['off0']}, "
              f"plant off {R['plant0']}, read-not-hidden {R['off1']}, "
              f"{classify(R)}")
    check("P0 the character criterion is the count, every base",
          all(R["off0"] == 0 for R in results))
    check("P0 plant (beta-bar^k) is caught somewhere",
          sum(R["plant0"] for R in results) > 0,
          f"{sum(R['plant0'] for R in results)} off")
    check("P1 read implies hidden, every base",
          all(R["off1"] == 0 for R in results))

    section("P2  the law L at the prime norms")
    for R in results:
        print(f"  {R['name']:20s} L off {R['offL']}, "
              f"plant (beta-bar^2) off {R['plantL']}")
    prime = [R for R in results if R["name"] in PRIME_NORM]
    check("P2 L is the count at every prime-norm base",
          all(R["offL"] == 0 for R in prime))
    check("P2 plant (beta-bar^2) is caught",
          sum(R["plantL"] for R in prime) > 0,
          f"{sum(R['plantL'] for R in prime)} off")

    section("P3  the conjugate prime at the split bases")
    ok3 = True
    for R in prime:
        p, q = R["p"], R["q"]
        Pbar = hnf(R["bar"], mul(R["bar"], (0, 1), p, q))
        P = hnf((0, 1), mul((0, 1), (0, 1), p, q))
        split = (P.a, P.c, P.d) != (Pbar.a, Pbar.c, Pbar.d)
        row = next(r for r in R["rows"]
                   if (r[0].a, r[0].c, r[0].d) == (Pbar.a, Pbar.c, Pbar.d))
        print(f"  {R['name']:20s} {'split' if split else 'ramified'}: "
              f"(beta-bar) = {(Pbar.a, Pbar.c, Pbar.d)} span {row[1]} "
              f"depth {row[2]}")
        if split:
            ok3 &= row[1] == 1 and row[2] is None
        else:
            ok3 &= row[1] == 1 and row[2] == 1
    check("P3 (beta-bar) hidden from span 1, read at no depth, split",
          ok3)

    section("P4, P5  the transplants")
    R = next(r for r in results if r["name"] == "2i")
    split = sum(1 for L, s, t, _ in R["rows"]
                if s is not None and (t is None or s < t))
    print(f"  2i: L off {R['offL']} (ideal, j) pairs; split ideals "
          f"{split}")
    lawread = all((mul(R["pw"][j - 1], R["bar"], R["p"], R["q"]) in L)
                  == (R["pw"][j] in L)
                  for L, _, _, _ in R["rows"] for j in range(1, R["J"] + 1))
    check("P4 2i: L is the read at every ideal", lawread)
    R = next(r for r in results if r["name"] == "-3 + i")
    print(f"  -3 + i: L off {R['offL']} (ideal, j) pairs")
    for L, s, t, bad in R["rows"]:
        for j, cnt, law in bad:
            print(f"    {(L.a, L.c, L.d)} index {L.n} span {j}: "
                  f"hidden {cnt}, L {law}; least span {s}, depth {t}")

    section("K  the kill")
    killed = not (all(R["offL"] == 0 for R in prime)
                  and all(R["off0"] == 0 for R in results))
    print("  L at the prime norms and the criterion: "
          + ("one FIRES" if killed else "both HOLD"))

    section("P6  the sweep, q = 2..12, -1 <= p <= q")
    sweep = []
    for q in range(2, 13):
        for p in range(-1, q + 1):
            sweep.append(run_base(f"x^2 + {p}x + {q}", p, q))
    good = [R for R in sweep if not nonmaximal_at_q(R["p"], R["q"])]
    rest = [R for R in sweep if nonmaximal_at_q(R["p"], R["q"])]
    for R in rest:
        print(f"  l | p, l^2 | q: {R['name']:16s} L off {R['offL']}, "
              f"{classify(R)}")
    worst = max(good, key=lambda R: R["offL"])
    print(f"  the other {len(good)}: L misses at "
          f"{sum(1 for R in good if R['offL'])} of them")
    print(f"  the other {len(good)}: L off at most {worst['offL']} "
          f"({worst['name']}), ideals {sum(len(R['Is']) for R in good)}")
    check("P6 L misses where no l | p has l^2 | q (the kill fired)",
          any(R["offL"] for R in good))
    check("P6 the criterion is the count at every sweep base",
          all(R["off0"] == 0 for R in sweep))
    print("  P6: " + ("FIRES" if worst["offL"] else "HOLDS"))

    section("P7  L', beta^(j-1) J_r for every root r of r(r + p) mod q")
    worst = max(sweep + results, key=lambda R: R["offR"])
    print(f"  {len(sweep)} sweep bases and the {len(results)} named: L' off"
          f" at most {worst['offR']} ({worst['name']}, roots "
          f"{worst['roots']})")
    for R in sweep:
        if len(R["roots"]) > 2 or R["offR"]:
            print(f"  {R['name']:16s} roots {R['roots']}, L' off "
                  f"{R['offR']}, L off {R['offL']}"
                  + (f", misses (ideal, j) {R['badR']}" if R["offR"]
                     else ""))
    check("P7 L' misses at q = 8, 4 | p, alone (the kill fired)",
          {(R["p"], R["q"]) for R in sweep + results if R["offR"]}
          == {(0, 8), (4, 8), (8, 8)})
    print("  P7: " + ("FIRES" if worst["offR"] else "HOLDS"))

    section("P8  the colon law: span j hides I iff span 1 hides "
            "(I : beta^(j-1))")
    worst = max(sweep + results, key=lambda R: R["offC"])
    distinct = {}
    for R in results + sweep:
        distinct.setdefault((R["p"], R["q"]), R)
    live = sum(R["live"] for R in distinct.values())
    hid = sum(R["hid"] for R in distinct.values())
    print(f"  {len(distinct)} distinct bases; {live} (ideal, j >= 1) pairs"
          f" with the index dividing n^j, {hid} of them hidden: off at "
          f"most {worst['offC']} ({worst['name']})")
    check("P8 the colon law is the count at every base",
          all(R["offC"] == 0 for R in sweep + results))
    plant = sum(R["plantC"] for R in distinct.values())
    check("P8 plant (beta^(j-2) for beta^(j-1)) is caught", plant > 0,
          f"{plant} off")
    print("  P8: " + ("FIRES" if worst["offC"] else "HOLDS"))

    section("S  the step: a live character lifts outside B")
    dead = {(R["p"], R["q"]): step_dead_ends(R)
            for R in distinct.values()}
    prim = sum(d[0] for d in dead.values())
    mixed = {k: d[1] for k, d in dead.items() if d[1]}
    print(f"  dead ends at one prime {prim}; at two or more "
          f"{sum(mixed.values())}, at q = "
          f"{sorted({q for _, q in mixed})}")
    check("S the step holds at every ideal of one prime's power",
          prim == 0)
    check("S the step dead-ends at some ideal of mixed index",
          sum(mixed.values()) > 0)

    section("G  the gap: least depth s - 1 + e, r^e = 0 mod m")
    off = sum(gap_offs(R, 0) for R in distinct.values())
    plant = sum(gap_offs(R, 1) for R in distinct.values())
    gaps = {}
    for R in distinct.values():
        for _, s, t, _ in R["rows"]:
            if s and t is not None:
                gaps[t - s] = gaps.get(t - s, 0) + 1
    print(f"  off {off}, plant (s + e) off {plant}; gaps printed "
          f"{dict(sorted(gaps.items()))}")
    check("G the gap law is the least depth at every base", off == 0)
    check("G plant (s + e) is caught", plant > 0, f"{plant} off")

    print()
    print(f"{sum(CHECKS)}/{len(CHECKS)} checks pass")
    return 0 if all(CHECKS) else 1


if __name__ == "__main__":
    sys.exit(main())