A-statistics code

// computes the 3D Moser sets and sorts them by a-statistic.

1. include <stdio.h>
2. include <stdlib.h>
3. include <string.h>
4. include <time.h>

1. define NUM_LINES_3D 25 // 25 non-horizontal lines in [3]^3

// the bitmasks of the 25 lines in [3]^3 const long line_bitmasks[NUM_LINES_3D] = {262657l,263172l,266304l,270592l,525314l,532608l,1049601l,1050628l,1056832l,1065216l,2101256l,2105376l,4202512l,8396808l,8405024l,16781313l,16785412l,16810048l,16843008l,33562626l,33620096l,67117057l,67125252l,67174464l,67240192l};

/* stat[i] is 0,1,2,3 depending on whether i (represented base 3) is of type a,b,c,d. */

const int stat[27] = {0,1,0, 1,2,1, 0,1,0,

                     1,2,1, 2,3,2, 1,2,1,

0,1,0, 1,2,1, 0,1,0};

1. define NUM_PERMS 384

int perms[NUM_PERMS][16]; // The 48 symmetries of the 4D cube long pow2[27]; // The 27 powers of 2

void buildperms(void) {

  int i,j,k,l,p,x,y,z,w,a,b,c,d;
  int n=0;
  long power=1;

  for (i=0; i<27; i++) { pow2[i] = power; power *= 2; }

  for (x=0;x<1;x++)

for (y=0;y<2;y++) for (z=0;z<2;z++)

 		      for (w=0;w<2;w++)

for (p=0;p<24;p++) { for (i=0;i<2;i++) for (j=0;j<2;j++) for (k=0;k<2;k++)

                                for (l=0;l<2;l++)

{ switch (p) { case 0: a=i; b=j; c=k; d=l; break; case 1: a=i; b=j; c=l; d=k; break; case 2: a=i; b=k; c=j; d=l; break; case 3: a=i; b=k; c=l; d=j; break; case 4: a=i; b=l; c=k; d=j; break; case 5: a=i; b=l; c=j; d=k; break;

case 6: a=j; b=i; c=k; d=l; break; case 7: a=j; b=i; c=l; d=k; break; case 8: a=j; b=k; c=i; d=l; break; case 9: a=j; b=k; c=l; d=i; break; case 10: a=j; b=l; c=k; d=i; break; case 11: a=j; b=l; c=i; d=k; break;

case 12: a=k; b=j; c=i; d=l; break; case 13: a=k; b=j; c=l; d=i; break; case 14: a=k; b=i; c=l; d=j; break; case 15: a=k; b=i; c=j; d=l; break; case 16: a=k; b=l; c=i; d=j; break; case 17: a=k; b=l; c=j; d=i; break;

case 18: a=l; b=j; c=k; d=i; break; case 19: a=l; b=j; c=i; d=k; break; case 20: a=l; b=k; c=j; d=i; break; case 21: a=l; b=k; c=i; d=j; break; case 22: a=l; b=i; c=k; d=j; break; case 23: a=l; b=i; c=j; d=k; break; } if (x) a = 1-a; if (y) b = 1-b; if (z) c = 1-c; if (w) d = 1-d; perms[n][i+2*j+4*k+8*l] = a+2*b+4*c+8*d; } n++; } }

long perm(long set, int i) // permutes set by permutation i { long in,out; int j;

in = set; out = 0; for (j=0; j < 16; j++) { if (in & 1) out |= pow2[perms[i][j]]; in >>= 1; }

return out; }

const int amask[8] = {0, 2, 6, 8, 18, 20, 24, 26}; // the 8 a-points of [3]^3

void init_data(void) {

int i,j,k;

   /* enumerate the 230 line-free sets in [3]^2 by brute force */

   int moser_2d[512];

int lincount=0;

   for (i=0; i < 512; i++)
    {
     moser_2d[i] = 0;

if ((i & 7) == 7) continue; if ((i & 56) == 56) continue; if ((i & 73) == 73) continue; if ((i & 84) == 84) continue; if ((i & 146) == 146) continue; if ((i & 273) == 273) continue; if ((i & 292) == 292) continue; if ((i & 448) == 448) continue; moser_2d[i] = 1; lincount++; }

printf("Number of 2D Moser sets: %d\n", lincount);

   /* Now enumerate all line-free sets in [3]^3 by brute force. */
   long set, tmpset;

long slice1, slice2, slice3; int s[4]; long num_linefree=0; int a;

   printf("Computing 3D Moser sets...\n");

   long statcounts[256];

for (i=0; i < 256; i++) statcounts[i] = 0;

for (slice1 = 0; slice1 < 512; slice1++) if (moser_2d[slice1]) for (slice2 = 0; slice2 < 512; slice2++) if (moser_2d[slice2]) for (slice3 = 0; slice3 < 512; slice3++) if (moser_2d[slice3]) { set = slice1 | (slice2 << 9) | (slice3 << 18);

if ((set & line_bitmasks[0]) == line_bitmasks[0]) continue; if ((set & line_bitmasks[1]) == line_bitmasks[1]) continue; if ((set & line_bitmasks[2]) == line_bitmasks[2]) continue; if ((set & line_bitmasks[3]) == line_bitmasks[3]) continue; if ((set & line_bitmasks[4]) == line_bitmasks[4]) continue; if ((set & line_bitmasks[5]) == line_bitmasks[5]) continue; if ((set & line_bitmasks[6]) == line_bitmasks[6]) continue; if ((set & line_bitmasks[7]) == line_bitmasks[7]) continue; if ((set & line_bitmasks[8]) == line_bitmasks[8]) continue; if ((set & line_bitmasks[9]) == line_bitmasks[9]) continue; if ((set & line_bitmasks[10]) == line_bitmasks[10]) continue; if ((set & line_bitmasks[11]) == line_bitmasks[11]) continue; if ((set & line_bitmasks[12]) == line_bitmasks[12]) continue; if ((set & line_bitmasks[13]) == line_bitmasks[13]) continue; if ((set & line_bitmasks[14]) == line_bitmasks[14]) continue; if ((set & line_bitmasks[15]) == line_bitmasks[15]) continue; if ((set & line_bitmasks[16]) == line_bitmasks[16]) continue; if ((set & line_bitmasks[17]) == line_bitmasks[17]) continue; if ((set & line_bitmasks[18]) == line_bitmasks[18]) continue; if ((set & line_bitmasks[19]) == line_bitmasks[19]) continue; if ((set & line_bitmasks[20]) == line_bitmasks[20]) continue; if ((set & line_bitmasks[21]) == line_bitmasks[21]) continue; if ((set & line_bitmasks[22]) == line_bitmasks[22]) continue; if ((set & line_bitmasks[23]) == line_bitmasks[23]) continue; if ((set & line_bitmasks[24]) == line_bitmasks[24]) continue;

/* compute the a mask */ a = 0; for (j=0; j<8; j++) if (set & pow2[amask[j]]) a |= pow2[j];

/* now compute the a,b,c,d stats */ s[0]=s[1]=s[2]=s[3]=0; tmpset = set; for (j=0; j < 27; j++) { if (tmpset & 1l) s[stat[j]]++; tmpset /= 2; }

statcounts[a]++; num_linefree++; if (num_linefree % 10000 == 0) printf("*");

}

for (j=0; j<256; j++) printf("%ld sets with a-statistic %o\n",statcounts[j],j); printf("\nNumber of 3D Moser sets: %ld\n", num_linefree);

}

int main() {

  int i,j,k,l,m;

  buildperms();
  init_data();

}