dualdeg.py
Checks: DUAL: The two ends of a quadratic base criterion.
Cited on: DUAL.
To run it: save dualdeg.py; standard library only.
"""dualdeg.py -- the colon law past degree 2: do the two ends of a base
of degree n part as they do at a quadratic one?
QUESTION. f monic of degree n, f(x) = x^n + a_(n-1) x^(n-1) + ... +
a_0, q = |a_0| >= 2, O = Z[x]/(f) on the basis 1, beta, ..., beta^(n-1),
digits 0, ..., q - 1. The span-j leading tile is a translate of T_j =
{sum_(k<j) d_k beta^k}, the depth-j trailing tile a coset of beta^j O;
depth j reads an ideal I iff beta^j is in I, span j hides it iff T_j
equidistributes mod I. At degree 2 span j >= 1 hides I iff span 1 hides
(I : beta^(j-1)), iff I contains beta^(j-1) (m, beta - r) for some
m | q, f(r) = 0 mod m (dualquad.py). Does that colon law hold at
degree 3?
THE ARGUMENT (frozen before this file). a_0 != 0 makes beta a
non-zero-divisor, so [O : beta^j O] = q^j and T_j is complete mod
beta^j O.
(a) The sum of e(chi) over T_j, e(x) = exp(2 pi i x), is
prod_(k<j) sum_(d<q) e(d chi(beta^k)), factor k zero iff
beta^k chi is dead (chi(beta^k) != 0, q chi(beta^k) = 0). Reads
no degree.
(b) Span 1 hides J iff every nonzero chi of O/J is dead iff O/J = Z/m
generated by 1, m | q, iff J = (m, beta - r), f(r) = 0 mod m.
Reads no degree.
(c) The if half: T_j = T_(j-1) + beta^(j-1) D, T_(j-1) complete mod
beta^(j-1) O, and y -> beta^(j-1) y carries O/(I : beta^(j-1))
onto (beta^(j-1) O + I)/I. Reads no degree, and for every k <= j
gives: span j hides I if span k hides (I : beta^(j-k)).
(d) TRANSPLANT, the only if's step. At degree 2 (beta-bar psi)(beta)
= q psi(1) = 0 is the one coordinate past the first; at degree n,
beta* = -g(beta), g(x) = (f(x) - a_0)/x, so beta times beta* =
a_0 and (beta* psi)(beta^i) = a_0 psi(beta^(i-1)), known 0 only at
i = 1.
ENGINE. Ideals of index dividing q^J as the beta-stable row-HNF
sublattices of Z^n; hidden by counting T_j's keys mod I; the colon law
as "I contains beta^(j-1) m and beta^(j-1) (beta - r), m | q, f(r) = 0
mod m"; characters of O/I by back-substitution on the HNF.
PREDICTIONS (frozen). C0 at n = 2 the engine reproduces dualquad.py's
sweep: 100 bases, q^J <= 10^4, 17144 live pairs (index | q^j), 9085
hidden, the colon law 0 off, gaps {0: 1615, 1: 213, 2: 16}. C1 T_j
complete mod beta^j O at every cubic base. P1 span 1 as (b), 0 off. P2
the colon law 0 off at every cubic base; plant beta^(j-2) > 0 off. P3
the gap law 0 off; plant s + e > 0 off. S one-step dead ends, printed.
R P2's offs by reducibility. Bases x^3 + a x^2 + b x + c, c = +/-q, q
in 2..6, a, b in -2..2: 250, q^J <= 128.
KILL. C0 or C1 red voids the run; the colon law at degree 3 dies if P2
prints one off.
H1 (frozen after the kill fired, before its run). By (c), and since a
q-killed O/K is hidden at span k iff 1, ..., beta^(k-1) span it: span j
>= 1 hides I iff for some 1 <= k <= j, q beta^(j-k) is in I and 1,
..., beta^(k-1) span O/(I : beta^(j-k)). At degree 2, k = 1 is the
colon law and k = 2 follows from it. KILL: one off at one base.
FINDINGS (entered after the run, from its printed output).
C0 live 17144, hidden 9085, colon law 0 off, gap law 0 off, one-step
dead ends [0, 1574] (none at one prime, dualquad.py's S), H1 0 off; the
first run's gap tally read {0: 1715, 1: 213, 2: 16}, the 100 unit
ideals that dualquad.py's tally skips (s = 0); aligned, {0: 1615,
1: 213, 2: 16}. C1 green at all 250 cubic bases.
P1 0 off: (b) holds at degree 3.
P2 FIRES: 1492 off at 190 of 250 bases, 412 at the 44 reducible, 1080
at the 206 irreducible; plant 3576 off. One witness, x^3 + x^2 + x
+ 2: Q = (2, beta^2 + beta + 1), O/Q = F_4, beta a unit there, is
hidden at span 2 ({0, 1} + {0, beta} is F_4) and read at no depth,
while (Q : beta) = Q is not cyclic. At degree 2, f = x(x + p) mod
l | q has the root 0, so no prime of residue degree 2 lies over q.
Not the only kind: x^3 - 2x^2 - x + 2 = x(x + 1)^2 mod 2 has no
such prime and is off 21, 6 of them at the quotient F_2[e]/(e^2),
beta a unit, hidden from span 2 on, the other 15 at index 8 to 128
from span 3 on, 2O among them (T_3 is the whole digit cube).
W the witness as printed: O/Q has 4 classes, 1 generates 2 of them,
hidden at spans 2 to 7, read at no depth to 7, (Q : beta) = Q.
P3 944 off (it leans on P2); plant 2514 off.
S one-step dead ends 2474 at one prime, 454 mixed: the step fails
at degree 3, not only its proof.
H1 FIRES: 84 off at 36 bases, 32 of the 84 at the reducible ones (12
of the 36), every one hidden with H1 false (its if half is
proved); plant 1478 off. The freeze's split of the offs by k was
not printed. The tally: by (index, j) (36, 2) 40, (64, 3) 16,
(16, 2) 8, (32, 3) 8, (16, 6) 4, (16, 7) 4, (32, 7) 4; by q 40 at
6, 32 at 4, 12 at 2; the 40 at q = 6 all (36, 2), a mixed index;
the 12 with q beta^i in I at no i < j all at q = 2, (16, 6) 4,
(16, 7) 4, (32, 7) 4, at x^3 +/- x^2 + x +/- 2.
K The colon law is a degree-2 law. H1 holds at 12542 of the 12626
live pairs; its misses fall in three groups by q: 12 at q = 2
with no q beta^i in I (index 16 and 32, spans 6 and 7), 32 at
q = 4 (index 16 to 64, spans 2 and 3), 40 at q = 6, mixed index
36, where the CRT may be what H1 lacks. Toy scale: q^J <= 128.
RUN RECORD. The first C0 ran without S and H1; a run of C0 with them reads
dead ends [0, 1574] and H1 0 off, dualquad.py's S. 9 checks, 6 pass
(P2, P3, H1 fire as recorded), 48.8 s, 10.4 MB under a memory guard.
With C0's S and H1, the witness W and P2, P3, H1 asserting the recorded
misses: 13/13, 90.1 s, 13.5 MB.
"""
import sys
import time
CHECKS = []
def check(name, ok, detail=""):
CHECKS.append(bool(ok))
print(f" {'PASS' if ok else 'FAIL'} {name}" +
(f" ({detail})" if detail else ""))
def section(title):
print()
print(title)
def divisors(m):
return [k for k in range(1, m + 1) if m % k == 0]
def times_beta(v, a):
"""beta v, beta^n = -sum a_i beta^i."""
top = v[-1]
w = [0] + list(v[:-1])
return tuple(w[i] - top * a[i] for i in range(len(a)))
def mul(u, v, a):
"""u v in O = Z[x]/(f)."""
out = [0] * len(a)
w = tuple(v)
for ui in u:
if ui:
out = [o + ui * x for o, x in zip(out, w)]
w = times_beta(w, a)
return tuple(out)
def powers(a, kmax):
n = len(a)
out = [tuple(1 if i == 0 else 0 for i in range(n))]
for _ in range(kmax):
out.append(times_beta(out[-1], a))
return out
class Lat:
"""The row-HNF lattice: row i is d_i at column i and H[i][k] at
k > i, 0 <= H[i][k] < d_k; every full-rank sublattice once."""
def __init__(self, d, H):
self.d, self.H = d, H
self.N = 1
for x in d:
self.N *= x
def key(self, v):
v = list(v)
n = len(v)
out = []
for i in range(n):
qq, r = divmod(v[i], self.d[i])
out.append(r)
if qq:
Hi = self.H[i]
for k in range(i + 1, n):
v[k] -= qq * Hi[k]
return tuple(out)
def __contains__(self, v):
return not any(self.key(v))
def rows(self):
n = len(self.d)
return [tuple(self.d[i] if k == i else (self.H[i][k] if k > i
else 0)
for k in range(n)) for i in range(n)]
def characters(self):
"""c mod N with row . c = 0 mod N: chi(v) = c . v / N."""
n, N = len(self.d), self.N
out = [()]
for i in range(n - 1, -1, -1):
nxt = []
for tail in out:
R = -sum(self.H[i][k] * tail[k - i - 1]
for k in range(i + 1, n)) % N
di = self.d[i]
if R % di:
continue
for t in range(di):
nxt.append(((R // di + t * (N // di)) % N,) + tail)
out = nxt
return out
def factorizations(m, n):
if n == 1:
yield (m,)
return
for x in divisors(m):
for rest in factorizations(m // x, n - 1):
yield (x,) + rest
def uppers(d):
"""Every upper part of a row HNF with diagonal d."""
n = len(d)
cells = [(i, k) for i in range(n) for k in range(i + 1, n)]
def rec(c, H):
if c == len(cells):
yield [row[:] for row in H]
return
i, k = cells[c]
for x in range(d[k]):
H[i][k] = x
yield from rec(c + 1, H)
H[i][k] = 0
yield from rec(0, [[0] * n for _ in range(n)])
def ideals(a, top):
"""The ideals of O of index dividing top: beta-stable lattices."""
n = len(a)
out = []
for N in divisors(top):
for d in factorizations(N, n):
for H in uppers(d):
L = Lat(d, H)
if all(times_beta(r, a) in L for r in L.rows()):
out.append(L)
return out
def tiles(a, q, J, pw):
T = [[tuple([0] * len(a))]]
for j in range(1, J + 1):
b = pw[j - 1]
T.append([tuple(x + dd * y for x, y in zip(t, b))
for t in T[-1] for dd in range(q)])
return T
def hidden_by_count(L, Tj):
if len(Tj) % L.N:
return False
cnt = {}
for v in Tj:
k = L.key(v)
cnt[k] = cnt.get(k, 0) + 1
return len(cnt) == L.N and len(set(cnt.values())) == 1
def fval(a, r):
return r ** len(a) + sum(c * r ** i for i, c in enumerate(a))
def reducible(a):
"""A rational root, which for n <= 3 is reducibility over Q."""
a0 = abs(a[0])
return any(fval(a, s * r) == 0 for r in divisors(a0) for s in (1, -1))
def prime_power(m):
if m == 1:
return False
l = next(k for k in range(2, m + 1) if m % k == 0)
while m % l == 0:
m //= l
return m == 1
def run_base(a, points):
"""One base, a = (a_0, ..., a_(n-1)): the counts, the laws, the step."""
n, q = len(a), abs(a[0])
J = 0
while q ** (J + 1) <= points:
J += 1
pw = powers(a, J + 1)
Is = ideals(a, q ** J)
T = tiles(a, q, J, pw)
one = pw[0]
# C1: T_j complete mod beta^j O
okc = len(set(T[J])) == q ** J
for j in range(J + 1):
B = [L for L in Is if L.N == q ** j and pw[j] in L]
okc &= len(B) == 1 and len({B[0].key(v) for v in T[j]}) == q ** j
mroots = [(m, r) for m in divisors(q) for r in range(m)
if fval(a, r) % m == 0]
cst = lambda m: tuple(m * x for x in one)
br = {r: tuple(x - (r if i == 0 else 0) for i, x in
enumerate(pw[1])) for r in range(q)}
def law(L, j, shift=1):
if j == 0:
return L.N == 1
p = pw[max(j - shift, 0)]
return any(mul(p, cst(m), a) in L and mul(p, br[r], a) in L
for m, r in mroots)
live = hid = off = off1 = plant = gapoff = gapplant = 0
gaps = {}
for L in Is:
cnt = [hidden_by_count(L, T[j]) for j in range(J + 1)]
read = [pw[j] in L for j in range(J + 1)]
for j in range(J + 1):
lw = law(L, j)
off += lw != cnt[j]
if j == 1:
off1 += lw != cnt[j]
if j >= 1:
plant += law(L, j, 2) != cnt[j]
if (q ** j) % L.N == 0:
live += 1
hid += cnt[j]
s = next((j for j in range(J + 1) if cnt[j]), None)
t = next((j for j in range(J + 1) if read[j]), None)
if s and s >= 1:
mr = [(m, r) for m, r in mroots
if mul(pw[s - 1], cst(m), a) in L
and mul(pw[s - 1], br[r], a) in L]
m, r = min(mr) if mr else (1, 1)
e = next((e for e in range(1, m + 1) if pow(r, e, m) == 0),
None) if mr else -J - 9
for shift, bad in ((0, "g"), (1, "p")):
want = None if e is None else s - 1 + e + shift
want = want if want is not None and want <= J else None
if bad == "g":
gapoff += want != t
else:
gapplant += want != t
if s and t is not None:
gaps[t - s] = gaps.get(t - s, 0) + 1
def colon(L, i):
"""(L : beta^i), the largest listed K with beta^i K in L."""
ks = [K for K in Is if all(mul(pw[i], r_, a) in L
for r_ in K.rows())]
return min(ks, key=lambda K: K.N)
def spans(K, k):
"""Do 1, beta, ..., beta^(k-1) span O/K over Z?"""
gens = [K.key(pw[i]) for i in range(k)]
seen = {tuple([0] * n)}
frontier = list(seen)
while frontier:
nxt = []
for v in frontier:
for g in gens:
w = K.key(tuple(x + y for x, y in zip(v, g)))
if w not in seen:
seen.add(w)
nxt.append(w)
frontier = nxt
return len(seen) == K.N
def h1(L, j, kmax):
for k in range(1, kmax + 1):
i = j - k
if mul(pw[i], cst(q), a) in L and spans(colon(L, i), k):
return True
return False
h1off = h1plant = h1q = 0
h1dir = [0, 0]
h1bad = []
for L in Is:
for j in range(1, J + 1):
hid_ = hidden_by_count(L, T[j])
v = h1(L, j, j)
h1off += v != hid_
if v != hid_:
h1q += not any(mul(pw[i], cst(q), a) in L
for i in range(j))
h1dir[hid_] += 1
h1bad.append((L.N, j, hid_))
h1plant += h1(L, j, j - 1) != hid_
dead = [0, 0]
for L in Is:
N = L.N
if N == 1:
continue
inB = lambda c: c[0] != 0 and (q * c[0]) % N == 0
pre = {}
for c in L.characters():
img = tuple(c[1:]) + (-sum(x * y for x, y in zip(a, c)) % N,)
pre.setdefault(img, []).append(c)
z = tuple([0] * n)
bad = sum(1 for t_, cs in pre.items() if t_ != z and not inB(t_)
and all(inB(c) for c in cs))
dead[0 if prime_power(N) else 1] += bad
return dict(a=a, n=n, q=q, J=J, ideals=len(Is), okc=okc, live=live,
hid=hid, off=off, off1=off1, plant=plant, gapoff=gapoff,
gapplant=gapplant, gaps=gaps, dead=dead, h1off=h1off,
h1plant=h1plant, h1bad=h1bad, h1q=h1q, h1dir=h1dir,
red=reducible(a))
def fmt(a):
n = len(a)
terms = [f"x^{n}"]
for i in range(n - 1, -1, -1):
c = a[i]
if c == 0:
continue
mon = "" if i == 0 else ("x" if i == 1 else f"x^{i}")
coef = str(abs(c)) if (abs(c) != 1 or i == 0) else ""
terms.append(("- " if c < 0 else "+ ") + coef + mon)
return " ".join(terms)
def section_control():
section("C0 the positive control: dualquad.py's sweep at n = 2")
bases = [(q, p) for q in range(2, 13) for p in range(-1, q + 1)]
bases.append((17, 8))
tot = dict(live=0, hid=0, off=0, gapoff=0, h1off=0)
gaps = {}
dead = [0, 0]
for q, p in bases:
R = run_base((q, p), 10 ** 4)
for k in tot:
tot[k] += R[k]
for g, c in R["gaps"].items():
gaps[g] = gaps.get(g, 0) + c
dead = [x + y for x, y in zip(dead, R["dead"])]
gaps = dict(sorted(gaps.items()))
print(f" {len(bases)} bases: live {tot['live']}, hidden {tot['hid']}, "
f"colon law off {tot['off']}, gap law off {tot['gapoff']}, "
f"gaps {gaps}, dead ends {dead}, H1 off {tot['h1off']}")
check("C0 the quadratic sweep reproduced",
len(bases) == 100 and tot["live"] == 17144 and tot["hid"] == 9085
and tot["off"] == 0 and tot["gapoff"] == 0
and gaps == {0: 1615, 1: 213, 2: 16})
check("C0 dualquad.py's S: one-step dead ends only at mixed index",
dead == [0, 1574], f"{dead}")
check("C0 H1 holds at every quadratic pair", tot["h1off"] == 0)
def section_witness():
section("W the witness: x^3 + x^2 + x + 2, Q = (2, beta^2 + beta + 1)")
a, q, J = (2, 1, 1), 2, 7
pw = powers(a, J + 1)
T = tiles(a, q, J, pw)
Is = ideals(a, q ** J)
Q = [L for L in Is if L.N == 4 and (2, 0, 0) in L and (1, 1, 1) in L]
if len(Q) != 1:
check("W Q is one ideal of index 4", False, f"{len(Q)} found")
return
Q = Q[0]
hid = [hidden_by_count(Q, T[j]) for j in range(J + 1)]
read = [pw[j] in Q for j in range(J + 1)]
colon = [K for K in Is if all(mul(pw[1], r_, a) in Q for r_ in K.rows())]
colon = min(colon, key=lambda K: K.N)
by1 = len({Q.key((k, 0, 0)) for k in range(Q.N)})
spans = [j for j in range(J + 1) if hid[j]]
depths = [j for j in range(J + 1) if read[j]]
print(f" O/Q has {Q.N} elements; hidden at spans {spans}, read at "
f"depths {depths}; "
f"(Q : beta) = Q {colon.rows() == Q.rows()}, 1 generates {by1} "
f"of its {Q.N} classes")
check("W Q is hidden from span 2 on and read at no depth to 7",
not hid[1] and all(hid[2:]) and not any(read))
check("W (Q : beta) = Q, 1 does not generate O/Q: the colon law "
"misses Q at span 2",
colon.rows() == Q.rows() and by1 < Q.N)
def section_cubic():
section("P the cubic bases: x^3 + a x^2 + b x + c, c = +/-q")
rows = []
for q in range(2, 7):
for c in (q, -q):
for a2 in range(-2, 3):
for a1 in range(-2, 3):
rows.append(run_base((c, a1, a2), 128))
okc = all(R["okc"] for R in rows)
check("C1 T_j complete mod beta^j O at every cubic base", okc,
f"{sum(not R['okc'] for R in rows)} red")
agg = lambda k, sel=lambda R: True: sum(R[k] for R in rows if sel(R))
print(f" {len(rows)} bases, {agg('ideals')} ideals, "
f"{agg('live')} live pairs, {agg('hid')} hidden")
for q in range(2, 7):
sel = lambda R, q=q: R["q"] == q
print(f" q = {q}: J = {next(R['J'] for R in rows if sel(R))}, "
f"ideals {agg('ideals', sel)}, live {agg('live', sel)}, "
f"hidden {agg('hid', sel)}, off {agg('off', sel)}, "
f"dead ends one prime {sum(R['dead'][0] for R in rows if sel(R))}"
f", mixed {sum(R['dead'][1] for R in rows if sel(R))}")
bad = [R for R in rows if R["off"]]
for R in bad[:12]:
print(f" OFF {fmt(R['a'])}: off {R['off']}, "
f"span 1 off {R['off1']}, reducible {R['red']}")
check("P1 span 1 hides I iff I contains some (m, beta - r)",
agg("off1") == 0, f"{agg('off1')} off")
check("P2 the colon law misses at degree 3: 1492 off at 190 bases",
agg("off") == 1492 and len(bad) == 190,
f"{agg('off')} off at {len(bad)} bases")
check("P2 plant beta^(j-2) is off", agg("plant") > 0,
f"{agg('plant')} off")
gaps = {}
for R in rows:
for g, c in R["gaps"].items():
gaps[g] = gaps.get(g, 0) + c
check("P3 the gap law, leaning on P2, misses 944", agg("gapoff") == 944,
f"{agg('gapoff')} off; gaps {dict(sorted(gaps.items()))}")
check("P3 plant s + e is off", agg("gapplant") > 0,
f"{agg('gapplant')} off")
hb = [R for R in rows if R["h1off"]]
for R in hb[:6]:
print(f" H1 OFF {fmt(R['a'])}: {R['h1off']} off, e.g. "
f"{R['h1bad'][:2]}")
byij, byq = {}, {}
for R in hb:
for N, j, _ in R["h1bad"]:
byij[(N, j)] = byij.get((N, j), 0) + 1
byq[R["q"]] = byq.get(R["q"], 0) + 1
print(f" H1's misses by (index, j) "
f"{dict(sorted(byij.items(), key=lambda x: -x[1]))}; by q "
f"{dict(sorted(byq.items()))}")
check("H1 misses at degree 3: 84 off at 36 bases", agg("h1off") == 84
and len(hb) == 36, f"{agg('h1off')} off at {len(hb)} bases; "
f"reducible {agg('h1off', lambda R: R['red'])}")
print(f" H1's misses: {agg('h1q')} with q beta^i in I at no i < j; "
f"hidden and H1 false {sum(R['h1dir'][1] for R in rows)}, "
f"unhidden and H1 true {sum(R['h1dir'][0] for R in rows)}")
check("H1 plant k <= j - 1 is off", agg("h1plant") > 0,
f"{agg('h1plant')} off")
nred = sum(R["red"] for R in rows)
print(f" R reducible bases {nred}: off {agg('off', lambda R: R['red'])};"
f" irreducible {len(rows) - nred}: "
f"off {agg('off', lambda R: not R['red'])}")
d1 = sum(R["dead"][0] for R in rows)
d2 = sum(R["dead"][1] for R in rows)
print(f" S one-step dead ends: one prime {d1}, mixed {d2}")
return agg("off")
def main():
t0 = time.time()
if "--cubic" not in sys.argv:
section_control()
section_witness()
off = section_cubic()
section("K the kill")
print(f" the colon law at degree 3: {'DIES' if off else 'STANDS'}")
print()
print(f"{sum(CHECKS)}/{len(CHECKS)} checks pass, "
f"{time.time() - t0:.1f} s")
return 0 if all(CHECKS) else 1
if __name__ == "__main__":
sys.exit(main())