towerlemma.py
Checks: THERMAL: The 2-power tower raises the sill theorem; rule.
Cited on: THERMAL.
To run it: save towerlemma.py, with fieldtop.py, growth.py, memory.py, sillclass.py, spectrum.py, splitend.py, splitloss.py, splitsharp.py, tower.py, towerpast.py, towerscan.py, typeloss.py beside it; needs mpmath, numpy, scipy.
"""towerlemma.py -- the 2-power tower's loss side at k = 4 to 15: one
lemma from a d0 near d*, and the certificates between D(k) and d0.
QUESTION. With Q(zeta_(2^k)) = C in K, K_ab = C, type (1, 0), K != C, d
= D m', D = 2^(k-1): (A) at k = 4, does some (c, X) at Y = 10^30 give
the lemma of towerpast.py (x1 the root of d B^-b = c lam, B the least
norm of W_k, w_B = 1) its three bounds g_up, M, M' below 0 with a d0
that the certificates reach? (B) Does one (c, X) close every k = 4 to
15, B = 257 at k = 4..7 (W_k = {257, 65537}) and B = 65537 at k = 8..15
(W_k = {65537}), with the certificate negative at every d = 0 mod D from
D(k) below d0, so that a field of the shape loses at every d = 0 mod D
>= D(k) at each k = 4..15?
THE DERIVATION (by hand, before this code).
The lemma as towerpast.py's (A), its law fixed below Y. With d = c lam
B^b at x1 and y = x1^a, a = log2 B: the places of norm B give log zeta
at most d y/(1 - y) = c lam/(1 - y), those of each other norm n of W_k
c lam (B/n)^b/(1 - n^-b), those past Y (d per prime, each of norm >=
p) c lam (B/Y)^b (1 + Y/(b - 1))/(1 - Y^-b): c_z = coeffs' at X, and
c_s the same with towerpast's coeffs. E, the places off W_k below Y:
a rational p of order f mod 2^k gives D/f primes of C, each with one
place of K over it of norm (p^f)^(m'), p^f = 1 mod 2^k, so p^f >=
2^k + 1 = n_c (the places of norm n_c^(m') among them), at most D per
p, the p^f distinct across p: by the integral test E <= D n_c^-s (1 +
n_c/(s - 1))/(1 - n_c^-s), s = m' b, m' >= d0/D, times (1 + 1/X) for
E + E/x, since (2/q^(m'))^b falls in b and in m' (towerpast's e3 is
this at k = 3: D = 4, n_c = 9). This one E covers the places between
10^5 and Y that tower.BigLaw bounds by its mid term (D per p, norm >=
p^(m')). Then g_up(x1) <= l(X) - (1 - c_z - c_s) lam0 + E', M(x1) <=
l(X) - c lam0 (1 - y0/2) (l_c >= 0 and the other W_k terms dropped),
M'(x1) <= 1/(1 - X) - a c lam0/(X (1 + y0)), lam0 = log(1/X), y0 = X^a,
and every d >= d0 = ceil(c lam0 X^-a) has x1 <= X.
Where the corner is: c_z ~ c and c_s, E', the tail are small, so g_up
< 0 needs (1 - c) lam0 > l(X) and M < 0 needs c lam0 > l(X): both hold
for some c iff l(X) < lam0/2, i.e. 2l + log x < 0, i.e. X < x0 = (3 -
sqrt 5)/2 = 0.381966; at c = 1/2 the two margins are about lam0/2 -
l(X), (l/lam)' = 3.04 at x0, and d0 -> l(x0) x0^-a ~ d* as X -> x0. At
X = 0.381, c = 1/2: l = 0.47965, lam0 = 0.96498, c lam0 = 0.48249;
y0 = 4.4e-4 (B = 257), so M <= -0.0027; c_s ~ c X^(a-1) = 5.8e-4, so
g_up <= 0.47965 - (1 - 0.5002 - 0.0006)(0.96498) ~ -0.0019; M' <=
1.616 - 8.0056 (0.48249)/0.381 ~ -8.5. d0 = 0.48249 e^(8.005625 x
0.964981) = 0.48249 x 2264.4 ~ 1092.5, so 1093 at B = 257; at B =
65537 (a = 16.000022), 0.48249 e^(15.43972) ~ 2.448e6, M <= -0.0028,
g_up ~ -0.0028, M' ~ -19. The tail at Y = 10^30, B = 257, b = 1.392:
c 257^b Y^(1 - b)/(b - 1) ~ 5e-9; at 10^60, B = 65537, ~1e-16: far
below the margins. E' at k = 15 (D = 16384, m' >= 149, s >= 207):
below 10^-900.
The certificates between D(k) and d0: d0 - D(k) is about 21 at k = 4,
so 1072, 1080, 1088 at k = 4, 1072 and 1088 at 5, 1088 at 6 and 7;
about 1.04e5 at k = 8..15: 815 degrees at k = 8 (step 128), 408, 204,
102, 51, 25, 13, 6 at k = 9..15. At k >= 8 the certificate at D(k) is
-5e-6 to -7.5e-4 (towerscan, towerpast) and falls in d. A
TRANSPLANT: that it stays negative over the whole gap, read here.
PREDICTIONS (frozen before the engine).
P1 (A) on the grid X in {0.20, 0.21, ..., 0.38, 0.381}, c in {0.05,
0.10, ..., 0.95} at k = 4, Y = 10^30, the least d0 with g_up, M, M'
< 0 and slope > 0 is at most 1130, at X = 0.381 or 0.38 with c in
[0.45, 0.55]; at (0.381, 1/2): d0 in [1085, 1100], g_up in [-0.003,
-0.001], M in [-0.0035, -0.002], M' <= -8.
P2 (B) at (X, c) = (0.381, 1/2), Y = 10^30 at k = 4..7 and 10^60 at
8..15, the three bounds are below 0 at every k, d0 = 1093 at k =
4..7 and in [2.44e6, 2.46e6] at 8..15.
P3 (B) the certificate (tower.BigLaw at that Y) is negative with M' <
0 and x_c < 1/2 at every d = 0 mod D from D(k) below d0, at every k
= 4..15.
P4 the lemma's bounds hold on the law itself (towerpast's P2) at the
first d = 0 mod D at or past d0, the next, and the first past 2 d0,
at every k = 4..15: x1 <= X, x1 <= x_low, C <= M(x1) <= the M
bound, g_up(x1) and M'(x1) below theirs.
KILL of the k = 4 lemma (frozen beforehand): every (c, X) on the grid
at Y = 10^30 has d0 past 10^5 (or no bounds all negative). KILL of
"a loss at every d = 0 mod D >= D(k)" at a k: P2's bounds not all
negative at that k, or a certificate in [D(k), d0) not negative, or
P4 failing there.
CONTROLS (frozen).
C1 POSITIVE: the lemma code here at k = 3 (B = 17, norms (257, 65537),
c = 1/3, X = 1/4, Y = 10^5, E with D = 4, n_c = 9) reproduces
towerpast.lemma3() (c_z, c_s, g_up, M, M', d0, E') to 1e-12.
C2 PLANT at the derivation's corner: at k = 4, X = 0.39 > x0, no c on
the grid gives all three bounds below 0; and at (0.381, c = 1),
the slope 1 - c_z - c_s is at most 0.
C3 the coefficients against the law's own sums (tower.BigLaw's parts)
at x = x1 and the first d past d0, at k = 4 (Y = 10^30) and k = 8
(10^60): lz + lt <= c_z lam + E and (lz - sw + lt)/x <= c_s lam +
E/x.
C4 POSITIVE: the certificate engine reproduces the printed values
-0.000483 at (4, 1072, 10^30) (tower.arm_y), -4.827e-04 at (5,
1072, 10^30) and -1.128e-05 at (9, 2344192, 10^60) (towerscan),
-4.973e-06 at (8, 2344064, 10^60) (towerpast), to the printed
digits.
C5 towerscan's order2 replaces tower.order_mod (equal at k = 5, 9 on
3000 primes there), checked again here at k = 15 on 500 primes.
FINDINGS (entered after the run, from its printed output).
CONTROLS. C1: 0.0, d0 134, E' 2.644e-63. C2: no c on the grid closes
at X = 0.39; at (0.381, c = 1) the slope is -0.002047, so the lemma
can say no at both corners. C3 FAILED on its first run, at k = 8, x1
only, by the float wall: lz - sw subtracts two sums near 0.48 (error
about 1e-16), while the norm-B term's bound c lam x^(a-1)/(1 - y^2)
exceeds its value d(-log(1 - y^2))/x by about 5e-21 there (40-digit
check, a gap of order y^4, y = x^16 ~ 2e-7); a diagnostic print read
lhs - rhs +8.5e-17 there. lzsw, the same sum place by place with no cancellation,
replaced the subtraction (a change to the frozen control, named here)
and C3 passes: at k = 4, d = 1096, lz + lt 0.48297431 <= 0.48308184
(at X) and <= 0.48308064 (at x1), the c_s side 5.5808e-04 against
5.5945e-04 and 5.5808e-04; at k = 8, d = 2447360, 0.48247873 <=
0.48247877, and 2.4965e-07 against 2.4966e-07 and 2.4965e-07 (the
norm-B term all but exact at x1, as at k = 3). C4: -4.8282e-04,
-4.8274e-04, -1.1279e-05, -4.9729e-06. C5: no mismatch. At review:
C1 passes k = 3's arguments to towerpast.lemma by hand and e_off at
(4, 9) is e3's formula, so its 0.0 reads e_off and lemma, never
lemma_k's choice of B and norms (read instead by (A)'s d0 1093 against
the hand 1092.5); C2 is caught by construction (at X > x0, M < 0 and
g_up < 0 sum to 2 l(X) < lam0, false; at c = 1, c_z > 1), a wiring
check, not a plant; C3 allows E + E/x on both sides where the frozen
text split E and E/x, and its c_s half at k = 8, x1 passes by 5e-21,
inside brentq's rounding of x1: the lemma reads the coefficients at
X, where the margin is about 1e-11.
(A) P1 HOLDS: the least d0 on the grid with all three bounds below 0
is 1093, at (X, c) = (0.381, 0.5); there c_z 0.500444, c_s 0.000580,
slope 0.498976, E' 2.934e-233, g_up <= -0.001840, M <= -0.002721, M'
<= -8.5179. KILL of the k = 4 lemma MISSED.
(B) P2 HOLDS: at (0.381, 1/2) the three bounds are below 0 at every k
= 4..15, d0 = 1093 at k = 4..7 (g_up, M, M' as at k = 4; E' 3.0e-143,
1.7e-84, 1.0e-47 at k = 5, 6, 7) and 2447300 at k = 8..15 (g_up and
M printed -0.002828, rounded away from 0: at full
precision lemma_k reads -0.0028276009 and -0.0028278981, so g_up and M
< -0.002827; M' <= -18.646). E' prints 0.0 at k = 8..15,
an underflow: by hand, n_c^-s with s = m' b >= (2447300/D)(1.392) is
257^-26600 at k = 8 and 32769^-208 = 10^-939 at k = 15, so E' < 10^-900
at every k = 8..15, no print moved. P3 HOLDS: the certificate is
negative with M' < 0 and x_c < 1/2 at every d = 0 mod D from D(k) below
d0: 3 degrees at k = 4 (1072, 1080, 1088), 2 at 5, 1 at 6 and 7, 807,
403, 201, 100, 50, 25, 12, 6 at k = 8..15, its largest at D(k) each
(-4.828e-04, -4.827e-04, -3.398e-03, -3.398e-03, -4.973e-06,
-1.128e-05, -2.389e-05, -4.911e-05, -4.911e-05, -1.499e-04, -3.513e-04,
-7.530e-04). P4 HOLDS at the three degrees read at every k. KILL of
"a loss at every d = 0 mod D >= D(k)" MISSED at every k = 4..15: with
the lemma from d0 and the certificates below it, a field of the shape
loses at every d = 0 mod D >= D(k), at each k = 4..15, so with the
no-loss half (towerscan) the losing degrees at each k = 2..15 are
exactly d >= D(k).
RUN RECORD. `python towerlemma.py` under a 2048 MB cap: the first run
stopped at C3's guard (4 of 5 checks, 0.6 s, 59.2 MB; C1, C2, C4, C5
printing as here). After lzsw: 5 of 5 checks, 94.8 s, peak commit 59.1
MB. Moved in from the lab unchanged, run bare through the gate: every
line printed as before, 5 of 5 checks, 95.8 s, 59.9 MB.
"""
import math
import os
import sys
import time
HERE = os.path.dirname(os.path.abspath(__file__))
sys.path.insert(0, os.path.join(HERE, "..", "code"))
os.environ.setdefault("OPENBLAS_NUM_THREADS", "1")
import numpy as np # noqa: E402
import splitsharp as ss # noqa: E402
import tower as tw # noqa: E402
import towerpast as tp # noqa: E402
import towerscan as ts # noqa: E402
CHECKS = []
XS, CS = 0.381, 0.5
GRID_X = [round(0.20 + 0.01 * i, 2) for i in range(19)] + [0.381]
GRID_C = [round(0.05 * i, 2) for i in range(1, 20)]
DK = {4: 1072, 5: 1072, 6: 1088, 7: 1088, 8: 2344064, 9: 2344192,
10: 2344448, 11: 2344960, 12: 2344960, 13: 2347008, 14: 2351104,
15: 2359296}
def check(name, ok, detail=""):
CHECKS.append(bool(ok))
print(f" [{'PASS' if ok else 'FAIL'}] {name}" + (f" ({detail})"
if detail else ""))
def guard():
if not all(CHECKS):
print("a control failed: no verdict")
sys.exit(1)
def ybig(k):
return 10 ** 30 if k <= 7 else 10 ** 60
def wprimes(k):
"""The primes of W_k: the Fermat primes past n_c, 1 mod 2^k."""
return [n for n in (17, 257, 65537)
if n > tw.nc(k) and (n - 1) % 2 ** k == 0]
def e_off(D, n, d0, X):
"""E + E/x: at most D places per p off W_k, of norm (p^f)^(m'), p^f
>= n distinct across p, m' >= d0/D, at x <= X."""
s = (d0 / D) * -math.log2(X)
return D * n ** -s * (1 + n / (s - 1)) / (1 - n ** -s) * (1 + 1 / X)
def lemma_k(k, X, c, Y):
W = wprimes(k)
B, norms = W[0], tuple(W[1:])
a = math.log2(B)
D = 2 ** (k - 1)
L = tp.lemma(X, c, 1.0, a, B, norms, Y, 0.0)
E = e_off(D, tw.nc(k), L["d0"], X)
L = tp.lemma(X, c, 1.0, a, B, norms, Y, E)
L.update(B=B, a=a, norms=norms, D=D)
return L
def good(L):
return L["G"] < 0 and L["M"] < 0 and L["Mp"] < 0 and L["slope"] > 0
def tlaw(k, d):
return tw.TowerLaw(k, tw.BigLaw(tw.law(k, d), ybig(k), 2 ** (k - 1)))
def lzsw(base, x):
"""lz - sw of base.parts(x), place by place with no cancellation:
-log(1 - u^2) at a column place, -log(1 - u) elsewhere (added after
C3's first run, the float wall at k = 8)."""
u = np.exp(math.log(x) * base.a)
col = base.k > 0
return float((base.c[col] * -np.log1p(-u[col] ** 2)).sum()
+ (base.c[~col] * -np.log1p(-u[~col])).sum())
# ------------------------------------------------------------- controls
def controls():
print("=== CONTROLS")
want = tp.lemma3()
L0 = tp.lemma(tp.X3, tp.H3, 1.0, tp.A17, 17, (257, 65537), ss.Y, 0.0)
got = tp.lemma(tp.X3, tp.H3, 1.0, tp.A17, 17, (257, 65537), ss.Y,
e_off(4, 9, L0["d0"], tp.X3))
err = max(abs(got[q] - want[q]) for q in ("cz", "cs", "G", "M", "Mp",
"E"))
check("C1 POSITIVE: k = 3's lemma reproduces towerpast.lemma3()",
err < 1e-12 and got["d0"] == want["d0"],
f"{err:.1e}, d0 {got['d0']}, E' {got['E']:.3e}")
any39 = [c for c in GRID_C if good(lemma_k(4, 0.39, c, ybig(4)))]
P = lemma_k(4, XS, 1.0, ybig(4))
check("C2 PLANT: at k = 4, X = 0.39 > x0 no c on the grid closes; at "
"(0.381, c = 1) the slope is <= 0",
not any39 and P["slope"] <= 0,
f"closing c at 0.39: {any39}; slope at c = 1 {P['slope']:+.6f}")
ok3 = True
for k in (4, 8):
L = lemma_k(k, XS, CS, ybig(k))
D = L["D"]
d = D * math.ceil(L["d0"] / D)
q = tlaw(k, d)
x = tp.x1(d, CS, L["a"])
lam = -math.log(x)
lz, sw, w, lt = q.base.parts(x)
ls = lzsw(q.base, x)
E = e_off(D, tw.nc(k), d, x)
for tag, xx in (("X", XS), ("x1", x)):
cz, cs = tp.coeffs(xx, CS, 1.0, L["a"], L["B"], L["norms"],
ybig(k))
g1 = lz + lt <= cz * lam + E
g2 = (ls + lt) / x <= cs * lam + E
ok3 = ok3 and g1 and g2
print(f" k = {k}, d = {d}, coefficients at {tag}: lz + lt "
f"{lz + lt:.8f} <= {cz * lam + E:.8f}; (lz - sw + lt)/x "
f"{(ls + lt) / x:.4e} <= {cs * lam + E:.4e}")
check("C3 the coefficients bound the law's own sums at the first d "
"past d0, k = 4 and 8", ok3)
want4 = ((4, 1072, -0.000483, 6), (5, 1072, -4.827e-04, 7),
(9, 2344192, -1.128e-05, 8), (8, 2344064, -4.973e-06, 9))
bad, seen = [], []
for k, d, v, digits in want4:
C = tlaw(k, d).cert()[0]
seen.append(f"{C:.4e}")
if round(C, digits) != round(v, digits):
bad.append((k, d, C))
check("C4 POSITIVE: the certificate engine reproduces four printed "
"values", not bad, f"{', '.join(seen)}")
bad = [p for p in tw.PR_Y[:500]
if ts.order2(p, 2 ** 15) != tw.order_mod(p, 2 ** 15)]
check("C5 order2 is tower.order_mod at k = 15 on 500 primes", not bad,
f"{bad[:4]}")
guard()
tw.order_mod = ts.order2
# ----------------------------------------------------------------- arms
def part_a():
print("\n=== (A) k = 4: the grid at Y = 10^30")
best = None
for X in GRID_X:
for c in GRID_C:
L = lemma_k(4, X, c, ybig(4))
if good(L) and (best is None or L["d0"] < best[0]):
best = (L["d0"], X, c)
print(f" least d0 with g_up, M, M' < 0: {best}")
L = lemma_k(4, XS, CS, ybig(4))
print(f" at (X, c) = ({XS}, {CS}): d0 {L['d0']}, c_z {L['cz']:.6f}, "
f"c_s {L['cs']:.6f}, slope {L['slope']:.6f}, E' {L['E']:.3e}; "
f"g_up <= {L['G']:+.6f}, M <= {L['M']:+.6f}, M' <= "
f"{L['Mp']:+.4f}")
p1 = (best is not None and best[0] <= 1130 and best[1] in (0.38, 0.381)
and 0.45 <= best[2] <= 0.55 and 1085 <= L["d0"] <= 1100
and -0.003 <= L["G"] <= -0.001 and -0.0035 <= L["M"] <= -0.002
and L["Mp"] <= -8)
kill = best is None or best[0] > 10 ** 5
print(f" P1 {'holds' if p1 else 'FAILS'}; KILL of the k = 4 lemma "
f"{'fires' if kill else 'missed'}")
def part_b():
print(f"\n=== (B) k = 4..15 at (X, c) = ({XS}, {CS})")
p2 = p3 = p4 = True
closed = []
for k in range(4, 16):
L = lemma_k(k, XS, CS, ybig(k))
D = L["D"]
Dk = ts.first_past(ts.dh(k), D)
assert Dk == DK[k], (k, Dk)
ok2 = good(L)
if k >= 8:
ok2 = ok2 and 2.44e6 <= L["d0"] <= 2.46e6
else:
ok2 = ok2 and L["d0"] == 1093
p2 &= ok2
bad, worst, n = [], None, 0
for d in range(Dk, L["d0"], D):
C, dm, xl, xc = tlaw(k, d).cert()
n += 1
if not (C < 0 and dm < 0 and xc < 0.5):
bad.append(d)
if worst is None or C > worst[0]:
worst = (C, d)
p3 &= not bad
d1 = D * math.ceil(L["d0"] / D)
d3 = D * (2 * L["d0"] // D + 1)
bad4 = []
for d in (d1, d1 + D, d3):
q = tlaw(k, d)
C, dm, xl, xc = q.cert()
x = tp.x1(d, CS, L["a"])
if not (x <= XS and x <= xl and C <= q.m(x) <= L["M"]
and q.g_up(x) <= L["G"] and q.dm(x) <= L["Mp"]):
bad4.append(d)
p4 &= not bad4
if ok2 and not bad and not bad4:
closed.append(k)
top = f"{worst[0]:+.3e} at {worst[1]}" if worst else "none"
print(f" k = {k:2d}: B = {L['B']}, Y = 10^{30 if k <= 7 else 60}, "
f"d0 {L['d0']}, g_up <= {L['G']:+.6f}, M <= {L['M']:+.6f}, "
f"M' <= {L['Mp']:+.3f}, E' {L['E']:.1e}; certificates at "
f"the {n} d = 0 mod {D} from {Dk} below d0: not negative at "
f"{bad[:4]}, largest {top}; P4 at {d1}, {d1 + D}, {d3}: "
f"failing at {bad4}")
print(f" P2 {'holds' if p2 else 'FAILS'}; P3 "
f"{'holds' if p3 else 'FAILS'}; P4 {'holds' if p4 else 'FAILS'}")
print(f" a loss at every d = 0 mod D >= D(k) at k = {closed}")
def main():
t0 = time.time()
controls()
part_a()
part_b()
print(f"\n{sum(CHECKS)} of {len(CHECKS)} checks, {time.time() - t0:.1f} s")
return 0 if all(CHECKS) else 1
if __name__ == "__main__":
sys.exit(main())