#include<stdio.h>
#include<stdlib.h>
#include<time.h>

/* this trivial program computes A225371(n) for n=1,2,3,4,5,6.
   
   There are 2^(N^2) binary matrices, which can be numbered from
   0 to 2^(N^2)-1 by simple concatenating the bits in the matrix.

   I allocate an array Q of 2^(N^2) bytes (or bits, for N=6) initially
   set to zero. I generate each binary matrix, I square it, I convert
   the square to a number x and I set the w-th element of Q to 1.
    At the end, I count the number of ones in Q, which is the total number
    of matrices which are squares.
    Since the running time in the largest case (N=6) was satisfactory
    (3.5 hours on my i7-3930k) I did not try to further optimize the
    program, which, by the way, should not be considered an example of 
    good programming. Actually it took my more time to write these comments
    than to write the program itself. Giovanni Resta 2013. */

// note, I use the std99 extensions of C that allow to 
// define local variables in for : like for(int i=0;.... )

// note, for N =  6 the programs needs about 8 GB of RAM.

#define N 6


// this is used to speed up the count of bits in Q
#if N  == 6
unsigned ByteCnt[256];
inline unsigned BitSum(unsigned int x)
{
  return ByteCnt[x&255]+ByteCnt[(x>>8)&255]+
        ByteCnt[(x>>16)&255]+ByteCnt[x>>24];
}
#endif

int main()
{
  // the matrix to be squared
  unsigned short int a[N][N];
  // Q[x]=1 tells me that the x-th matrix is a square
  // if N=6 I use  single bits (packed in 32 bit unsigned integers)
  // to store 
#if N == 6
  unsigned int *Q;
  Q = calloc(1LLU<<(N*N-5),sizeof(unsigned int));
  if (!Q) return printf("I need 8GB !\n");
#else
  unsigned char *Q;
  Q = calloc(1U<<(N*N),sizeof(unsigned char));
#endif

  // init of BitCnt array
  // used to count bits in Q when N=6
#if N == 6
  time_t tt=time(0);
  for (int i=0; i<256; i++)
    { 
      int j;
      ByteCnt[i]=0;
      for (j=0; j<8; j++) ByteCnt[i]+=(i&(1<<j)?1:0); 
    }
#endif

  // I repeat for each of the 2^(N*N) matrices
  for (unsigned long u=0; u<1LLU<<(N*N); u++)
    { 
      unsigned long v=u, w=0;
      // if N = 6 I print some stats otherwise I become impatient
#if N == 6
      if ((u & ((1LLU<<28)-1))==0 && u)
	printf(">> %.2f%% in %d sec. \n",(100.0*u)/(1LLU<<36),
	       (int)(time(0)-tt));
#endif 
      //============================
      // I build the matrix a corresponding to index u
      for (int i=0; i<N; i++)
	for (int j=0; j<N; j++) { a[i][j] = v & 1;  v>>=1; }
      //============================
      // simple matrix product, with And instead of * and
      // Xor instead of +, to respect GF(2) math.
     for (int i=0; i<N; i++)
	for (int j=0; j<N; j++)
	  {
	    unsigned short int Acc=0;
	    for (int k=0; k<N; k++) Acc ^= a[i][k] & a[k][j];
	    // as soon as an entry of the result is computed,
	    // I update w which at the end will contain
	    // a representation of the square
	    w = (w<<1) | (Acc&1);
	  }
      //============================
     // I set to 1 the w-th byte (or bit) of Q
#if N == 6
     Q[w>>5]|=1u<<(w&31);
#else 
     Q[w]=1;
#endif
  }
  // I count the bytes (or bit) of Q set to 1 and that is the result.
  long int Tot = 0;
#if N == 6
  for (long int k=0; k<1LLU<<(N*N-5); k++)  Tot+=BitSum(Q[k]);
#else
  for (int k=0; k<1LU<<(N*N); k++)  Tot+=Q[k];
#endif
  printf("\nA225371(%d) =  %ld\n",N,Tot);

  return 0;
}