"""dam.py -- a slack max(ct, t log2(ln t) / 8) is dammed from every seed;
the covering's exponent L_J divides 2 lcm(1..J) <= 2 4^(J-1).

QUESTION. Every c-greedy walk over Z is dammed, and a nondecreasing
slack at least p(t) on the prime powers lets the record run escape from
the void. Is a slack a little past every linear one,
f(t) = max(ct, kappa t ln(ln t)), dammed from every seed, and at which
kappa does the covering prove it?

DESIGN. The f-greedy walk, f nondecreasing with f(t) >= t, takes any
lambda-raising m <= f(m0). (1) THE DAM. If f(tau(P)) < P+, every run
from a seed whose primes are <= P keeps them <= P (every move is at most
f(m0) <= f(tau(P))); at infinitely many P, every seed is dammed.
g(t) = max(ct, kappa t ln(ln t)), c >= 1, kappa > 0, is increasing on
t > 1, since (t ln(ln t))' = ln u + 1/u >= 1 with u = ln t.
(2) THE COVERING (cgreedy.py's L2 and T2). At e = 1 + L_J, J >= 4,
L_J the lcm of the orders of 2 mod 2j + 1 over j <= J, every j 2^e + 1
with j <= J is composite, so at the largest prime P <= J 2^e,
tau(P) <= 2^(e+2) and P+ > J 2^e; for f <= g past some t the dam holds
there, J large, when g(2^(e+2)) <= J 2^e, that is when 4c <= J and
    (D)  kappa ln((L_J + 3) ln 2) <= J/4.
K(J) = J / (4 ln((L_J + 3) ln 2)) is the largest kappa (D) admits.
(3) THE BOUND. ord_n(2) | lambda(n), the Carmichael function, even with
lambda(n)/2 <= (n - 1)/2 = j at n = 2j + 1, so L_J | 2 lcm(1..J). And
lcm(1..x) <= 4^(x-1) at every x >= 1: each prime power p^a in
(ceil(x/2), x] has Legendre term 1 at its own exponent a in the
exponent of p in C(x, floor(x/2)) (every term 0 or 1), and no two share
a prime, so the product of their primes, exp(psi(x) - psi(ceil(x/2))),
divides C(x, floor(x/2)) <= 4^floor(x/2); the steps
x_i - ceil(x_i/2) = floor(x_i/2) telescope to x - 1. So
(L_J + 3) ln 2 < 2 4^(J-1) + 3 <= 4^J at J >= 2, and (D) holds at every
J >= 4 when kappa <= 1/(4 ln 4) = 1/(8 ln 2): the dam at P_J, the
largest prime <= J 2^(1 + L_J), for every J >= max(4, 4c) with
f <= g at 2^(L_J + 3) (every large J), P_J unbounded, so every seed is dammed, at every c.
(4) THE LIMIT. The same chain gives ln((L_J + 3) ln 2) < ln 3 + psi(J)
at J >= 3, and psi(J) ~ J (the prime number theorem, known: Hadamard and
de la Vallee Poussin, 1896), so K(J) >= (1 - o(1))/4: every kappa < 1/4
meets (D) at every large J. (4) is proof only; no run reaches it.

BY HAND. L_4 = 12, K(4) = 0.427; L_7 = 60, K(7) = 0.463; L_11 = 3960,
K(11) = 0.347; L_29 = 240,360,120 (cgreedy.py's print), K(29) = 0.383.

PREDICTIONS (frozen before the code).
Q1 (forced) L_J | 2 lcm(1..J) at every 1 <= J <= 40.
Q2 (forced) lcm(1..x) <= 4^(x-1) at every 1 <= x <= 2000.
Q3 (forced) ln((L_J + 3) ln 2) < J ln 4, and (D) holds at
    kappa = 1/(8 ln 2), at every 4 <= J <= 40.
Q4 (a reading) max K(J) over 4 <= J <= 40 lies in [0.40, 0.60], near
    J = 7; ln L_J / J at J = 40 in [0.5, 0.8].

KILL. "C5 FAIL" printed at a J or an x where Q1, Q2 or Q3 breaks.

CONTROLS, read first. C6: L_4 = 12 and L_29 = 240,360,120; 2j + 1 |
j 2^e + 1 at e = 1 + L_J for every j <= J <= 40. C7, planted breaks the
tests must flag: Q1's test on 59 L_29; Q3's two conjuncts apart, Q3a
(ln((L_J + 3) ln 2) < J ln 4) on L_20 = 2 4^20, where Q3b at
kappa = 0.1 passes, and Q3b ((D)) at kappa = 1 and at kappa = 0.34, just
above the least K(J), where it must flag J = 5 alone; Q2's test against
4^(x-1) - 1 and against 2^(x-1), which fires at large x.

FINDINGS (from the printed run; the controls passed first, C7 flagging
59 L_29 at J = 29, kappa = 1 at 37 J from J = 4 with Q3a silent there,
the planted L_20 = 2 4^20 under Q3a alone, kappa = 0.34 at J = 5
alone, 4^(x-1) - 1 at x = 1, and 2^(x-1) at 1998 x, 3 to 2000). KILL silent: Q1 holds at every
1 <= J <= 40, Q2 at every 1 <= x <= 2000, Q3 at every 4 <= J <= 40. Q4
held its bands but not its place: max K(J) = 0.5576 at J = 22, not near
7 (K(7) = 0.4634), the least 0.3310 at J = 5; ln L_J / J runs 0.51 to
0.82 over the printed J, 0.5099 at J = 40, where L_40 = 721,080,360 =
3 L_29 (the order 54 of 2 mod 81). The run checks the proof's steps at
J <= 40 and x <= 2000; the dam from every seed at kappa <= 1/(8 ln 2) =
0.18034, at every c, is (3)'s proof, and every kappa < 1/4 is (4)'s. A
kappa up to K(J) is dammed at the one P near J 2^(1 + L_J) for the
J <= 40 computed, and from every seed only below 1/4 as proved.
psi(x)/x is 0.9972 at x = 2000, against the bound's ln 4 = 1.3863.
RUN RECORD. 8/8 checks, 0.0 s, peak commit 12.8 MB. A second run, Q3
split into Q3a and Q3b and the plants on L_20, at kappa = 0.34 and
against 2^(x-1) added: 11/11, 0.0 s, peak commit 12.8 MB, every earlier
print unchanged.
"""

