# -*- coding: utf-8 -*-
# run_t2f.py  --  companion of verify_ca_type.py for Proposition 4.3 of
#   "Recycled roots under reduction and descent in the Casas-Alvero problem"
#   David Niedbala Giraudin, ORCID 0009-0009-1526-1178.
#   Technical assistance: Claude (Anthropic).
#
# The C program t2f.c below (embedded verbatim, since the archive accepts no .c file) checks,
# for each canonical sigma of degree n (first nonzero label 1), that q = 2**31 - 1 does not
# divide lc(E2) and that gcd(E1, E2) = 1 in F_q[y].  This script writes t2f.c to a temporary
# directory, compiles it with the system C compiler (cc or gcc) and runs it.
#
# Usage:  python3 run_t2f.py                 degrees 13..16 (a few seconds)
#         python3 run_t2f.py 17 18 19 20     any list of degrees (n = 20: minutes; 22: about
#                                            an hour; 24: many hours on one core)
#         python3 run_t2f.py 24 --part k m   only the k-th of m slices (0 <= k < m), so that
#                                            degree 24 can be spread over several machines.
#         python3 run_t2f.py --write         only write t2f.c to the current directory.
# For each degree it checks that the number of canonical sigma equals (3**(n-2) - 2**(n-1) + 1)/2
# (over all slices) and prints the FLAG lines.  The expected FLAG lines are none for n <= 21,
# three for n = 22 and five for n = 24; they are settled exactly by part [E] of verify_ca_type.py.
import hashlib, os, subprocess, sys, tempfile

T2F_SHA256 = "74828581f8bfd7d03093803696b9fe8621b2b55f7bb55ea07bab8adabfa7d581"
T2F_C = r'''/* t2f.c -- check that q = 2^31-1 does not divide N_n (type <= 2 scenarios), degree n.
 * Companion of verify_ca_type.py, same algorithm (Lemma 2.2 recursion, gcd in F_q[y]),
 * for the larger degrees.  Usage:  gcc -O3 -o t2f t2f.c ;  ./t2f n [part nparts]
 * Output: one FLAG line per sigma (positions m = 1..n-1, labels 0,1,2) for which
 * q divides lc(E2) or gcd(E1,E2) != 1 mod q, then a RESULT line with the count of
 * canonical sigma (first nonzero label = 1).  A FLAG is settled by an exact gcd over Q
 * (Res(E1,E2) != 0); it then only means q | Res(E1,E2), not N_n = 0.
 */
// Casas-Alvero, scenarios a racines recyclees {0,1,y} (t=2), caracteristique 0 via F_Q.
// sigma(m) in {0,1,2} pour m=1..d-2 ; sigma(d-1)=0 (racine multiple a l'origine) ; a_d=0.
// Canonique : le premier sigma non nul vaut 1 (le miroir 1<->y couvre l'autre moitie).
// Certificat : Q ne divise pas lc(E2) (verifie explicitement) et pgcd(E1,E2) trivial dans F_Q[y].
#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include <string.h>
#ifndef Q
#define Q 2147483647ULL
#endif
typedef uint64_t u64;
static int d; static u64 C[64][64];
static u64 a[64][64];      // a[m][k], degre <= m
static int sg[64];
static long long nleaf=0, nflag=0, nlc=0;
static int prefix_lo, prefix_hi, prefix_len; // partage du travail
static inline u64 red(u64 x){ if(Q==2147483647ULL){ x=(x&Q)+(x>>31); x=(x&Q)+(x>>31); return x>=Q?x-Q:x; } return x%Q; }
static inline u64 mulm(u64 x,u64 y){ return red(x*y); }
static u64 powm(u64 b,u64 e){u64 r=1;while(e){if(e&1)r=mulm(r,b);b=mulm(b,b);e>>=1;}return r;}
static int trimdeg(u64*p,int n){ while(n>0 && p[n]==0) n--; return n; } // degre (0 si constante)
static int gcddeg(u64*A0,int da,u64*B0,int db){
  // Euclide sans inverse : A <- lc(B)*A - lc(A)*x^sh*B (meme pgcd sur un corps)
  u64 A[64],B[64]; memcpy(A,A0,8*(da+1)); memcpy(B,B0,8*(db+1));
  da=trimdeg(A,da); db=trimdeg(B,db);
  int za=(da==0&&A[0]==0), zb=(db==0&&B[0]==0);
  if(za&&zb) return 99; if(za) return db; if(zb) return da;
  u64 *X=A,*Y=B; int dx=da, dy=db;
  for(;;){
    if(dx<dy){ u64*T=X;X=Y;Y=T; int t=dx;dx=dy;dy=t; }
    if(dy==0) return 0;               // Y constante non nulle
    u64 ly=Y[dy];
    while(dx>=dy){
      u64 lx=X[dx]; int sh=dx-dy;
      for(int k=0;k<sh;k++) X[k]=mulm(X[k],ly);
      for(int k=0;k<dy;k++) X[sh+k]=red(mulm(X[sh+k],ly)+Q-mulm(lx,Y[k]));
      X[dx]=0;
      dx--; while(dx>0 && X[dx]==0) dx--;
      if(dx==0 && X[0]==0) return dy;   // Y divise X : pgcd = Y
    }
  }
}
static void leaf(void){
  nleaf++;
  u64 E1[64],E2[130]; memset(E1,0,sizeof E1); memset(E2,0,sizeof E2);
  for(int j=0;j<=d;j++) for(int k=0;k<=j;k++){ u64 v=a[j][k]; if(!v) continue;
     E1[k]=(E1[k]+v)%Q; E2[k+d-j]=(E2[k+d-j]+v)%Q; }
  if(E2[d]==0){ nlc++; goto flag; }
  { int g=gcddeg(E1,d,E2,d); if(g==0) return; }
 flag:
  nflag++;
  printf("FLAG d=%d sigma=",d); for(int m=1;m<d;m++) printf("%d",sg[m]); printf("\n"); fflush(stdout);
}
static void compute(int m,int s){
  for(int k=0;k<=m;k++) a[m][k]=0;
  if(s==0) return;
  for(int j=0;j<m;j++){ u64 c=C[d-j][m-j]; int sh=(s==2)? m-j:0;
    for(int k=0;k<=j;k++){ u64 v=a[j][k]; if(v) a[m][k+sh]=(a[m][k+sh]+Q-mulm(c,v))%Q; } }
}
static void dfs(int m,int used1,int used2){
  if(m==d-1){ if(!(used1&&used2)) return; compute(d-1,0); sg[d-1]=0; for(int k=0;k<=d;k++)a[d][k]=0; leaf(); return; }
  for(int s=0;s<3;s++){
    if(s==2 && !used1) continue;           // canonique
    if(m<=prefix_len){ // filtre de partage : code du prefixe en base 3
    }
    sg[m]=s; compute(m,s); dfs(m+1,used1||s==1,used2||s==2);
  }
}
int main(int argc,char**argv){
  d=atoi(argv[1]);
  int part=argc>2?atoi(argv[2]):0, nparts=argc>3?atoi(argv[3]):1;
  for(int n=0;n<64;n++){C[n][0]=1;for(int k=1;k<=n;k++)C[n][k]=(C[n-1][k-1]+(k<n?C[n-1][k]:0))%Q;}
  for(int k=0;k<=d;k++)a[0][k]=0; a[0][0]=1;
  // partage : on enumere les 3 premiers choix (m=1..L) et on garde ceux d'indice = part mod nparts
  int L = (d-2>=6)?6:d-2; long idx=0;
  int tot=1; for(int i=0;i<L;i++) tot*=3;
  for(int code=0;code<tot;code++){
    int c=code, ok=1, u1=0,u2=0;
    for(int m=1;m<=L;m++){ int s=c%3; c/=3; if(s==2&&!u1){ok=0;break;} sg[m]=s; compute(m,s); u1|=(s==1); u2|=(s==2);}
    if(!ok) continue;
    if((idx++)%nparts!=part) continue;
    if(L==d-2){ if(u1&&u2){ compute(d-1,0); sg[d-1]=0; for(int k=0;k<=d;k++)a[d][k]=0; leaf(); } }
    else dfs(L+1,u1,u2);
  }
  printf("RESULT d=%d part=%d/%d Q=%llu leaves=%lld flagged=%lld lc0=%lld\n",d,part,nparts,(unsigned long long)Q,nleaf,nflag,nlc);
  return 0;
}
'''

