localtype.py
Cited on: PRINCIPAL.
To run it: save localtype.py, with module_law.py, principal.py beside it; needs numpy.
"""localtype.py -- does the local type at the small primes carry the
rise of PI_1 - PI_k0 with |D|?
QUESTION. Over the real fields with the class k0 (their fundamental
unit of norm +1), |D| <= 48000, d = PI_1 - PI_k0, the trivial class's corrected
prime-power count less k0's, rises with log D: +0.29 per unit at the cut
250 and +0.88 at 10^4 (weighted by h+/2, the odd characters). Is the rise
the population's local type at the small primes shifting with |D|, or
does it survive with the local type held fixed?
THE DERIVATION (before code).
(a) d splits exactly by the kind of p: d = d_I + d_R + d_S. An inert
p gives (p), of norm p^2, principal and totally positive, so every
power lies in the trivial class, weight 1/e at p^(2e) < x: at 250,
w(2) = 11/6, w(3) = 3/2, w(p) = 1 at 5, 7, 11, 13, and nothing from
17, so d_I at 250 is a function of the local type at 2..13 alone.
A ramified R has R^2 = (p) totally positive: its even powers go to
the trivial class, its odd powers to R's class, +-1 to d only when
that is 1 or k0. A split P, Pbar in c, c^-1 reach d only through
powers landing in 1 or k0.
(b) Unfiltered, the fundamental discriminants have 2 inert (D = 5 mod
8) at density 1/3, an odd p inert at p/(2(p + 1)) and ramified at
1/(p + 1); so the unfiltered mean of d_I at 250 is 11/18 + 9/16 +
5/12 + 7/16 + 11/24 + 13/28 = 2.950, the size of d itself.
(c) Every D with a negative prime discriminant factor (-4, -8, -q, q
= 3 mod 4) has k0. At small |D| the factors are small, so the
filter enriches the ramification of 2, 3, 7, 11 and lowers d_I; as
|D| grows a large factor can carry the sign and the enrichment
relaxes. The guess: the filtered d_I rises and carries the slope.
(d) The weight h+/2 grows with the 2-rank: both weightings print.
(e) The within slope s_w: weighted OLS of y on log D with an
intercept per stratum, the stratum (D/p) at p = 2, 3, 5, 7, 11,
13. Its denominator, the within variance of log D, prints beside
the total; its se is HC0 per field. It is linear in y, so it
splits over the parts, and s_w(d_I) = 0 exactly at 250.
PREDICTIONS, frozen before the engine.
P1 the parts sum to the counts' d per field to 1e-12 at every cut.
P2 s(d_I) >= 0.5 s(d) at 250, character-weighted.
P3 the share of fields with 3 ramified falls from the first |D| band
to the last, filtered; unfiltered it sits within 0.03 of 1/4 in
every band.
P4 the unfiltered band means of d_I sit within 0.15 of 2.950; the
filtered first band's below 2.80.
KILL LOCAL (printed: "within-stratum slope of d at 250"), character-
weighted: fires if s_w >= 0.5 s, the local type then not carrying
the rise. At 10^4 the strata do not fix d_I: printed, no verdict.
C1 read first: s(d) at 250 and 10^4 within 0.005 of +0.29 and +0.88,
and P1.
C2 read before the kill: s_w(d_I) at 250 is 0 to 1e-12; s_w(d_I + c
log D), c = s(d), reads c to 1e-9; s_w(d) with the strata shuffled
across fields (one seed) within 3 se of the pooled s.
ARMS: filtered (the fields with k0) and unfiltered (every fundamental
D to 48000, the local type only), the filter the ingredient
removed.
RESOLVED BEFORE THE RUN, nothing read: P3's shares and P4's means read
field-weighted on both arms, like for like; the kill and P2 stay
character-weighted, as frozen. ADDED BEFORE THE RUN: a check that the
inert part at 250 is (a)'s hand function of the local type, touching no
kill.
FINDINGS, the 10323 real fields with k0 to |D| <= 48000 (32585 odd
characters); unfiltered, every fundamental D to 48000.
F0 CONTROLS. P1 held: the parts sum to d to 1.14e-13, and the inert
part at 250 is the hand function of the local type exactly. C1:
s(d) +0.2899 at 250 and +0.8788 at 10^4, parity.py's. C2: s_w of
the inert part 5.5e-33; plus c log D it reads c to the printed
digit; with the strata shuffled +0.3074 +- 0.0409 against
+0.2899. 685 strata, 625 of two or more fields holding 10263; the
within sum of squares of log D 20526.8 of the total 22370.5.
F1 KILL LOCAL FIRED: within the strata of fixed local type at 2..13
the slope of d at 250 reads +0.2740 +- 0.0410 of the pooled
+0.2899 (field-weighted +0.3993 +- 0.0461). The local type at the
small primes does not carry the rise.
F2 P2 MISSED, the hand guess wrong: the split primes carry it. At 250
s(d) = inert +0.0243 +- 0.0190 + ramified +0.0076 +- 0.0164 +
split +0.2581 +- 0.0392 (within the strata: +0.0000, +0.0123 +-
0.0133, +0.2617 +- 0.0403); at 10^4 +0.8788 = +0.0541 + 0.0169 +
0.8078 +- 0.1804 (within: +0.9399 +- 0.1713, split +0.8862).
Field-weighted at 250, s(d) +0.4101 = inert +0.0219 + ramified
-0.0065 + split +0.3947.
F3 P3 and P4 HELD, the confound real but flat: the share with 3
ramified reads 0.383, 0.364, 0.355, 0.351, 0.347 by band (to
4000, 8000, 16000, 32000, 48000) filtered, 0.252 to 0.250
unfiltered; the inert part at 250 2.672, 2.729, 2.741, 2.743,
2.744 filtered against 2.962, 2.959, 2.945, 2.948, 2.954
unfiltered (2.950 by hand). d at 250, field-weighted, +2.140,
+3.269, +3.308, +3.474, +3.543: most of the rise is from the first
band to the second.
TIERS. (a) a property, (b) known densities; F1 to F3 observations over
this population. Why the split part rises is not derived.
RUN RECORD. A rehearsal at --d 3000 (C1 fails there by design, its
references being the full population's). Full: 7/7, 41.5 s; peak commit
55.8 MB under a memory guard (512 MB) at the first full run. P2, P3,
P4 and the kill stay prints, since a verdict is read, never asserted, a missed one is the
finding, and a check before the controls' gate would block the kill.
"""
import os
os.environ.setdefault("OPENBLAS_NUM_THREADS", "1")
import random
import sys
import time
from math import log
HERE = os.path.dirname(os.path.abspath(__file__))
sys.path.insert(0, HERE)
import numpy as np # noqa: E402
import principal as P # noqa: E402
from module_law import check, section, CHECKS # noqa: E402
D_MAX = 48000
if "--d" in sys.argv:
D_MAX = int(sys.argv[sys.argv.index("--d") + 1])
STRATA_P = (2, 3, 5, 7, 11, 13)
BANDS = ((0, 4000), (4000, 8000), (8000, 16000), (16000, 32000),
(32000, 10 ** 9))
J250, J10K = 0, len(P.CUTS) - 1
D_I_UNFILTERED = 11 / 18 + 9 / 16 + 5 / 12 + 7 / 16 + 11 / 24 + 13 / 28
def parts(G):
"""Per cut, d's inert, ramified and split parts (cumulative)."""
nc = len(P.CUTS)
out = np.zeros((3, nc)) # rows: I, R, S
def add(kind, n, c, wt):
j = P.bucket(n)
if j is None:
return
if c == G.e:
out[kind, j:] += wt
elif c == G.k0:
out[kind, j:] -= wt
for p in P.SMALL:
k = P.kron(G.D, p)
if k == 1:
c = G.place(p)
cb = G.inv(c)
e, pe = 1, p
while pe < P.TOP:
add(2, pe, G.power(c, e), 1 / e)
add(2, pe, G.power(cb, e), 1 / e)
e, pe = e + 1, pe * p
elif k == -1:
e, pe = 1, p * p
while pe < P.TOP:
add(0, pe, G.e, 1 / e)
e, pe = e + 1, pe * p * p
else:
r = G.place(p)
e, pe = 1, p
while pe < P.TOP:
add(1, pe, G.power(r, e), 1 / e)
e, pe = e + 1, pe * p
return out
def stratum(D):
return tuple(P.kron(D, p) for p in STRATA_P)
def wslope(x, y, w):
"""Weighted OLS slope of y on x and its field-level (HC0) se."""
xm = np.sum(w * x) / np.sum(w)
ym = np.sum(w * y) / np.sum(w)
sxx = np.sum(w * (x - xm) ** 2)
s = np.sum(w * (x - xm) * (y - ym)) / sxx
e = y - ym - s * (x - xm)
return s, np.sqrt(np.sum((w * (x - xm) * e) ** 2)) / sxx
def within(x, y, w, g):
"""The within-stratum slope, its se, and the within and total sums
of squares of x (weighted)."""
xd, yd = np.empty_like(x), np.empty_like(y)
for key in set(g):
m = np.array([gi == key for gi in g])
xd[m] = x[m] - np.sum(w[m] * x[m]) / np.sum(w[m])
yd[m] = y[m] - np.sum(w[m] * y[m]) / np.sum(w[m])
sxx = np.sum(w * xd ** 2)
s = np.sum(w * xd * yd) / sxx
e = yd - s * xd
se = np.sqrt(np.sum((w * xd * e) ** 2)) / sxx
tot = np.sum(w * (x - np.sum(w * x) / np.sum(w)) ** 2)
return s, se, sxx, tot
def read_filtered():
rows, worst = [], 0.0
for D in range(5, D_MAX + 1):
if not P.fundamental(D):
continue
G = P.Group(D)
if G.h < 2 or G.k0 == G.e:
continue
_, cor = P.counts(G)
cor = np.array(cor)
d = cor[:, G.e] - cor[:, G.k0]
pt = parts(G)
worst = max(worst, float(np.max(np.abs(pt.sum(axis=0) - d))))
rows.append((D, G.h // 2, d, pt))
return rows, worst
def band_of(D):
for i, (lo, hi) in enumerate(BANDS):
if lo < D <= hi:
return i
def d_inert_250(D):
wt = {2: 11 / 6, 3: 3 / 2, 5: 1, 7: 1, 11: 1, 13: 1}
return sum(v for p, v in wt.items() if P.kron(D, p) == -1)
def main():
t0 = time.time()
section("A THE FILTERED POPULATION AND THE SPLIT")
rows, worst = read_filtered()
n = len(rows)
Ds = np.array([r[0] for r in rows])
x = np.log(Ds.astype(float))
wc = np.array([r[1] for r in rows], float)
wf = np.ones(n)
g = [stratum(int(D)) for D in Ds]
print(f" {n} real fields with k0, |D| <= {D_MAX}, "
f"{int(wc.sum())} odd characters [{time.time() - t0:.0f} s]")
check("P1 the parts sum to d per field at every cut", worst < 1e-12,
f"worst {worst:.2e}")
dI250 = np.array([r[3][0, J250] for r in rows])
hand = np.array([d_inert_250(int(D)) for D in Ds])
check("the inert part at 250 is the hand function of the local type",
np.max(np.abs(dI250 - hand)) < 1e-12,
f"worst {np.max(np.abs(dI250 - hand)):.2e}")
pooled = {}
for j, ref in ((J250, 0.29), (J10K, 0.88)):
d = np.array([r[2][j] for r in rows])
s, se = wslope(x, d, wc)
pooled[j] = (s, se)
ps = [wslope(x, np.array([r[3][k, j] for r in rows]), wc)
for k in range(3)]
print(f" cut {P.CUTS[j]:5d}, character-weighted: s(d) {s:+.4f}"
f" +- {se:.4f} = inert {ps[0][0]:+.4f} +- {ps[0][1]:.4f}"
f" + ramified {ps[1][0]:+.4f} +- {ps[1][1]:.4f}"
f" + split {ps[2][0]:+.4f} +- {ps[2][1]:.4f}")
sf = [wslope(x, np.array([r[3][k, j] for r in rows]), wf)[0]
for k in range(3)]
print(f" field-weighted: s(d) "
f"{wslope(x, d, wf)[0]:+.4f} = inert {sf[0]:+.4f}"
f" + ramified {sf[1]:+.4f} + split {sf[2]:+.4f}")
check(f"C1 s(d) at {P.CUTS[j]} within 0.005 of parity.py's "
f"{ref:+.2f}", abs(s - ref) <= 0.005, f"{s:+.4f}")
if j == J250:
p2 = ps[0][0] >= 0.5 * s
print(f" P2 the inert part carries at least half at 250: "
f"{'HELD' if p2 else 'MISSED'} ({ps[0][0]:+.4f} of"
f" {s:+.4f})")
section("B THE SHARES AND MEANS BY |D| BAND, FIELD-WEIGHTED")
un = [D for D in range(5, D_MAX + 1) if P.fundamental(D)]
ub = [band_of(D) for D in un]
fb = [band_of(int(D)) for D in Ds]
print(" band filtered: n, 3 ramified, d_I(250), d(250)"
" | unfiltered: n, 3 ramified, d_I(250)")
sh_f, sh_u, mi_u, mi_f = [], [], [], []
for i, (lo, hi) in enumerate(BANDS):
fi = [k for k in range(n) if fb[k] == i]
ui = [D for D, b in zip(un, ub) if b == i]
r3f = np.mean([Ds[k] % 3 == 0 for k in fi])
r3u = np.mean([D % 3 == 0 for D in ui])
dif = np.mean(dI250[fi])
df = np.mean([rows[k][2][J250] for k in fi])
diu = np.mean([d_inert_250(D) for D in ui])
sh_f.append(r3f)
sh_u.append(r3u)
mi_u.append(diu)
mi_f.append(dif)
print(f" {lo:5d}-{min(hi, D_MAX):5d} {len(fi):5d} {r3f:.3f}"
f" {dif:.3f} {df:+.3f} | {len(ui):5d} {r3u:.3f}"
f" {diu:.3f}")
p3 = sh_f[-1] < sh_f[0] and all(abs(v - 0.25) <= 0.03 for v in sh_u)
print(f" P3 the filtered 3-ramified share falls, the unfiltered near"
f" 1/4: {'HELD' if p3 else 'MISSED'}")
p4 = (all(abs(v - D_I_UNFILTERED) <= 0.15 for v in mi_u)
and mi_f[0] < 2.80)
print(f" P4 the unfiltered d_I within 0.15 of {D_I_UNFILTERED:.3f},"
f" the filtered first band below 2.80: "
f"{'HELD' if p4 else 'MISSED'}")
section("C THE WITHIN-STRATUM SLOPE")
sizes = {}
for key in g:
sizes[key] = sizes.get(key, 0) + 1
multi = sum(1 for v in sizes.values() if v >= 2)
inmulti = sum(v for v in sizes.values() if v >= 2)
print(f" {len(sizes)} strata, {multi} with two or more fields,"
f" holding {inmulti} of {n}")
s250, se250 = pooled[J250]
z0 = within(x, dI250, wc, g)
zc = within(x, dI250 + s250 * x, wc, g)
rnd = random.Random(1)
gs = list(g)
rnd.shuffle(gs)
d250 = np.array([r[2][J250] for r in rows])
zs = within(x, d250, wc, gs)
print(f" within/total sum of squares of log D: {z0[2]:.1f} / "
f"{z0[3]:.1f}")
check("C2 s_w of the inert part at 250 is 0", abs(z0[0]) < 1e-12,
f"{z0[0]:.2e}")
check("C2 s_w of the inert part plus c log D reads c",
abs(zc[0] - s250) < 1e-9, f"{zc[0]:+.6f} against {s250:+.6f}")
check("C2 s_w with the strata shuffled sits within 3 se of s",
abs(zs[0] - s250) <= 3 * zs[1],
f"{zs[0]:+.4f} +- {zs[1]:.4f} against {s250:+.4f}")
if not all(CHECKS):
print("CONTROL FAILED: the kill is not read.")
raise SystemExit(1)
for j in (J250, J10K):
d = np.array([r[2][j] for r in rows])
sw = within(x, d, wc, g)
pr = [within(x, np.array([r[3][k, j] for r in rows]), wc, g)
for k in range(3)]
print(f" within-stratum slope of d at {P.CUTS[j]}: {sw[0]:+.4f}"
f" +- {sw[1]:.4f} of pooled {pooled[j][0]:+.4f}; inert"
f" {pr[0][0]:+.4f}, ramified {pr[1][0]:+.4f} +-"
f" {pr[1][1]:.4f}, split {pr[2][0]:+.4f} +- {pr[2][1]:.4f}")
swf = within(x, d, wf, g)
print(f" field-weighted: {swf[0]:+.4f} +- {swf[1]:.4f}")
if j == J250:
fired = sw[0] >= 0.5 * pooled[j][0]
print(f" KILL LOCAL (s_w >= half the pooled): "
f"{'FIRED' if fired else 'MISSED'}")
print(f"\n{sum(CHECKS)}/{len(CHECKS)} checks PASS,"
f" {time.time() - t0:.1f} s")
raise SystemExit(0 if all(CHECKS) else 1)
if __name__ == "__main__":
main()