cubrace.py
Checks: THERMAL: The race over cubic fields theorem; rule.
Cited on: THERMAL.
To run it: save cubrace.py, with cubic.py, cubtop.py, fieldtop.py, growth.py, memory.py, quadtop.py, spectrum.py beside it; needs mpmath, numpy, scipy.
"""cubrace.py -- the race between the F_2 place and the column of norm 3
read at a cubic field, where cubtop.py's bound cannot decide it.
QUESTION. cubtop.py's bound fails only at the column of norm 3, where 2
is 1^3 or 1^2 1 and the F_2 column holds two or three places of norm
5. Ten complex cubic fields with |d_K| <= 3000 meet a failing
configuration. At each, does that column's point lie above the F_2
place's?
THE ARGUMENT. By fieldtop.py's (1), beta_P > beta_F2 iff sigma_P <
sigma_2 at beta_F2; both points are read too. The column P of norm 3
(N - 1 = 2, p = 3) holds, by spectrum.py's (5), exactly: a place over
2 of residue degree 1 and e = 2 or 3 to depth 2 (lambda(Q^2) = 2;
1 + pi squares to valuation exactly 2, so lambda(Q^3) = 4), one with
e = 1 to depth 3 ((Z/8)^* has exponent 2, (Z/16)^* 4), a place of norm
4 to depth 2 (lambda = 3 2^(j-1)); every other place of norm 3 without
bound, none of norm 9 or 27; a place over an odd r != 3 of norm M to
depth 1 where M - 1 | 2 3^oo (lambda(R^2) carries r), else 0. Below
10^6 such M are 7, 19, 163, 487, 1459 and the primes 2 3^a + 1 past
them; past 10^6 at most three places each, a mass printed. The F_2
column: T2 = (1 - x)^-(g - 1) times 1 + N^-b at each place of Fermat
norm (fieldtop.py's (2)); the Fermat primes past 65537 are past 10^6
and at least 2^32 + 1, a factor below exp(6 2^-32b), carried by the F_2
column's lower end.
DESIGN. The splitting of every p < 10^6: at p | a d_K, cubic.py's
places; else the roots of g(x) = a x^3 - b x^2 + c x - d mod p,
deg gcd(x^p - x, g) = 3, 1, 0 for 111, 12, 3. zeta_K is the Euler
product below 10^6 times T^m, T = zeta(b) prod_(p < 10^6) (1 - p^-b)
the tail of zeta's product: each local factor past 10^6 lies in
[1, (1 - p^-b)^-3] (its places have sum e f = 3, at every cubic
field, ramified or not), so m = 0 and m = 3 bracket zeta_K (m = 1 the
central reading, every split law averaging one place of norm p). Each
point is read at m = 0 and m = 3 by bisection to 1e-9. VERDICT per
field: LOSS where N=3's bracket lies above F2's, TOP where it lies
below, else OPEN; (R) checks the tally as found. HUNT (H): every
complex cubic field with |d_K| <= 3000 (cubic.py's enumerate_fields), its types at 2, 3, 5, 17, 257,
65537 read off its places, run through cubtop.py's bracket at every
type of 2 with an F_2 place; the fields in a failing configuration
must be FIELDS, exactly.
CONTROL (C) the same engine with chi_D splitting at D = -8 (2
ramified, 3 split): its brackets hold fieldtop.py's F2 1.63177 and
N=3 1.49122, the central reading within 2e-4 of each.
PREDICTIONS, fixed before the first run. (r1) TOP at -460. (r2) LOSS
at one or more of -1836, -2284, -2764, -2956. (r3) no OPEN.
FINDINGS. C: F2 [1.63175, 1.63179] central 1.63177, N=3 [1.49110,
1.49135] central 1.49122. The mass past 10^6 of the column of norm 3,
1.40e-08 at b = 1.4. H: 419 fields read, 10 in a failing
configuration, each failing at the column of norm 3 alone, = FIELDS.
RACE: LOSS 2 [-1420, -2764], TOP 8, OPEN 0.
d_K -1420 (g = 3x^3 + 7x^2 + 7x + 5, 2 = Q^3): F2 [1.59661, 1.59670],
N=3 [1.59889, 1.59897], central gap -0.00228. The F_2 column holds
places of norm 3, 5, 5, 9, 17, 17, 17, 65537; the column of norm 3
holds Q to depth 2 and one place of norm 1459.
d_K -2764 (g = 2x^3 + 6x^2 + 8x + 9, 2 = Q^3): F2 [1.63685, 1.63690],
N=3 [1.65592, 1.65596], central gap -0.01907; F_2 column 3, 5, 5,
5, 9, 17, 257; held 19, 163, 487, 1459, 39367.
TOP, central gap: -460 +0.02012, -940 +0.04082, -1836 +0.02276,
-2284 +0.00432, -2380 +0.06134, -2484 +0.00846 (2 = Q1^2 Q2), -2860
+0.04357, -2956 +0.01071.
K: both losing fields pass cubic.py's class group check against the
analytic class number formula (h = 1 at -1420, h = 2 at -2764). N=3's
lower end carries the bound on its places past 10^6.
(r1) holds. (r2) holds at -2764 and misses -1420, a loss with the
bound's ratio 0.960832, above three that top. (r3) holds.
READING. The place with residue field F_2 loses the top over a cubic
field; the least |d_K| complex cubic where it does is -1420, since
cubtop.py's bound decides every column of every configuration but the
column of norm 3 at its 675 failing ones, and of the fields below 1420
only -460 and -940 meet one, both topping. A loss needs what the
bound's failure needs, 2 of type 1^2 1 or 1^3 and two or three places
of norm 5 in the F_2 column; both losses read have 2 totally ramified,
places of norm 5 and 9 in the F_2 column, and few places of norm
2 3^a + 1 to feed the column of norm 3 (one at -1420).
RUN RECORD. `python cubrace.py`: 5 of 5 checks, 33.1 s, peak commit
79.3 MB. With the F_2 lower end carrying the Fermat factor: 5 of 5,
33.1 s, 79.1 MB, every bracket and verdict reprinted.
"""
import os
import sys
from math import exp, log
os.environ.setdefault("OPENBLAS_NUM_THREADS", "1")
sys.path.insert(0, os.path.dirname(os.path.abspath(__file__)))
import numpy as np # noqa: E402
from scipy.special import zeta # noqa: E402
from cubic import Field, enumerate_fields, primes_below # noqa: E402
from fieldtop import kronecker # noqa: E402
X = 10 ** 6
PRIMES = primes_below(X)
FERMAT = {3, 5, 9, 17, 257, 65537}
FIELDS = {-460: (1, -1, 5, -3), -940: (1, -3, 1, -5),
-1420: (3, -7, 7, -5), -1836: (1, -6, 6, -6),
-2284: (2, -2, 6, -5), -2380: (1, -4, 4, -10),
-2484: (4, -3, 6, -2), -2764: (2, -6, 8, -9),
-2860: (2, -2, 0, -5), -2956: (1, -3, 7, -15)}
CHECKS = []
def check(name, ok, detail=""):
CHECKS.append(bool(ok))
print(f" [{'PASS' if ok else 'FAIL'}] {name}" + (f" ({detail})"
if detail else ""))
def roots_mod(p, a, b, c, d):
"""deg gcd(x^p - x, g) mod p, g = a x^3 - b x^2 + c x - d, p | a
excluded: the number of roots of g in F_p."""
inv = pow(a, p - 2, p)
c2, c1, c0 = (-b * inv) % p, (c * inv) % p, (-d * inv) % p
def mul(u, v):
w = [0] * 5
for i in range(3):
for j in range(3):
w[i + j] += u[i] * v[j]
for k in (4, 3): # x^3 = -c2 x^2 - c1 x - c0
t = w[k] % p
w[k] = 0
w[k - 1] -= t * c2
w[k - 2] -= t * c1
w[k - 3] -= t * c0
return [w[0] % p, w[1] % p, w[2] % p]
r, base, e = [1, 0, 0], [0, 1, 0], p
while e:
if e & 1:
r = mul(r, base)
base = mul(base, base)
e >>= 1
h = [r[0], (r[1] - 1) % p, r[2]]
f = [c0, c1, c2, 1]
while any(h):
while h and h[-1] == 0:
h.pop()
if not h:
break
# f mod h
lead = pow(h[-1], p - 2, p)
f = f[:]
while len(f) >= len(h):
q = f[-1] * lead % p
s = len(f) - len(h)
for i, hc in enumerate(h):
f[s + i] = (f[s + i] - q * hc) % p
f.pop()
while f and f[-1] == 0:
f.pop()
f, h = h, f
while f and f[-1] == 0:
f.pop()
return len(f) - 1
def cubic_places(form, dk):
"""[(p, f, e)] for every place over p < X."""
K = Field(form, dk)
a, b, c, d = form
out = []
for p in PRIMES:
if a % p == 0 or dk % p == 0:
out += [(p, P.deg, P.e) for P in K.places(p)]
continue
n = roots_mod(p, a, b, c, d)
if p < 500: # the gcd against brute roots
assert n == sum(1 for r in range(p)
if (((a * r - b) * r + c) * r - d) % p == 0)
out += {3: [(p, 1, 1)] * 3, 1: [(p, 1, 1), (p, 2, 1)],
0: [(p, 3, 1)]}[n]
return out
def quad_places(D):
out = []
for p in PRIMES:
k = kronecker(D, p)
out += {1: [(p, 1, 1)] * 2, -1: [(p, 2, 1)], 0: [(p, 1, 2)]}[k]
return out
LOGP = np.log(np.array(PRIMES, dtype=float))
def tail(b):
"""T(b) = zeta(b) prod_(p < X) (1 - p^-b)."""
return float(np.exp(log(zeta(b)) + np.log1p(-np.exp(-b * LOGP)).sum()))
class Race:
def __init__(self, places, n):
self.n = n
self.places = places
self.logN = np.array([f * log(p) for p, f, e in places])
self.two = [(f, e) for p, f, e in places if p == 2]
self.g = sum(1 for f, e in self.two if f == 1)
self.fermat = [p ** f for p, f, e in places if p != 2
and p ** f in FERMAT]
self.n3 = sum(1 for p, f, e in places if p == 3 and f == 1)
held = []
for p, f, e in places:
M = p ** f
if p in (2, 3):
continue
m = M - 1
while m % 2 == 0:
m //= 2
if (M - 1) % 4 == 0:
continue # 4 | M - 1: not | 2 3^oo
while m % 3 == 0:
m //= 3
if m == 1:
held.append(M)
self.held = held
def zeta_K(self, b, m):
return float(np.exp(-np.log1p(-np.exp(-b * self.logN)).sum())) \
* tail(b) ** m
def sigma2(self, b):
x = 2.0 ** -b
s = (1 - x) ** -(self.g - 1)
for N in self.fermat:
s *= 1 + N ** -b
return s
def sigma3(self, b):
"""The column at one place of norm 3."""
x = 2.0 ** -b
s = 1.0
for f, e in self.two:
if f == 1:
d = 3 if e == 1 else 2
s *= sum(x ** j for j in range(d + 1))
elif f == 2:
s *= 1 + 4.0 ** -b + 16.0 ** -b
s /= (1 - 3.0 ** -b) ** (self.n3 - 1)
for M in self.held:
s *= 1 + M ** -b
return s
def root(self, sigma, m):
lo, hi = 1.05, 4.0
f = lambda b: self.zeta_K(b, m) - sigma(b) - 1 # noqa: E731
assert f(lo) > 0 > f(hi)
while hi - lo > 1e-9:
mid = (lo + hi) / 2
if f(mid) > 0:
lo = mid
else:
hi = mid
return (lo + hi) / 2
def read(self):
out = {}
# F_2's lower end: sigma2 times the Fermat places past 65537
lo2 = lambda b: (self.sigma2(b) # noqa: E731
* exp(6 * 2.0 ** (-32 * b)))
out["F2"] = (self.root(lo2, 0), self.root(self.sigma2, 1),
self.root(self.sigma2, self.n))
# N=3's lower end: sigma3 times the bound on its places past X
hi3 = lambda b: self.sigma3(b) * exp(missed_mass(b)) # noqa: E731
out["N3"] = (self.root(hi3, 0), self.root(self.sigma3, 1),
self.root(self.sigma3, self.n))
return out
def missed_mass(b):
"""Three places at each norm 2 3^a + 1 past X: above the true
missed mass of the column of norm 3."""
s, a = 0.0, 0
while True:
M = 2 * 3 ** a + 1
if M > X:
s += 3 * M ** -b
if M ** -b < 1e-30:
return s
a += 1
def control():
print("CONTROL C: D = -8 through the race engine")
R = Race(quad_places(-8), 2)
r = R.read()
f2, n3 = r["F2"], r["N3"]
print(f" F2 [{f2[0]:.5f}, {f2[2]:.5f}] central {f2[1]:.5f} "
f"(fieldtop.py 1.63177)")
print(f" N3 [{n3[0]:.5f}, {n3[2]:.5f}] central {n3[1]:.5f} "
f"(fieldtop.py 1.49122)")
check("C each bracket holds fieldtop.py's point, the central reading "
"within 2e-4",
f2[0] <= 1.63177 <= f2[2] and n3[0] <= 1.49122 <= n3[2]
and abs(f2[1] - 1.63177) < 2e-4 and abs(n3[1] - 1.49122) < 2e-4)
def ptype(K, p):
"""The splitting type of p in K, in cubtop.py's codes."""
pl = tuple(sorted((P.deg, P.e) for P in K.places(p)))
return {((1, 1), (1, 1), (1, 1)): "s", ((1, 1), (2, 1)): "q",
((3, 1),): "i", ((1, 1), (1, 2)): "r", ((1, 3),): "t"}[pl]
def hunt(cap=3000):
"""(H) every complex cubic field with |d_K| <= cap whose
configuration fails cubtop.py's bracket: FIELDS, exactly."""
import cubtop
print(f"\nHUNT H: the complex cubic fields |d_K| <= {cap} in a "
f"failing configuration")
G = cubtop.Grid(1.0001 + cubtop.STEP * np.arange(29999))
top = cubtop.top_index(G, 3)
fields = enumerate_fields(cap)
found = {}
for form, dk in fields.items():
K = Field(form, dk)
t2 = ptype(K, 2)
if t2 not in "srqt":
continue
ty = "".join(ptype(K, p) for p in cubtop.REL)
bad, _, least = cubtop.sweep(
G, [((t2, ty), cubtop.cubic_profile(t2, ty))], top)
if bad:
found[(dk, form)] = ty
print(f" d_K {dk}: 2 {t2}, config {ty}, least ratio "
f"{least[0]:.6f}, failing {[b[1] for b in bad]}")
print(f" {len(fields)} fields read, {len(found)} in a failing "
f"configuration")
check("H the fields of a failing configuration are FIELDS",
set(found) == set(FIELDS.items()))
def main():
control()
hunt()
if not all(CHECKS):
print("a control failed: no verdict")
sys.exit(1)
print(f"\nmissed mass of the column of norm 3 past 10^6 at b = 1.4: "
f"{missed_mass(1.4):.2e}")
tally = {"LOSS": [], "TOP": [], "OPEN": []}
for dk, form in FIELDS.items():
R = Race(cubic_places(form, dk), 3)
r = R.read()
f2, n3 = r["F2"], r["N3"]
v = ("LOSS" if n3[0] > f2[2] else "TOP" if n3[2] < f2[0]
else "OPEN")
tally[v].append(dk)
if v == "LOSS":
K = Field(form, dk)
ok = K.class_group()
check(f"K d_K {dk}: the class group against the analytic "
f"class number formula", ok, f"h = {K.h}")
two = "".join(f"f{f}e{e} " for f, e in R.two)
print(f" d_K {dk}: 2 = {two}| F2 [{f2[0]:.5f}, {f2[2]:.5f}] | "
f"N=3 [{n3[0]:.5f}, {n3[2]:.5f}] | central gap "
f"{f2[1] - n3[1]:+.5f} | {v}")
print(f" Fermat places {sorted(R.fermat)}, norm-3 places "
f"{R.n3}, held {sorted(R.held)}")
print("RACE " + ", ".join(f"{k} {len(v)} {v}" for k, v in tally.items()))
check("R LOSS at exactly -1420 and -2764, OPEN at none",
sorted(tally["LOSS"]) == [-2764, -1420] and not tally["OPEN"])
print(f"{sum(CHECKS)} of {len(CHECKS)} checks pass")
sys.exit(0 if all(CHECKS) else 1)
if __name__ == "__main__":
main()