import os
import time
from math import gcd, log

os.environ.setdefault("OPENBLAS_NUM_THREADS", "1")
import cgreedy  # noqa: E402

CHECKS = []
JMAX = 40
XMAX = 2000
KAPPA = 1 / (8 * log(2))


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


def lcm_table(n):
    """lcm(1..x) for 0 <= x <= n (lcm of the empty range is 1)."""
    out, L = [1], 1
    for x in range(1, n + 1):
        L = L * x // gcd(L, x)
        out.append(L)
    return out


def l_table():
    """L_J for 1 <= J <= JMAX, by cgreedy.order2."""
    out, L = {}, 1
    for J in range(1, JMAX + 1):
        o = cgreedy.order2(2 * J + 1)
        L = L * o // gcd(L, o)
        out[J] = L
    return out


def lhs(L):
    return log((L + 3) * log(2))


def q1_bad(LJ, lcm):
    return [J for J in LJ if (2 * lcm[J]) % LJ[J]]


def q2_bad(lcm, bound=lambda x: 4 ** (x - 1)):
    return [x for x in range(1, XMAX + 1) if lcm[x] > bound(x)]


def q3a_bad(LJ):
    return [J for J in range(4, JMAX + 1) if not lhs(LJ[J]) < J * log(4)]


def q3b_bad(LJ, kappa):
    return [J for J in range(4, JMAX + 1) if not kappa * lhs(LJ[J]) <= J / 4]


