""" pairs.py -- the pair automaton G_d: every two codings of one point of P_m, over every m at once, as the paths of one finite graph, and its live part as a carry-window automaton K_d. QUESTION. Two codings y, y' of one point start at s_0, s'_0 with equal first coordinate and the rest equal mod m; E_k = (s_k - s'_k)/m runs E_(k+1) = M^(-1)(E_k - (y_k - y'_k) v_0) from E_0 in {0} x Z^(d-1), integral and bounded, m gone. z is multiply coded iff one coding of z is the y-label of a non-diagonal infinite path of G_d (nodes (E, a, a'), a, a' the runs) from a start node ((0, e), 0, 0). What is G_d's live part, in one d-shape, and does it read sigma single and -sigma multiple, as fibre.py does? DESIGN. fibre.py's eigen-coordinates: E's contracting coordinate l lies in the polygon of b_l times a difference of legal tail sums, support h(u) + h(-u); the beta-coordinate is bounded by max(S0, |b_beta|/(beta - 1)), S0 the bound a first coordinate 0 puts on it at the start. The start nodes: the lattice points (0, e) inside, enumerated over the bounding box of the polygon bounds. G_d by BFS, its live part by peeling sinks (fibre.peel). A product with a lasso word u V^inf (y forced) adds the phase; a path is non-diagonal once it leaves E = 0 or starts off it, read as: the live product holds a start node with E != 0, or a node with E = 0 and y != y' with a live successor. A product run after G_d is seeded from G_d's live starts and kept in its live set. The carry coordinates, derived by hand from the step E' = (E_2 - delta, E_3, ..., E_d, E_1 - E_2 - ... - E_d), delta = y - y'. With w = (E_2, ..., E_d, E_1 - E_2 - ... - E_d), so that E = (sum w, w_1, ..., w_(d-1)), the step is a shift: w' = (w_2, ..., w_d, c), c = 2 w_1 - sum w - delta, and sum w' = w_1 - delta. So a pair of codings is one integer carry sequence x, w its length-d windows, with delta_m = x_m - (x_(m+1) + ... + x_(m+d)): y - y' is, up to a term at the start, the carry 1^d 0 <-> 0^d 1 applied x_m times at position m. A start (0, e) is the window (e, -sum e), of sum 0. Setting y = y' = 0 wherever delta = 0 is never worse, so a continuation is legal iff delta's runs of +1 and of -1 stay below d, the first counted on from a or a'. K_d, the carry graph: nodes (w, a, a'), w in {-1, 0, 1}^d with sum in {-1, 0, 1}, those edges, starts (w, 0, 0) with sum w = 0: G_d's step on the E in {-1, 0, 1}^d with |E_1 - E_2 - ... - E_d| <= 1, with no box. Its live part reachable from its starts is G_d's iff every live E of G_d lies in {-1, 0, 1}^d (P1; the second bound then holds, an entry of each successor) and the box cuts no infinite path of K_d (a bounded path's contracting coordinates are tail sums, its beta-coordinate within the bound the recursion keeps). Every zero-sum w is a live start of K_d, at every d: delta = w_1 keeps the sum 0 and appends c = w_1, so the window rotates and the pair is y = [w = 1]^inf, y' = [w = -1]^inf, period d, legal as a nonzero zero-sum w holds both signs. The two words code one point by periodic.py's criterion: the E-path returns after d steps, so E_0 is the periodic solution, E_0 = -iota(X_y - X_y'), first coordinate 0. SK: the products of K_d with a lasso word u V^inf (y forced, y' free), in exact integers, a depth-first search for the longest non-diagonal path, one arm at a time with its own memo: the START arm from the nonzero starts, the BRANCH arm off the diagonal at E = 0 at every time up to len(u) + 2 len(V); a node met again on the stack is an infinite path. Words: S^inf and A^inf; R = (1^(d-1) 0)^inf, a coding of -alpha, which has d codings, its partners found on the start arm; Q, read off K_d's live set (branch_word: an edge y != y' from the zero start, then live edges to a repeat), whose partner branches at time 0. PREDICTIONS, CONTROLS AND KILL (fixed before the runs). C1 CONTROL, read first, d = 4..8: the S-product has no non-diagonal live path from a start node (sigma single, as fibre.py reads), and the A-product has one, its start E_0 = +-e_3 (by hand: s_A - s_B = -4 iota(beta^2 - beta - 1) = -4 e_3). A mismatch voids the run. C2 CONTROL, d = 3, 4: every census class fibre.py labels multiple at m = 2..4 has its first coding paired, every single one not. PLANT: the box halved loses the A/B pair at some d. P1 the live E-set's entries lie in {-1, 0, 1} at every d read. P2 (TRANSPLANT from tilings, no confidence) 2(2^(d-1) - 1) nonzero E's among the live start nodes. KILL of the d-uniform route: P1 fails at d = 8 with an entry of size >= 3, or the live E count at d = 8 exceeds 10^4. Fixed after d = 3..7 were read as a fit, before d = 8 was run: K1 at d = 8, G_8's live set read in w is K_8's live part reachable from its starts, symmetric difference 0. K2 at d = 8, G_8's live starts are the zero-sum windows. CONTROL, read first: K_d with c's delta sign flipped misses, a nonzero symmetric difference at d = 3, 4, 5. KILL of the rule: a nonzero symmetric difference at d = 8. SK's control: the A-product does not die, -sigma being multiply coded. Added after the S-product's read, before the run that ships: R lives on the start arm, Q on the branch arm. FINDINGS (entered after the runs; prints copied from them). - C1 holds at d = 4..8: the S-product has no non-diagonal live path (live 10, 12, 14, 16, 18 nodes, the cycle alone), the A-product's starts at E_0 = (0, 0, -1, 0, ..., 0) = -e_3. C2 holds at d = 3, 4, m = 2..4: paired iff labelled multiple, no exception. The PLANT loses the A pair at every d = 4..8. - KILL misses; P1 holds at d = 3..8: every live E has entries in {-1, 0, 1}. G_d: nodes 33, 232, 1225, 13428, 92777, 1144068; live 33, 130, 463, 1566, 5117, 16402; live E's 13, 35, 93, 253, 697, 1945; live nonzero starts 6, 18, 50, 140, 392, 1106. P2 holds at d = 3 alone (6) and fails from d = 4 (18 against 14). - The rule's KILL misses; K1 and K2 hold at d = 3..8, d = 8 the frozen test: G_d's live part is K_d's reachable live part, 33, 130, 463, 1566, 5117, 16402 nodes in both; the flipped control misses by 30, 126, 458. So the live starts are the zero-sum windows of {-1, 0, 1}^d, T_d - 1 nonzero ones, T_d the central trinomial coefficient: the 6, 18, 50, 140, 392, 1106 above. - SK: inside K_d the S-product has no infinite non-diagonal path at d = 4..11 on either arm, its longest from a nonzero start 7, 13, 17, 26, 33, 40, 56, 63 steps and off the diagonal 15, 23, 29, 36, 44, 55, 65, 82; the A-product has one on the start arm at each d (C1's start at d = 4..8). The controls hold at d = 4..11: R lives on the start arm, Q on the branch arm. Inside K_d is all of G_d's live part where K1 holds (d <= 8); past it the read covers the second codings whose carries and window sums lie in {-1, 0, 1}. RUN RECORD: the bare run, d = 3..8, 204.1 s, peak commit 326.5 MB (G_8 184.9 s; the products 0.3 s seeded from G_d's live starts and kept in its live set, the box check skipped as that set lies in the box, against 316.1 s seeded from the box). SK at d = 12 passes 512 MB in its memo, so it stops at 11 (9.1 s with its controls). --plant d = 4..8 4.5 s, 15.3 MB. """ import cmath import itertools from array import array import math import sys import time from fibre import DIRS, Minv_col, Num, census, eigen, support CHECKS = [] def check(name, ok, detail=""): CHECKS.append(bool(ok)) print(f" [{'PASS' if ok else 'FAIL'}] {name}" + (f" -- {detail}" if detail else ""), flush=True) class Box: """The box of E: polygons for the contracting coordinates, a bound on the beta-coordinate.""" def __init__(self, N, shrink=1.0): d = self.d = N.dim lam, b, rows = eigen(N) self.rows = rows ib = self.ib = max(range(d), key=lambda i: abs(lam[i])) beta = lam[ib].real dirs = [cmath.exp(2j * math.pi * t / DIRS) for t in range(DIRS)] self.poly = {i: [(u, shrink * (support(lam[i], b[i], d, u) + support(lam[i], b[i], d, -u)) + 1e-7) for u in dirs] for i in range(d) if i != ib} bnd = [0.0] * d for i, hs in self.poly.items(): bnd[i] = max(h for _, h in hs) / math.cos(math.pi / DIRS) S0 = sum(bnd[i] * abs(lam[i]) ** (d - 1) for i in self.poly) \ / beta ** (d - 1) bnd[ib] = shrink * max(S0, abs(b[ib]) / (beta - 1)) self.bnd = [x * (1 + 1e-9) + 1e-9 for x in bnd] self.lam, self.b = lam, b R = [[l ** (d - 1 - i) for l in lam] for i in range(d)] self.top = [int(sum(abs(R[j][l]) * self.bnd[l] for l in range(d))) + 1 for j in range(d)] def inside(self, s): for i in range(self.d): z = sum(self.rows[i][j] * s[j] for j in range(self.d)) if abs(z) > self.bnd[i]: return False if i != self.ib and any((u * z).real > h for u, h in self.poly[i]): return False return True def starts(self): """The lattice points (0, e) inside: e over the coordinate box that the eigen-bounds give, |s_j| <= sum_l |r_l,j| bnd_l.""" out = [] for e in itertools.product(*[range(-t, t + 1) for t in self.top[1:]]): s = (0,) + e if self.inside(s): out.append(s) return out class Nodes: """The nodes (E, a, a', phase) packed as integers, decoded on read: G_d at d = 8 holds about 10^6 of them.""" def __init__(self, d, rad, P): self.d, self.rad, self.P, self.codes = d, rad, P, [] def pack(self, nd): E, a, a2, ph = nd c = 0 for e in E: c = c * self.rad + e + self.rad // 2 return ((c * self.d + a) * self.d + a2) * self.P + ph def __getitem__(self, i): c = self.codes[i] c, ph = divmod(c, self.P) c, a2 = divmod(c, self.d) c, a = divmod(c, self.d) E = [] for _ in range(self.d): c, e = divmod(c, self.rad) E.append(e - self.rad // 2) return tuple(E[::-1]), a, a2, ph def __len__(self): return len(self.codes) def build(N, box, word=None, allow=None): """G_d, or its product with the lasso word (u, V): nodes (E, a, a', phase); an edge labelled (y, y') to node t stored as 4t + 2y + y'. With allow, only nodes whose (E, a, a') is in it: a live product path projects to a live G_d path. Returns nodes, succ, the number of starts.""" d, v0 = N.dim, N.v0 seq = None if word is not None: u, V = word seq = u + V def nxt(ph): return ph + 1 if ph + 1 < len(seq) else len(u) nodes = Nodes(d, 2 * max(box.top) + 1, len(seq) if seq else 1) idx = {} seeds = (box.starts() if allow is None else sorted(E for E, a, a2 in allow if (a, a2, E[0]) == (0, 0, 0))) for E in seeds: if allow is None or (E, 0, 0) in allow: c = nodes.pack((E, 0, 0, 0)) idx[c] = len(nodes.codes) nodes.codes.append(c) nst = len(nodes.codes) succ = [] i = 0 while i < len(nodes.codes): E, a, a2, ph = nodes[i] out = [] ys = (0, 1) if seq is None else (seq[ph],) for y in ys: for y2 in (0, 1): r1, r2 = (a + 1 if y else 0), (a2 + 1 if y2 else 0) if r1 >= d or r2 >= d: continue E1 = Minv_col(N, tuple(e - (y - y2) * w for e, w in zip(E, v0))) if allow is not None: if (E1, r1, r2) not in allow: continue elif not box.inside(E1): continue c = nodes.pack((E1, r1, r2, 0 if seq is None else nxt(ph))) t = idx.get(c) if t is None: t = idx[c] = len(nodes.codes) nodes.codes.append(c) out.append(4 * t + 2 * y + y2) succ.append(out) i += 1 return nodes, succ, nst def live_part(succ): """fibre.peel on flat arrays: the nodes with an infinite continuation (the predecessor lists of G_d at d = 8 pass 512 MB as lists).""" nn = len(succ) cnt = array("i", (len(o) for o in succ)) off = array("i", [0]) * (nn + 1) for o in succ: for x in o: off[(x >> 2) + 1] += 1 for u in range(nn): off[u + 1] += off[u] pred, fill = array("i", [0]) * off[nn], array("i", off[:nn]) for u, o in enumerate(succ): for x in o: t = x >> 2 pred[fill[t]] = u fill[t] += 1 del fill live = bytearray([1]) * nn stack = [u for u in range(nn) if cnt[u] == 0] while stack: u = stack.pop() if not live[u]: continue live[u] = 0 for k in range(off[u], off[u + 1]): q = pred[k] cnt[q] -= 1 if cnt[q] == 0 and live[q]: stack.append(q) return live def paired(nodes, succ, nst, live): """A non-diagonal live path from a start: a live start off E = 0, or a live node on E = 0 reachable along the diagonal with a live edge y != y'.""" zero = tuple([0] * len(nodes[0][0])) off = [i for i in range(nst) if live[i] and nodes[i][0] != zero] if off: return ("start", nodes[off[0]][0]) seen = [i for i in range(nst) if live[i] and nodes[i][0] == zero] stack, mark = list(seen), set(seen) while stack: u = stack.pop() for x in succ[u]: t, y, y2 = x >> 2, (x >> 1) & 1, x & 1 if not live[t]: continue if y != y2: return ("branch", nodes[u]) if t not in mark: mark.add(t) stack.append(t) return None def words(d): S = [1] + [0] * (d + 1) + [1] * (d - 3) + [0, 1, 0] A = ([0, 0] + [1] * (d - 1) + [0] * (d - 2) + [1, 0, 0, 0] + [1] * (d - 2) + [0, 1, 1, 0] + [1] * (d - 4) + [0, 0, 1]) return S, A def lasso(y, top=150): """y (400 digits) as (u, V) with y = u V^inf on its digits.""" for p in range(1, top): s = max((i + 1 for i in range(len(y) - p) if y[i] != y[i + p]), default=0) if s < len(y) - 2 * p - 50: return y[:s], y[s:s + p] return None def w_of(E): """The carry window of E: (E_2, ..., E_d, E_1 - E_2 - ... - E_d).""" return tuple(E[1:]) + (E[0] - sum(E[1:]),) def carry(d, flip=False): """K_d's live part reachable from its starts: nodes (w, a, a'), w in {-1, 0, 1}^d with sum in {-1, 0, 1}, the step w -> (w_2, ..., w_d, c), c = 2 w_1 - sum w - delta, sum w' = w_1 - delta; starts (w, 0, 0) with sum w = 0. With flip, c = 2 w_1 - sum w + delta (the control).""" idx, nodes, succ = {}, [], [] for w in itertools.product((-1, 0, 1), repeat=d): if sum(w) == 0: idx[(w, 0, 0)] = len(nodes) nodes.append((w, 0, 0)) i = 0 while i < len(nodes): w, a, a2 = nodes[i] s, out = sum(w), [] for y in (0, 1): for y2 in (0, 1): r1, r2 = (a + 1 if y else 0), (a2 + 1 if y2 else 0) dl = y - y2 c = 2 * w[0] - s + (dl if flip else -dl) if r1 >= d or r2 >= d or abs(c) > 1 or abs(w[0] - dl) > 1: continue nd = (w[1:] + (c,), r1, r2) t = idx.get(nd) if t is None: t = idx[nd] = len(nodes) nodes.append(nd) out.append(4 * t + 2 * y + y2) succ.append(out) i += 1 live = live_part(succ) return {nodes[i] for i in range(len(nodes)) if live[i]} def section_k(alive, starts): print("K G_d's live part against K_d's, read in w", flush=True) for d in sorted(alive): if d <= 5: bad = carry(d, flip=True) ^ {(w_of(E), a, a2) for E, a, a2 in alive[d]} check(f"control d = {d}: the flipped step misses", bad, f"symmetric difference {len(bad)}") for d in sorted(alive): G = {(w_of(E), a, a2) for E, a, a2 in alive[d]} K = carry(d) zs = {(0,) + w[:-1] for w in itertools.product((-1, 0, 1), repeat=d) if sum(w) == 0} check(f"K1 d = {d}: G_d live = K_d reachable live", G == K, f"G {len(G)}, K {len(K)}, symmetric difference {len(G ^ K)}") check(f"K2 d = {d}: live starts = the zero-sum windows", starts[d] == zs, f"{len(zs) - 1} nonzero, symmetric " f"difference {len(starts[d] ^ zs)}") def s_life(d, word, arm): """The product of K_d with the lasso word (u, V) (y = u V^inf forced, y' free): the longest non-diagonal path of one arm, "start" (from a nonzero start) or "branch" (off the diagonal at E = 0, at every time up to len(u) + 2 len(V), past which time repeats phase and run), None if one is infinite (a node met again on the stack). Each arm has its own memo, so a cycle one arm met cannot hide from the other.""" u, V = word seq, P, memo, inf = u + V, len(u) + len(V), {}, [] def nxt(ph): return ph + 1 if ph + 1 < P else len(u) def life(w, a2, ph): key = (sum((x + 1) * 3 ** i for i, x in enumerate(w)) * d + a2) * P + ph if key in memo: if memo[key] is None: inf.append(key) return 0 return memo[key] memo[key] = None s, out = sum(w), 0 for y2 in (0, 1): dl, r2 = seq[ph] - y2, (a2 + 1 if y2 else 0) c = 2 * w[0] - s - dl if r2 < d and abs(c) <= 1 and abs(w[0] - dl) <= 1: out = max(out, 1 + life(w[1:] + (c,), r2, nxt(ph))) memo[key] = out return out sys.setrecursionlimit(10000) worst = 0 if arm == "start": worst = max(life(w, 0, 0) for w in itertools.product((-1, 0, 1), repeat=d) if sum(w) == 0 and any(w)) else: ph, run = 0, 0 for _ in range(len(u) + 2 * len(V)): y2 = 1 - seq[ph] if not y2 or run + 1 < d: worst = max(worst, 1 + life((0,) * (d - 1) + (y2 - seq[ph],), run + 1 if y2 else 0, nxt(ph))) run = run + 1 if seq[ph] else 0 ph = nxt(ph) return None if inf else worst def branch_word(d): """A word whose second coding branches off the diagonal, read off K_d's live set: from the zero start an edge with y != y' into a live node, then the first live edge at each node until a node repeats, the labels y read as (u, V).""" K = carry(d) nd, seen, path = ((0,) * d, 0, 0), {}, [] while nd not in seen: seen[nd] = len(path) w, a, a2 = nd nexts = [(y, (w[1:] + (2 * w[0] - sum(w) - (y - y2),), a + 1 if y else 0, a2 + 1 if y2 else 0)) for y in (0, 1) for y2 in (0, 1) if path or y != y2] y, nd = next((y, t) for y, t in nexts if t in K) path.append(y) return path[:seen[nd]], path[seen[nd]:] def section_sk(top): print("SK sigma's S-product inside K_d, exact integers, by arm", flush=True) for d in range(4, top + 1): t0 = time.time() S, A = words(d) R = [1] * (d - 1) + [0] got = {(nm, arm): s_life(d, V, arm) for nm, V in (("S", ([], S)), ("A", ([], A)), ("R", ([], R)), ("Q", branch_word(d))) for arm in ("start", "branch")} txt = ", ".join(f"{nm} {arm} {'inf' if v is None else v}" for (nm, arm), v in got.items()) check(f"d = {d}: S dies on both arms; A and R live on the start " f"arm, Q on the branch arm", got[("S", "start")] is not None and got[("S", "branch")] is not None and got[("A", "start")] is None and got[("R", "start")] is None and got[("Q", "branch")] is None, f"{txt} [{time.time() - t0:.1f} s]") def section_g(top): print("G_d the live part", flush=True) p1, alive, starts = [], {}, {} for d in range(3, top + 1): N, t0 = Num(d), time.time() box = Box(N) nodes, succ, nst = build(N, box) live = live_part(succ) alive[d] = {nodes[i][:3] for i in range(len(nodes)) if live[i]} Es = sorted({nd[0] for nd in alive[d]}) starts[d] = {nodes[i][0] for i in range(nst) if live[i]} st = sorted(starts[d] - {tuple([0] * d)}) big = max(max(abs(x) for x in E) for E in Es) if big > 1: p1.append((d, big)) print(f" d = {d}: nodes {len(nodes)}, live {sum(live)}, live E's " f"{len(Es)}, live nonzero starts {len(st)} (P2 " f"{2 * (2 ** (d - 1) - 1)}), max |entry| {big} " f"[{time.time() - t0:.1f} s]", flush=True) if d <= 5: print(f" starts {st}", flush=True) print(f" E's {Es}", flush=True) del nodes, succ, live print(f"P1: {'holds' if not p1 else 'fails ' + str(p1)}") section_k(alive, starts) return alive def main(): shrink = 0.5 if "--plant" in sys.argv else 1.0 top = (int(sys.argv[sys.argv.index("--d") + 1]) if "--d" in sys.argv else 8) alive = section_g(top) if shrink == 1.0 else {} if shrink == 1.0: section_sk(11) if "--rule" in sys.argv: print(f"{sum(CHECKS)}/{len(CHECKS)} checks passed") sys.exit(0 if all(CHECKS) else 1) print("C1 sigma's S and -sigma's A against G_d", flush=True) for d in range(4, top + 1): N, t0 = Num(d), time.time() box = Box(N, shrink) S, A = words(d) rs, ra = [], [] for w, acc in ((S, rs), (A, ra)): nodes, succ, nst = build(N, box, ([], w), alive.get(d)) live = live_part(succ) acc.append(paired(nodes, succ, nst, live)) acc.append((len(nodes), sum(live))) e3 = tuple(-1 if j == 2 else 0 for j in range(d)) ok = (rs[0] is None and ra[0] is not None and ra[0][0] == "start" and ra[0][1] in (e3, tuple(-x for x in e3))) check(f"d = {d}: S unpaired, A paired at E_0 = +-e_3", ok, f"S {rs[0]} {rs[1]}, A {ra[0]} {ra[1]} " f"[{time.time() - t0:.1f} s]") if shrink != 1.0: print(f"{sum(CHECKS)}/{len(CHECKS)} checks passed") sys.exit(0 if all(CHECKS) else 1) print("C2 census classes at d = 3, 4, m = 2..4", flush=True) for d in (3, 4): N = Num(d) box = Box(N) bad = [] for m in (2, 3, 4): out = census(N, m) for lab in ("single", "multiple"): for key, f in out[lab]: w = lasso(f["coding"]) nodes, succ, nst = build(N, box, w) pr = paired(nodes, succ, nst, live_part(succ)) if (pr is not None) != (lab == "multiple"): bad.append((m, key, lab, pr)) check(f"d = {d}: paired iff labelled multiple, m = 2..4", not bad, f"{bad[:4]}") print(f"{sum(CHECKS)}/{len(CHECKS)} checks passed") sys.exit(0 if all(CHECKS) else 1) if __name__ == "__main__": main()