def main():
    args = sys.argv[1:]
    src = T2F_C
    h = hashlib.sha256(src.encode()).hexdigest()
    print("t2f.c sha256 %s %s" % (h, "OK" if h == T2F_SHA256 else "MISMATCH"), flush=True)
    if h != T2F_SHA256:
        sys.exit(1)
    if args == ["--write"]:
        open("t2f.c", "w").write(src); print("written t2f.c"); return
    part, nparts = 0, 1
    if "--part" in args:
        i = args.index("--part"); part, nparts = int(args[i + 1]), int(args[i + 2]); del args[i:i + 3]
    degrees = [int(a) for a in args] or [13, 14, 15, 16]
    tmp = tempfile.mkdtemp()
    cfile, exe = os.path.join(tmp, "t2f.c"), os.path.join(tmp, "t2f")
    open(cfile, "w").write(src)
    for cc in ("cc", "gcc", "clang"):
        try:
            subprocess.run([cc, "-O3", "-o", exe, cfile], check=True); break
        except (OSError, subprocess.CalledProcessError):
            continue
    else:
        print("no working C compiler found"); sys.exit(1)
    for n in degrees:
        out = subprocess.run([exe, str(n), str(part), str(nparts)], capture_output=True, text=True, check=True).stdout
        flags = [l for l in out.splitlines() if l.startswith("FLAG")]
        res = [l for l in out.splitlines() if l.startswith("RESULT")][-1]
        leaves = int(res.split("leaves=")[1].split()[0])
        expected = (3 ** (n - 2) - 2 ** (n - 1) + 1) // 2
        for l in flags:
            print("  " + l)
        if nparts == 1:
            print("n=%2d : %d canonical sigma (expected %d: %s), flagged %d"
                  % (n, leaves, expected, "OK" if leaves == expected else "WRONG", len(flags)), flush=True)
        else:
            print("n=%2d slice %d/%d : %d canonical sigma, flagged %d (add the counts of all slices; total expected %d)"
                  % (n, part, nparts, leaves, len(flags), expected), flush=True)

if __name__ == "__main__":
    main()