def section_controls(LJ, lcm):
    print("CONTROLS")
    check("C6 L_4 = 12, L_29 = 240,360,120", LJ[4] == 12
          and LJ[29] == 240_360_120, f"L_4 = {LJ[4]}, L_29 = {LJ[29]}")
    bad = [(J, j) for J in LJ for j in range(1, J + 1)
           if (j * pow(2, 1 + LJ[J], 2 * j + 1) + 1) % (2 * j + 1)]
    check("C6 2j + 1 | j 2^e + 1, e = 1 + L_J, j <= J <= 40", not bad,
          f"{len(bad)} failures")
    planted = dict(LJ)
    planted[29] = 59 * LJ[29]
    check("C7 Q1's test flags 59 L_29", q1_bad(planted, lcm) == [29])
    check("C7 Q3b's test at kappa = 0.34 flags J = 5 alone",
          q3b_bad(LJ, 0.34) == [5])
    b = q3b_bad(LJ, 1.0)
    check("C7 Q3b's test flags kappa = 1, Q3a's passes the true L_J",
          4 in b and not q3a_bad(LJ), f"flags {len(b)} J, "
          f"from {b[0] if b else None}")
    planted = dict(LJ)
    planted[20] = 2 * 4 ** 20
    check("C7 Q3a's test flags L_20 = 2 4^20, Q3b's at kappa = 0.1 passes",
          q3a_bad(planted) == [20] and not q3b_bad(planted, 0.1))
    b = q2_bad(lcm, lambda x: 4 ** (x - 1) - 1)
    check("C7 Q2's test flags the bound 4^(x-1) - 1", 1 in b,
          f"flags {len(b)} x, from {b[0] if b else None}")
    b = q2_bad(lcm, lambda x: 2 ** (x - 1))
    check("C7 Q2's test flags the bound 2^(x-1) up to x = 2000",
          XMAX in b, f"flags {len(b)} x, {b[0] if b else None} to "
          f"{b[-1] if b else None}")


def section_dam(LJ, lcm):
    print("THE LOWER DAM")
    b1, b2 = q1_bad(LJ, lcm), q2_bad(lcm)
    b3 = sorted(set(q3a_bad(LJ)) | set(q3b_bad(LJ, KAPPA)))
    for J in b1[:3]:
        print(f"  C5 FAIL Q1 at J = {J}")
    for x in b2[:3]:
        print(f"  C5 FAIL Q2 at x = {x}")
    for J in b3[:3]:
        print(f"  C5 FAIL Q3 at J = {J}")
    check(f"Q1 L_J | 2 lcm(1..J) at every 1 <= J <= {JMAX}", not b1)
    check(f"Q2 lcm(1..x) <= 4^(x-1) at every 1 <= x <= {XMAX}", not b2)
    check(f"Q3 ln((L_J + 3) ln 2) < J ln 4 and (D) at kappa = 1/(8 ln 2)"
          f" = {KAPPA:.5f}, every 4 <= J <= {JMAX}", not b3)
    print("  J    L_J                     ln L_J   ln L_J/J  J ln 4   K(J)")
    for J in [4, 5, 7, 8, 10, 11, 16, 20, 29, 30, 35, 40]:
        L = LJ[J]
        print(f"  {J:2d}   {L:<22d}  {log(L):7.3f}   {log(L) / J:6.3f}"
              f"   {J * log(4):6.2f}   {J / (4 * lhs(L)):.4f}")
    ks = {J: J / (4 * lhs(LJ[J])) for J in range(4, JMAX + 1)}
    top = max(ks, key=ks.get)
    low = min(ks, key=ks.get)
    print(f"  Q4 reading: max K(J) = {ks[top]:.4f} at J = {top}; least "
          f"{ks[low]:.4f} at J = {low}; ln L_40 / 40 = "
          f"{log(LJ[JMAX]) / JMAX:.4f}")
    psi = [log(lcm[x]) / x for x in (100, 500, 1000, 2000)]
    print("  psi(x)/x at x = 100, 500, 1000, 2000: "
          + ", ".join(f"{v:.4f}" for v in psi) + f" (ln 4 = {log(4):.4f})")
    return bool(b1 or b2 or b3)


def main():
    t0 = time.time()
    lcm = lcm_table(XMAX)
    LJ = l_table()
    section_controls(LJ, lcm)
    if not all(ok for _, ok in CHECKS):
        print("CONTROL FAILED: nothing below is read.")
        raise SystemExit(1)
    fired = section_dam(LJ, lcm)
    print("KILL")
    print(f"  {'FIRED' if fired else 'silent'}")
    print(f"{sum(ok for _, ok in CHECKS)}/{len(CHECKS)} checks passed, "
          f"{time.time() - t0:.1f} s")
    raise SystemExit(0 if all(ok for _, ok in CHECKS) else 1)


if __name__ == "__main__":
    main()