"""principal.py -- where a quadratic field's first principal split prime sits, and why principal split primes run short of their share at small p. QUESTION. K = Q(sqrt D), D a fundamental discriminant, O its ring. An odd prime p with chi_D(p) = +1 splits as P Pbar, and P is PRINCIPAL when P = (alpha); over a real field "principal" means principal as an ideal (the wide sense) unless the narrow group is named. L_1(D) is the least odd prime p with chi_D(p) = +1 and P principal. Two questions. (a) How far down can L_1 sit? Over an imaginary field is there a floor set by |D|, and does a real field have one? (b) By Chebotarev in the Hilbert class field the principal split primes have density 1/(2h) among all primes, 1/h of the split ones. At a finite cut they run SHORT, the more so as h grows. What is the shortfall? THE ARGUMENT (written before the engine). (1) THE TEST. The principal form of discriminant D is x^2 + b0 x y + c0 y^2, b0 = D mod 2, c0 = (b0 - D)/4, and N(x + y w) is its value at (x, y), w = (b0 + sqrt D)/2. A split P over p is principal iff p or -p is a value of it (the conjugate place is then principal too), iff 4p = |u^2 - D y^2| with u = 2x + b0 y. Over an imaginary field the form is definite and the sign is +. (2) THE FLOOR. Let D < 0 and P over p principal: 4p = u^2 + |D| y^2. y = 0 makes 4p = u^2, so p is a square, which it is not; hence y != 0 and 4p >= |D|. At odd D, u = b0 y = y mod 2 is odd with y, so 4p >= 1 + |D|. So L_1 >= |D|/4, and >= (|D| + 1)/4 at odd D. EQUALITY. At odd D the value at (0, 1) is c0 = (1 - D)/4 = (|D| + 1)/4, and D = 1 - 4 c0 gives gcd(c0, D) = 1; so if c0 is an odd prime it is unramified, represented, hence split and principal, and L_1 = c0. At D = 0 mod 4 the floor |D|/4 divides D, so a prime there is ramified and never counted, and 4p = |D| is impossible for a split p; the floor is not attained. Conversely L_1 is an odd prime, so L_1 = (|D| + 1)/4 needs that value to be one. So L_1 MEETS THE FLOOR EXACTLY WHEN D IS ODD AND (|D| + 1)/4 IS AN ODD PRIME. The floor never consults h. COROLLARY. No imaginary field with |D| > 4B has L_1 <= B. (3) NO FLOOR OVER A REAL FIELD. Take D = m^2 - 12 or m^2 + 12 with m odd and 3 not dividing m. Then D = 1 mod 4 (m^2 = 1 mod 8, 12 = 4 mod 8, so D = 5 mod 8 either way), and with w = (1 + sqrt D)/2 and x = (m - 1)/2, N(x + w) = x^2 + x - (D - 1)/4 = (m^2 - D)/4 = +-3. D = m^2 mod 3 is a nonzero square, so chi_D(3) = +1, and x + w generates a place over 3. So L_1(D) = 3 whenever such a D is fundamental, that is squarefree; m^2 -+ 12 is an irreducible quadratic in m, so its squarefree values are infinite in number (Estermann 1931), and L_1 = 3 at unboundedly large D. [CORRECTION, made at the audit and left in place because the slate is a frozen record: Estermann 1931 is the case z^2 + k. With m = 6n +- 1 the family is an irreducible quadratic in n with no fixed square divisor, and its squarefree values have positive density by Ricci, Rend. Circ. Mat. Palermo 57 (1933).] A floor in terms of D is impossible, and the reason is (1)'s sign: an indefinite form takes small values far out. (4) THE SHORTFALL (the explicit formula, known in shape). Chebotarev in the Hilbert class field H/K, whose group is Cl, is a statement about the count PI_C(x) of all prime-ideal POWERS P^k of K with N(P^k) < x and class [P^k] = C, each weighted 1/k: its main term is Li(x)/h for every C, and its remainder oscillates with no drift, so over a POPULATION of fields the level PI_C/(PI/h) averages toward 1. A count of split places omits three families, and each lands on known classes: the inert prime ideals (q), norm q^2, are principal: weight 1 per inert q < sqrt x, ALL on the trivial class; the squares of split places P^2, norm p^2: weight 1/2 at 2[P] and 2[Pbar], so on the squares 2Cl; the ramified prime ideals R, R^2 = (r) principal, so [R] in Cl[2]. Higher powers are smaller still. So the split-place count of a class is short by the weight its class receives, and the three families give a class weight only if it is TRIVIAL, a SQUARE, or 2-TORSION. [CORRECTION, made at the audit, the slate left as frozen: the odd powers k >= 3 of split places land on k[P], any class, at weight 1/k over p < x^(1/k). The engine counts them; "only if" holds for the squares and the inert and ramified primes, and G loses the sliver of odd-power weight.] A class that is none of these (a non-square outside Cl[2]) loses nothing and reads LONG after the renormalisation by the realised total. Sizes: over p < x the split places number about pi(x), so each class nominally holds pi(x)/h; the trivial class loses about pi(sqrt x)/2 inert primes plus the share |Cl[2]|/h of the pi(sqrt x)/2 weight of split squares, a relative deficit of (h + |Cl[2]|)/2 * pi(sqrt x)/pi(x), growing linearly in h. Over a real field the same holds in the NARROW group Cl+: (q) = q O has a generator of positive norm, and P Pbar = (p) and R^2 = (r) likewise. (5) THE CLASSICAL INSTANCE. Over Q with the residue classes mod an odd prime q: a square l^2 < x is a quadratic residue, so the residue classes' prime count is short of the non-residues' by about pi(sqrt x)/2 and the non-residue lead, (N - R)/((N + R)/2), is about pi(sqrt x)/pi(x); counting prime POWERS with weight 1/k puts it back. (6) THE STATISTIC. Per field and cut, X_C is the split-place count (RAW) or PI_C (CORRECTED); the level of a set S of classes is the mean of X_C over S divided by the field's total over h. The total is at least pi(250) = 53, so the denominator never nears zero. A field's levels average to 1 over its h classes, so a cell is read against its siblings, never alone. Cells: T the trivial class; A the nontrivial classes of Cl[2]; S the nontrivial squares outside Cl[2]; G the classes in neither 2Cl nor Cl[2], which (4) says receive no weight. Over a real field K0, the narrow class k0 of a principal ideal whose generators all have negative norm, is read alone and A excludes it. Pooled figures are the mean over fields with its standard error. THE ENGINE. Classes are binary quadratic forms (a, b, c), b^2 - 4ac = D: reduced forms over D < 0; over D > 0 the reduced forms fall into cycles under the step rho, and a proper (narrow) class is a cycle. A place over p at a root b of b^2 = D mod 4p is the form (p, b, (b^2 - D)/(4p)), its conjugate (p, -b, .). A power P^k is read by composition; the power form (p^2, b_2, .), b_2 the lift of b mod 4p^2, is built only to check the square (K3). Composition is Dirichlet's. PREDICTIONS, fixed before the engine. K CONTROLS, read before anything else. K1 h(D) = 1, 1, 3, 5, 7, 9, 4, 4, 8 at D = -3, -4, -23, -47, -71, -199, -56, -84, -420, and h+(D) = 1, 1, 2, 2, 4, 3 at D = 5, 8, 12, 40, 60, 229; nine imaginary D with h = 1 in the imaginary population's range. K2 GENUS: |Cl[2]| = 2^(t - 1) (Cl+ over a real field), t the number of primes dividing D, at every field of both populations. K3 COMPOSITION: at every field and every split p < 100, the square of [P] by composition equals the class of the power form P^2, and [P] + [Pbar] is trivial. [As built, the second clause composes [P] with the inverse form (a, -b, c) of its reduced representative, so it checks composition against inversion, not the conjugate place's own form.] K4 WEIGHTS: at every field and cut, the summed PI_C over classes equals the prime-power weight computed from chi_D alone. K5 THE TEST: over D < 0 with |D| <= 400 and p < 200, the class engine's "P trivial" agrees with a box search for 4p = u^2 + |D| y^2. K6 THE CLASSICAL RACE (5): over the odd primes q < 400 at x = 10^6, the raw non-residue lead pooled over q is positive at 3 standard errors or more and within a factor 2 of pi(sqrt x)/pi(x); the prime-power count's lead is within 2 standard errors of 0. F THE FLOOR. Over every imaginary fundamental D, -20000 <= D < 0: L_1 >= |D|/4 at every field, and L_1 = (|D| + 1)/4 at exactly the odd D where that value is an odd prime; no even D meets |D|/4. R NO FLOOR. Every squarefree D = m^2 +- 12 < 10^6, m odd, 3 not dividing m, has L_1 = 3, read off the class engine (the place over 3 wide-principal) and off the explicit element of norm +-3. P THE SHORTFALL, imaginary population: every fundamental D with -4000 <= D < 0 and h >= 2; cuts 250, 400, 630, 1000, 2500, 10000. P1 G, pooled: raw above 1 by 10 standard errors or more at every cut; corrected within 0.01 of 1 at 10000, and the correction removes at least 70% of the raw excess at every cut from 630 up. P2 At 10000 every (h, cell) with at least 30 fields reads, corrected, within 0.03 of 1 or within 3 standard errors of it. P3 T at 10000: the raw deficit 1 - level over the (4) prediction (h + |Cl[2]|)/2 * pi(sqrt x)/pi(x), pooled, lies in [0.7, 1.3]. Real population: every fundamental D, 5 <= D <= 16000, with h+ >= 2, in the narrow group. P4 G, pooled, at 10000: corrected within 0.01 of 1 and at least 90% of the raw excess removed. At cuts below 10000 the real cells are printed with no prediction. THE TRANSPLANTS. P1's 70%, P2's 0.03 and P4's 90% are fixed from a neighbouring reading at these populations and cuts, not derived; the classes in G and the corrections are (4)'s. [SETTLED ON A LATER READ, the argument and slate above left as frozen. (4)'s deficit is against the Chebotarev share Li(x)/h; a level divides by the realised total, short of Li(x) by about pi(sqrt x), so T's level falls short by (4)'s deficit less pi(sqrt x)/pi(x), ((h + |Cl[2]|)/2 - 1) pi(sqrt x)/pi(x) to leading order, printed beside P3. (6)'s raw total is twice the split primes below the cut, about pi(cut), with no floor of 53. K3 and K5 read odd p. These hold by construction and are not checks: K4 (the weight sum grows in the same branches as the classes), F's bound and even case (the search hits only where |D| y^2 <= 4p with y >= 1, y = 0 never hitting as 4p is no square, and chi_D filters the even case), and R's norm +-3 and chi_D(3) (algebra). F's equality holds so too and is kept as a printed count: at odd D, below (|D| + 1)/4 every 4p < |D| leaves only y = 0, and at it (u, y) = (1, 1) hits, with chi_D = +1 as D = 1 mod it. R also asserts L_1 < 10^4 over the real population, a check added after the slate. Every section runs whatever the K lines print: a K line off is read as the kill, not enforced by the code.] KILLS, as printed observables: any K line off (nothing below is read); an F violation or an equality off its criterion; an R member whose place over 3 is not wide-principal; P1's G corrected level off 1 by more than 0.01 at 10000; a P2 cell outside both bands; P3's ratio outside [0.7, 1.3]; P4's corrected G off 1 by more than 0.01. FINDINGS. Every prediction landed: 19/19 checks PASS at the first run, no kill fired (17/17 since, as the RUN RECORD says). K K1 all fifteen class numbers as named, h = 1 at exactly -3, -4, -7, -8, -11, -19, -43, -67, -163. K2 |Cl[2]| = 2^(t - 1) at all 1217 imaginary and all 4865 real fields. K3 and K5 0 off. K6 over 77 moduli at 10^6 the raw non-residue lead reads 0.00180 +- 0.00031 (5.9 se) against pi(sqrt x)/pi(x) = 0.00214, and the prime-power lead -0.00042 +- 0.00031 (-1.4 se). F The floor is met at 499 of 6079 fields, every one odd with (|D| + 1)/4 an odd prime, 0 off. The median L_1 runs 181, 587, 1163, 2269, 4483 over the |D| bands (0, 1000], (1000, 2500], ..., (10000, 20000], each inside the floor's range over its band. R 627 squarefree members up to 994021 (13 and 37 counted once each, 1 + 12 = 25 - 12 and 25 + 12 = 49 - 12), each with its place over 3 wide-principal, 0 off. Over all 4865 real fields to 16000, L_1 < 10000 at every one, and the median L_1 runs 7, 11, 11, 11, 11 over the D bands (0, 1000], (1000, 2500], ..., (10000, 16000]: 11 from the second band on. P Imaginary, 1208 fields; G on 880: cut 250 400 630 1000 2500 10000 raw 1.1862 1.1453 1.1060 1.0885 1.0567 1.0233 corr. 1.0128 1.0018 0.9960 1.0002 1.0032 0.9991 (se 0.0036 to 0.0005 raw), 42.7 se or more above 1 raw at every cut, 93% to 104% of the excess removed. At 10000 the 39 readable (h, cell) read 0.9715 to 1.0195 corrected, h = 4 to 28, T pooled over the 1208 fields 0.9873; raw T runs 0.950 at h = 4 down to 0.748 at h = 28. T's raw deficit at 10000 is 0.2080 against the predicted 0.2245 (ratio 0.926), and against (4)'s deficit less pi(sqrt x)/pi(x), 0.2041 (ratio 1.019). At low cuts T is the hard zero (raw 0.097 at 250: a field's principal places start at |D|/4) and its corrected level there is no prediction. Real narrow, 4107 fields; G on 1069: raw 1.2103 1.1716 1.1301 1.1043 1.0590 1.0255 corr. 1.0344 1.0241 1.0166 1.0135 1.0048 1.0009 84% to 97% removed. The rest is not uniform: k0 reads 0.824 corrected at 250 and 0.988 at 10000, T 1.068 and 1.005, the other cells within 0.035 of 1 at 250 (G 1.0344) and 0.005 at 10000. Tiers: (2) a theorem, the equality criterion with it; (3) a theorem, its infinitude by Ricci; (4) a property as to where the omitted weight lands, and a pattern as to the population means (imaginary whole, real to 84-97% with the k0 remainder left open); (5) the classical instance, read as the control. RUN RECORD. 19/19, 10.1 s wall, peak working set 58.6 MB under a memory guard. Later the checks that held by construction were cut (F's bound and even case, R's norm and chi_D(3), the two K4 weight sums), K1 asserts the nine h = 1 fields, R asserts L_1 < 10^4 over the real population, and the renormalised T prediction, the G ratios and K6's leads in standard errors are printed: 17/17, 10.3 s, peak commit 54 MB under a memory guard. """ import os os.environ.setdefault("OPENBLAS_NUM_THREADS", "1") import time from array import array from math import isqrt from statistics import mean, stdev from module_law import check, section, CHECKS RACE_X = 10 ** 6 # the classical race's cut RACE_Q = 400 # its moduli: the odd primes below this FLOOR_D = 20000 # the floor is read over -FLOOR_D <= D < 0 FAMILY_D = 10 ** 6 # the real family m^2 +- 12 is read below this IMAG_D = 4000 # the imaginary population: -IMAG_D <= D < 0 REAL_D = 16000 # the real population: 5 <= D <= REAL_D CUTS = (250, 400, 630, 1000, 2500, 10000) TOP = CUTS[-1] MIN_CELL = 30 # an (h, cell) is read with at least this many fields def primes_below(n): s = bytearray([1]) * n s[0:2] = b"\x00\x00" for i in range(2, isqrt(n - 1) + 1): if s[i]: s[i * i::i] = bytearray(len(range(i * i, n, i))) return [i for i in range(n) if s[i]] PRIMES = primes_below(RACE_X) SMALL = [p for p in PRIMES if p < TOP] SET_PRIMES = set(PRIMES) def pi_below(x): return sum(1 for p in SMALL if p < x) if x <= TOP else \ sum(1 for p in PRIMES if p < x) def squarefree(n): n = abs(n) for p in PRIMES: if p * p > n: return True if n % (p * p) == 0: return False return True def fundamental(D): if D in (0, 1): return False if D % 4 == 1: return squarefree(D) if D % 4 == 0: m = D // 4 return m % 4 in (2, 3) and squarefree(m) return False def omega(D): n, t = abs(D), 0 for p in PRIMES: if p * p > n: break if n % p == 0: t += 1 while n % p == 0: n //= p return t + (n > 1) def kron(D, p): if p == 2: return 0 if D % 2 == 0 else (1 if D % 8 == 1 else -1) r = pow(D % p, (p - 1) // 2, p) return 0 if r == 0 else (1 if r == 1 else -1) ROOT = {} for _p in SMALL[1:]: _t = array("i", [-1]) * _p for _r in range((_p + 1) // 2): _t[_r * _r % _p] = _r ROOT[_p] = _t def place_b(D, p): """b with b = D mod 2 and b^2 = D mod 4p: the place (p, b, .).""" if p == 2: return next(b for b in range(4) if (b * b - D) % 8 == 0) r = ROOT[p][D % p] if p in ROOT else \ next(r for r in range(p) if (r * r - D) % p == 0) return r + p if (r - D) % 2 else r def form(D, a, b): return (a, b, (b * b - D) // (4 * a)) # --------------------------------------------------------------------------- # classes as forms def xgcd(a, b): x0, x1, y0, y1 = 1, 0, 0, 1 while b: q, a, b = a // b, b, a % b x0, x1 = x1, x0 - q * x1 y0, y1 = y1, y0 - q * y1 return a, x0, y0 def compose(f1, f2, D): """Dirichlet composition of two primitive forms with a > 0.""" (a1, b1, c1), (a2, b2, c2) = f1, f2 if a1 > a2: (a1, b1, c1), (a2, b2, c2) = (a2, b2, c2), (a1, b1, c1) s = (b1 + b2) // 2 n = b2 - s if a2 % a1 == 0: y1, d = 0, a1 else: d, y1, _ = xgcd(a2, a1) if s % d == 0: y2, x2, d1 = -1, 0, d else: d1, x2, v = xgcd(s, d) y2 = -v v1, v2 = a1 // d1, a2 // d1 r = (y1 * y2 * n - x2 * c2) % v1 b3 = b2 + 2 * v2 * r return form(D, v1 * v2, b3) def reduce_imag(f): a, b, c = f D = b * b - 4 * a * c while True: k = (b + a - 1) // (2 * a) b -= 2 * a * k c = (b * b - D) // (4 * a) if a > c: a, b, c = c, -b, a continue if a == c and b < 0: b = -b return a, b, c class Real: """Reduction and the cycle step for one D > 0.""" def __init__(self, D): self.D, self.sq = D, isqrt(D) def reduced(self, f): a, b, _ = f s = self.sq return 0 < b <= s and 2 * abs(a) + b > s and 2 * abs(a) - b <= s def rho(self, f): _, b, c = f cc, s = abs(c), self.sq if cc > s: r = (-b) % (2 * cc) if r > cc: r -= 2 * cc else: r = s - (s + b) % (2 * cc) return form(self.D, c, r) def reduce(self, f): while not self.reduced(f): f = self.rho(f) return f def cycle(self, f): f = self.reduce(f) out = [f] g = self.rho(f) while g != f: out.append(g) g = self.rho(g) return out class Group: """The class group of D < 0, or the narrow class group of D > 0: a class is an index, `of` maps a form to it, `mul` composes.""" def __init__(self, D): self.D = D self.b0 = D % 2 self.one_form = form(D, 1, self.b0) self.idx, self.rep = {}, [] if D < 0: a = 1 while 3 * a * a <= -D: for b in range(-a + 1, a + 1): if (b * b - D) % (4 * a) == 0: f = form(D, a, b) if f[2] >= a and not (f[2] == a and b < 0): self.idx[f] = len(self.rep) self.rep.append(f) a += 1 self.red = reduce_imag else: R = Real(D) self.red, s = R.reduce, R.sq for b in range(1, s + 1): if (b - D) % 2: continue m = (D - b * b) // 4 for a in range((s - b) // 2 + 1, (s + b) // 2 + 1): if m % a: continue for f in (form(D, a, b), form(D, -a, b)): if f in self.idx: continue cyc = R.cycle(f) i = len(self.rep) for g in cyc: self.idx[g] = i self.rep.append(next(g for g in cyc if g[0] > 0)) self.h = len(self.rep) self.e = self.of(self.one_form) self._mul = {} self.dbl = [self.mul(i, i) for i in range(self.h)] self.two = {i for i in range(self.h) if self.dbl[i] == self.e} self.squares = set(self.dbl) self.k0 = self.of(form(D, -1, self.b0)) if D > 0 else self.e def of(self, f): return self.idx[self.red(f)] def mul(self, i, j): key = (i, j) if i <= j else (j, i) if key not in self._mul: self._mul[key] = self.of(compose(self.rep[i], self.rep[j], self.D)) return self._mul[key] def inv(self, i): a, b, c = self.rep[i] return self.of((a, -b, c)) def power(self, i, k): out = self.e for _ in range(k): out = self.mul(out, i) return out def place(self, p): return self.of(form(self.D, p, place_b(self.D, p))) # --------------------------------------------------------------------------- # the controls def section_controls(): section("K CONTROLS") want = {-3: 1, -4: 1, -23: 3, -47: 5, -71: 7, -199: 9, -56: 4, -84: 4, -420: 8, 5: 1, 8: 1, 12: 2, 40: 2, 60: 4, 229: 3} got = {D: Group(D).h for D in want} ones = [D for D in range(-IMAG_D, 0) if fundamental(D) and Group(D).h == 1] check("K1 class numbers at fifteen named fields; h = 1 at exactly the " "nine imaginary D", got == want and ones == [-163, -67, -43, -19, -11, -8, -7, -4, -3], f"off: {[(D, got[D]) for D in want if got[D] != want[D]]}, " f"h = 1 at {ones}") off = 0 for D in range(-400, 0): if not fundamental(D): continue G = Group(D) b0, c0 = D % 2, (D % 2 - D) // 4 for p in SMALL[1:]: if p >= 200: break if kron(D, p) != 1: continue box = any(x * x + b0 * x * y + c0 * y * y == p for x in range(-30, 31) for y in range(-30, 31)) off += box != (G.place(p) == G.e) check("K5 'P trivial' agrees with a box search, |D| <= 400, odd p < 200", off == 0, f"{off} off") def section_race(): section("K6 THE CLASSICAL RACE: residues mod q against non-residues") x = RACE_X raw, cor = [], [] for q in PRIMES[1:]: if q >= RACE_Q: break qr = {r * r % q for r in range(1, q)} N = R = 0 Nw = Rw = 0.0 for p in PRIMES: if p == q: continue if p % q in qr: R += 1 else: N += 1 k, pk = 2, p * p while pk < x: if pk % q in qr: Rw += 1 / k else: Nw += 1 / k k, pk = k + 1, pk * p raw.append((N - R) / ((N + R) / 2)) cor.append((N + Nw - R - Rw) / ((N + Nw + R + Rw) / 2)) pred = pi_below(isqrt(x) + 1) / len(PRIMES) mr, sr = mean(raw), stdev(raw) / len(raw) ** 0.5 mc, sc = mean(cor), stdev(cor) / len(cor) ** 0.5 print(f" {len(raw)} moduli, x = {x}: raw lead {mr:.5f} +- {sr:.5f} " f"({mr / sr:.1f} se), prediction pi(sqrt x)/pi(x) = {pred:.5f}; " f"prime-power lead {mc:.5f} +- {sc:.5f} ({mc / sc:.1f} se)") check("K6 raw lead >= 3 se and within a factor 2 of the prediction", mr >= 3 * sr and pred / 2 <= mr <= 2 * pred) check("K6 the prime-power lead is within 2 se of 0", abs(mc) <= 2 * sc) # --------------------------------------------------------------------------- # the floor, and its absence over a real field def least_principal_imag(D): """L_1 by the completed square, chi_D(p) = +1 checked beside it (a ramified p can be represented), searched from p = 3 up.""" b0, aD = D % 2, -D for p in PRIMES[1:]: y, hit = 0, False while aD * y * y <= 4 * p: t = 4 * p - aD * y * y u = isqrt(t) if u * u == t and (u - b0 * y) % 2 == 0: hit = True break y += 1 if hit and kron(D, p) == 1: return p return None def section_floor(): section("F THE FLOOR: L_1 >= |D|/4 over an imaginary field") fields = [D for D in range(-FLOOR_D, 0) if fundamental(D)] eq_off, tight = 0, 0 bands = {} for D in fields: L = least_principal_imag(D) prime_c0 = D % 2 == 1 and (1 - D) // 4 in SET_PRIMES and \ (1 - D) // 4 > 2 at = 4 * L == 1 - D eq_off += at != prime_c0 tight += at band = next(B for B in (1000, 2500, 5000, 10000, 20000) if -D <= B) bands.setdefault(band, []).append(L) check("F L_1 = (|D| + 1)/4 exactly when D is odd and it is an odd prime", eq_off == 0, f"{len(fields)} fields, {tight} meet the floor, " f"{eq_off} off") lo = 0 for B, Ls in sorted(bands.items()): Ls.sort() print(f" |D| in ({lo:5d}, {B:5d}]: {len(Ls):4d} fields, median L_1 " f"{Ls[len(Ls) // 2]:6d}, the floor |D|/4 runs {lo // 4} to " f"{B // 4}") lo = B def least_principal_real(G): for p in SMALL[1:]: if kron(G.D, p) == 1 and G.place(p) in (G.e, G.k0): return p return None def section_family(): section("R NO FLOOR: L_1 = 3 at every squarefree m^2 +- 12") n, off, top, seen = 0, 0, 0, set() m = -1 while m * m - 12 < FAMILY_D: m += 2 if m % 3 == 0: continue for D in (m * m - 12, m * m + 12): if not 0 < D < FAMILY_D or not squarefree(D) or D in seen: continue seen.add(D) n += 1 top = max(top, D) R = Real(D) cyc = set(R.cycle(form(D, 3, place_b(D, 3)))) wide = R.reduce(form(D, 1, 1)) in cyc or \ R.reduce(form(D, -1, 1)) in cyc off += not wide check("R every member has a wide-principal place over 3", off == 0 and n > 0, f"{n} fundamental D up to {top}, {off} off") # --------------------------------------------------------------------------- # the shortfall def bucket(n): """The first cut index with n below the cut, or None.""" for j, x in enumerate(CUTS): if n < x: return j return None def counts(G): """Per cut, the split-place count and the prime-power count of every class.""" h, nc = G.h, len(CUTS) raw = [[0.0] * h for _ in range(nc)] ext = [[0.0] * h for _ in range(nc)] for p in SMALL: k = kron(G.D, p) if k == 1: c = G.place(p) cb = G.inv(c) j = bucket(p) raw[j][c] += 1 raw[j][cb] += 1 e, pe = 2, p * p while pe < TOP: j = bucket(pe) ext[j][G.power(c, e)] += 1 / e ext[j][G.power(cb, e)] += 1 / e e, pe = e + 1, pe * p elif k == -1: e, pe = 1, p * p while pe < TOP: j = bucket(pe) ext[j][G.e] += 1 / e e, pe = e + 1, pe * p * p else: r = G.place(p) e, pe = 1, p while pe < TOP: j = bucket(pe) ext[j][G.power(r, e)] += 1 / e e, pe = e + 1, pe * p for j in range(1, nc): for c in range(h): raw[j][c] += raw[j - 1][c] ext[j][c] += ext[j - 1][c] cor = [[raw[j][c] + ext[j][c] for c in range(h)] for j in range(nc)] return raw, cor def cells(G): out = {"T": [G.e]} if G.D > 0 and G.k0 != G.e: out["K0"] = [G.k0] out["A"] = [c for c in G.two if c not in (G.e, G.k0)] out["S"] = [c for c in G.squares if c not in G.two] out["G"] = [c for c in range(G.h) if c not in G.squares and c not in G.two] return {k: v for k, v in out.items() if v} def level(X, S, h): return sum(X[c] for c in S) / len(S) / (sum(X) / h) def pooled(v): return mean(v), (stdev(v) / len(v) ** 0.5 if len(v) > 1 else 0.0) def population(sign): """Every field of the population: its group, cells and levels. The K2 and K3 controls run field by field.""" lo, hi = ((-IMAG_D, 0) if sign < 0 else (5, REAL_D + 1)) rows, genus_off, comp_off, seen = [], 0, 0, 0 for D in range(lo, hi): if not fundamental(D): continue G = Group(D) genus_off += len(G.two) != 2 ** (omega(D) - 1) for p in SMALL: if p >= 100: break if p > 2 and kron(D, p) == 1: b = place_b(D, p) b2 = next(b + 2 * p * t for t in range(p) if ((b + 2 * p * t) ** 2 - D) % (4 * p * p) == 0) c = G.place(p) comp_off += G.dbl[c] != G.of(form(D, p * p, b2)) comp_off += G.mul(c, G.inv(c)) != G.e seen += 1 if G.h < 2: continue raw, cor = counts(G) cs = cells(G) rows.append({ "D": D, "h": G.h, "two": len(G.two), "raw": [{k: level(raw[j], S, G.h) for k, S in cs.items()} for j in range(len(CUTS))], "cor": [{k: level(cor[j], S, G.h) for k, S in cs.items()} for j in range(len(CUTS))], }) name = "imaginary" if sign < 0 else "real narrow" check(f"K2 genus: |Cl[2]| = 2^(t - 1) at every {name} field", genus_off == 0, f"{seen} fields, {genus_off} off") check(f"K3 [P]^2 by composition = the power form; [P][Pbar] = 1 " f"({name})", comp_off == 0, f"{comp_off} off") return rows def table(rows, cellnames): print(" cut " + "".join(f"{k:>34s}" for k in cellnames)) for j, x in enumerate(CUTS): line = f" {x:5d} " for k in cellnames: r = [row["raw"][j][k] for row in rows if k in row["raw"][j]] c = [row["cor"][j][k] for row in rows if k in row["cor"][j]] if not r: line += f"{'-':>34s}" continue (mr, sr), (mc, sc) = pooled(r), pooled(c) line += f" {mr:.4f}+-{sr:.4f} -> {mc:.4f}+-{sc:.4f}" print(line) def g_ladder(rows): out = [] for j, x in enumerate(CUTS): r = [row["raw"][j]["G"] for row in rows if "G" in row["raw"][j]] c = [row["cor"][j]["G"] for row in rows if "G" in row["cor"][j]] (mr, sr), (mc, sc) = pooled(r), pooled(c) gone = (mr - mc) / (mr - 1) if mr != 1 else float("nan") print(f" G at {x:5d}: raw {mr:.4f} +- {sr:.4f} ({(mr - 1) / sr:5.1f}" f" se), corrected {mc:.4f} +- {sc:.4f}, removed " f"{100 * gone:5.1f}% ({len(r)} fields)") out.append((mr, sr, mc, sc, gone)) return out def section_imag(): section("P THE SHORTFALL, imaginary: RAW -> PRIME-POWER levels") rows = population(-1) print(f" {len(rows)} fields with h >= 2, |D| <= {IMAG_D}") table(rows, ["T", "A", "S", "G"]) lad = g_ladder(rows) check("P1 G raw above 1 by >= 10 se at every cut", all((mr - 1) >= 10 * sr for mr, sr, _, _, _ in lad)) check("P1 G corrected within 0.01 of 1 at 10000", abs(lad[-1][2] - 1) <= 0.01, f"{lad[-1][2]:.4f}") check("P1 >= 70% of the raw G excess removed from cut 630 up", all(g >= 0.70 for _, _, _, _, g in lad[2:])) lo, hi = CUTS[0], TOP lead = [pi_below(isqrt(x) + 1) / pi_below(x) for x in (lo, hi)] print(f" from cut {lo} to {hi}: G's raw excess falls by " f"{(lad[0][0] - 1) / (lad[-1][0] - 1):.1f}, pi(sqrt x)/pi(x) by " f"{lead[0] / lead[1]:.1f}, pi(x) rises by " f"{pi_below(hi) / pi_below(lo):.1f}") j = len(CUTS) - 1 worst, read = [], 0 print(f" (h, cell) at {TOP}, cells with >= {MIN_CELL} fields:") for h in sorted({row["h"] for row in rows}): rs = [row for row in rows if row["h"] == h] for k in ("T", "A", "S", "G"): c = [row["cor"][j][k] for row in rs if k in row["cor"][j]] if len(c) < MIN_CELL: continue r = [row["raw"][j][k] for row in rs if k in row["raw"][j]] (mr, _), (mc, sc) = pooled(r), pooled(c) read += 1 ok = abs(mc - 1) <= 0.03 or abs(mc - 1) <= 3 * sc if not ok: worst.append((h, k, round(mc, 4))) print(f" h = {h:2d} {k}: n = {len(c):3d} raw {mr:.4f} " f"corrected {mc:.4f} +- {sc:.4f}") check("P2 every read (h, cell) corrected within 0.03 or 3 se of 1", read > 0 and not worst, f"{read} cells read, off: {worst}") pr = pi_below(isqrt(TOP) + 1) / pi_below(TOP) d = [1 - row["raw"][j]["T"] for row in rows] pd = [(row["h"] + row["two"]) / 2 * pr for row in rows] ratio = mean(d) / mean(pd) print(f" T at {TOP}: raw deficit {mean(d):.4f}, predicted " f"{mean(pd):.4f}, ratio {ratio:.3f}; corrected T " f"{mean(row['cor'][j]['T'] for row in rows):.4f}") check("P3 T's raw deficit over (h + |Cl[2]|)/2 pi(sqrt x)/pi(x) in " "[0.7, 1.3]", 0.7 <= ratio <= 1.3, f"{ratio:.3f}") pl = [((row["h"] + row["two"]) / 2 - 1) * pr for row in rows] print(f" the level's prediction, renormalised by the realised total, " f"((h + |Cl[2]|)/2 - 1) pi(sqrt x)/pi(x): {mean(pl):.4f}, ratio " f"{mean(d) / mean(pl):.3f}") def section_real(): section("P THE SHORTFALL, real, narrow group: RAW -> PRIME-POWER") rows = population(+1) print(f" {len(rows)} fields with h+ >= 2, D <= {REAL_D}") table(rows, ["T", "K0", "A", "S", "G"]) lad = g_ladder(rows) mr, _, mc, _, gone = lad[-1] check("P4 G corrected within 0.01 of 1 at 10000, >= 90% removed", abs(mc - 1) <= 0.01 and gone >= 0.90, f"{mc:.4f}, {100 * gone:.1f}%") Ls, lo = {}, 0 for D in range(5, REAL_D + 1): if fundamental(D): band = next(B for B in (1000, 2500, 5000, 10000, 16000) if D <= B) Ls.setdefault(band, []).append(least_principal_real(Group(D))) misses = 0 for B, v in sorted(Ls.items()): miss = sum(1 for L in v if L is None) misses += miss v = sorted(L for L in v if L is not None) print(f" real D in ({lo:5d}, {B:5d}]: {len(v) + miss:4d} fields, " f"median L_1 {v[len(v) // 2]:4d}, none below {TOP} at {miss}") lo = B check(f"R every real field to {REAL_D} has L_1 < {TOP}", misses == 0, f"{sum(len(v) for v in Ls.values())} fields, {misses} without") def main(): t0 = time.time() section_controls() section_race() section_floor() section_family() section_imag() section_real() print() print(f"{sum(CHECKS)}/{len(CHECKS)} checks PASS, " f"{time.time() - t0:.1f} s") if not all(CHECKS): raise SystemExit(1) if __name__ == "__main__": main()