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

/* compute number of nxn normal binary matrices
 as listed in A055547 of http://research.att.com/~njas/sequences .
 n  0     1      2        3        4       5       6         7
 ---------------------------------------------------------------
          2      8       68     1124    36112   2263268          ?

 Note that
 A055547(n)>= 2^[n(n+1)/2] (because all symmetric matrices
 are normal)

 and A055547(n)>=A055547(n-1)*2^n because we can generate a new normal
 n by n matrix from a normal (n-1)x(n-1) matrix by attaching any
 binary vector of length n at the same to the right (new column) and
 below (new row) to the old matrix, and because there are 2^n binary
 vectors available.

 Any binary normal matrix can be split into a symmetric and a antisymmetric
 piece [ H=A^s + A^a where A^s=(A+A')/2 and A^a=(A-A')/2, and A' is the
 transposed matrix].
 The number of different antisymmetric matrices of this set is
 counted in A007081(n).
 Richard J. Mathar, 03/22/2006
*/


/* define this if the matrices are also to be shown, not just counted */
/* #define SHOWMATS */

/* this indicates that we're using a binary basis... don't change */
#define SIZ 2

/* maximum matrix side length */
#define NMAX 10

/* a flattened index of column 'c' and row 'r' in the table of
* size nVec by nVec
*/
#define FLAT(r,c) ((c)+(r)*nVec)

/**
* This function computes the dot product of vector no 'i' and vector no 'j'
* which are both stored as binV[i][0..n-1] and binV[j][0..n-1].
* Since these are binary vectors of length n, the outcome is anywhere from
* 0 to n (inclusive).
*/
short int overl(const short (*binV)[NMAX],const int n, const int i, const int j)
{
	short resul = 0 ;
	for(int c=0; c < n; c++)
		resul += binV[i][c]*binV[j][c] ;

#ifdef DEBUG
	printf("overlap ") ;
	for(int c=0; c < n; c++)
		printf("%2d ",binV[i][c]) ;
	printf(" with ") ;
	for(int c=0; c < n; c++)
		printf("%2d ",binV[j][c]) ;
	printf(" is %d\n",resul) ;
#endif
	return resul ;
}

/**
* This function initializes the list of all binary vectors of length n (with n bits)
* and accumulates them into a matrix which is returned.
* The implementation is simply to count from 0 to 2^n-1 and decompose these
* numbers into their binary version by consecutive division modulo 2.
* On return, the first element of the vector is resul[0][0..n-1]=000000..0,
* the second is resul[1][0..n-1]=0000000..01,
* the third is resul[2][0..n-1]=0000000..10,
* the first but last is resul[n-2][0..n-1] = 11111...110,
* and the last is resul[n-1][0..n-1] = 11111...111.
*/
short (*initBin(const int n))[NMAX]
{
	const int nVec = 1<< n ;
	short (*resul)[NMAX] = malloc( nVec*sizeof(short[NMAX])) ;

	for(int no=0; no < nVec; no++)
	{
		int noShiftd = no ;
		for(int bit=0; bit <n; bit++)
		{
			resul[no][bit] = noShiftd % 2 ;
			noShiftd >>= 1 ;
		}
	}

#ifdef DEBUG
	printf("%d vectors\n",nVec) ;
	for(int no=0; no < nVec ; no++)
	{
		for(int bit=0; bit<n; bit++)
			printf("%2d",resul[no][bit]) ;
		printf("\n") ;
	}
	printf("\n") ;
#endif

	return resul ;
}

/**
* This initializes the symmetric matrix of the dot products between any two of the
* binary vectors (of which we have nVec=2^n, already allocated in binV).
* The dot product between binary vector i and binary vector j is stored in
* resul[FLAT(i,j)], where FLAT is a macro to simplify handling of  the
* dynamically allocated array which is returned. 0<=i,j<nVec.
* The number of entries in this matrix for a specific dot product between 0
* and n is given by A027465.
*/
short * initCoupl(const int n, const int nVec, const short (*binV)[NMAX])
{
	short *resul = (short *)malloc( nVec*nVec*sizeof(short)) ;
	for(int row=0; row< nVec; row++)
	{
		for(int col=row; col< nVec; col++)
			resul[FLAT(row,col)] = resul[FLAT(col,row)] = overl(binV,n,row,col) ;
	}

#ifdef DEBUG
	printf("overlaps\n") ;
	for(int i=0; i < nVec ; i++)
	{
		for(int j=0; j<nVec; j++)
			printf("%d ",resul[FLAT(i,j)]) ;
		printf("\n") ;
	}
	printf("\n") ;
#endif

	return resul ;
}

/**
* This program attempts to fill the n by n matrix mat[][], which is a partially
* filled candidate for a binary normal matrix, at level 'levl', where the matrix is
* filled from the top to bottom and from the left to right, alternatively with rows
* and columns. [This allows early skipping of the cases where the dot product
* of a pair of columns does not match the dot product of the pair of rows with the
* same pair index.]
* The number of matrices generated so far is provided in 'had' and updated on return.
* 'n' is the side length of the binary matrix.
* 'nVec' is 2^n, the number of binary vectors stored in 'binvec'.
* 'coupl' is the overlap matrix between pairs of vectors stored in 'binvec'
* 'vecLis' contains the indices into 'binvec': vecLis[0] the index for the uppermost
* row, vecLis[1] for teh leftmost column, vecLis[2] for the 2nd row from above, etc.
* On return, vecLis[levl] has been set.
*/
long long fill(short mat[NMAX][NMAX], const int n, int levl, const int nVec, const short *coupl, const short (*binvec)[NMAX],
	int *vecLis, long long had)
{
	/* This boolean indicator says we're filling a column if 'levl' is odd, else a row */
	const int isCol = levl % 2 ;

	/* This runs from 0 to n-1 for rows and columns */
	const int RowCol = levl / 2 ;

	/* For this particular row or column we check in turn all binary vectors whether they
	* can be placed in  here.
	*/
	for(int v=0; v < nVec ; v++)
	{
		/* This boolean indicates whether the vector no 'v' can be placed into
		* the current matrix as a row or column without violating the constraints
		* set by creating a normal matrix.
		*/
		int works = 1 ;

		/* bookkeeping that we are now trying to insert binary vector with indes 'v' */
		vecLis[levl] = v ;

		/* Constraints exist unless this is the topmost row (which is inserted into the empty matrix) */
		if ( levl > 0)
		{
			/* For new columns check that the overlap with itself equals the
			* overlap of the row vector with the same index.
			*/
			if( isCol )
			{
				/* find out which binary vector had been placed in the matrix row
				* with the equivalent index as the current column. Since 'levl' switches
				* between rows and colums as these are inserted, we only have to 
				* move back by one index.
				*/
				const int matchr = levl-1 ;
				const int matchrV = vecLis[matchr] ;

				/* we read of the dot product of the associated row with itself */
				const int rowC = coupl[FLAT(matchrV,matchrV)] ;

				/* If this differs from the dot product of the attempted column,
				* the result would not be normal and we move on to the next 'v'.
				*/
				if (coupl[FLAT(v,v)] != rowC)
					continue ;
			}

			/* Check that the selected binary vector matches with the first digits
			* with what is already in the matrix. The vectors lower than 'levl' that
			* have already been inserted into mat[][] already occupy the upper rows
			* and left columns.
			*/
			if ( isCol)
			{
				/* If the new entry is a column, we need to compare its first
				* vector elements with the existing top rows at the equivalent place.
				* If there is a mismatch, the vector no 'v' cannot be inserted.
				* We set 'works' to 'false' to indicate that we should not continue 
				* with this setting.
				*/
				for(int row=0 ; row <= RowCol; row++)
				{
					if ( binvec[v][row] != mat[row][RowCol] )
					{
						works = 0 ;
						break ;
					}
				}
			}
			else
			{
				/* If the new entry is a row, we need to compare its first
				* vector elements with the existing left columns at the equivalent place.
				* If there is a mismatch, the vector no 'v' cannot be inserted.
				*/
				for(int col=0 ; col < RowCol; col++)
					if ( binvec[v][col] != mat[RowCol][col] )
					{
						works = 0 ;
						break ;
					}
			}

			if ( isCol)
			{
				/* check that the overlap with all previous rows/columns matches the ones
				* of the columns/rows. If this is for example the third column from the left,
				* we check that the dot product with the first column equals the dot
				* product between the third row and the first row, and that the dot product
				* with the 2nd column equals the dot product between the third row and the 2nd row.
				* The dot product between column i and column j must be the same as the
				* dot product between row i and row j for all pairs i and j (including i=j)
				* to create a normal matrix.
				*/
				for(int p =0 ; p < RowCol ; p++)
				{
					const int matchrVi = vecLis[2*p] ;
					const int matchrVj = vecLis[levl-1] ;
					const int matchcVi = vecLis[2*p+1] ;
					const int matchcVj = vecLis[levl] ;
					if ( coupl[FLAT(matchcVi,matchcVj)] != coupl[FLAT(matchrVi,matchrVj)] )
					{
						works = 0 ;
						break ;
					}
				}
			}
		}
		/* if one of the checks above failed to insert a row or column 
		* that matches the criteria, continue with the next binary vector
		* at the same row or column position.
		*/
		if ( ! works)
			continue ;

		/* other wise the new row or column fits, and we store it into the 
		* 'living' instance of the binary matrix.
		*/
		if ( isCol)
		{
			for(int row=0; row < n; row++)
				mat[row][RowCol] = binvec[v][row] ;
		}
		else
		{
			for(int col=0; col < n; col++)
				mat[RowCol][col] = binvec[v][col] ;
		}

		/* If we have not yet inserted a full number of n rows and n columns, we
		* call this function recursively to add another row or column.
		* The fine point is that this is also called for levl==2n-2, altthough
		* at this point the matrix mat[][] is already filled: This allows to have
		* the check for the rightmost column be done on the same level as
		* all the other columns. 
		*/
		if (levl < 2*n-1 )
		{
			had = fill(mat,n,levl+1,nVec,coupl,binvec,vecLis,had) ;
		}
		else
		{
			/* If otherwise we have a full matrix with no discrepancies with
			* the checks for a normal matrix, we raise the count by 1 and
			* optionally print it on output.
			*/
			had++ ;
#ifdef SHOWMATS
			for(int row=0 ;row < n; row++)
			{
				for(int col=0 ; col < n ; col++)
					printf("%2d",mat[row][col]) ;
				printf("\n") ;
			}
#endif
			/* for each 10 thousandths display what we've reached so far */
			if ( had % 10000 == 0)
				printf("%d\n",had) ;
		}
	}
	return had ;
}

/**
* Main entry. The program is called with a single parameter, 1<=n<=NMAX
*/
int main(int argc, char *argv[])
{
	/* the number of matrices we've found so far */
	long long have=0 ;

	/* the list of the binary vectors with n ditigs, 2^n of them */
	short (*binvec)[NMAX]=NULL ;

	/* the matrix with the 'overlaps" (=dot products) between pairs of the binary vectors */
	short *coupl=NULL ;

	/* an instance of a candidate for a normal binary matrix */
	short mat[NMAX][NMAX] ;

	/* vecLis[0] contains the index into binvec[] for the digits of the first (uppermost) row,
	* vecLis[1] the index for the digits of the first (leftmost) column, vecLis[2] the
	* index for the digits of the 2nd row and so on, up to vecLis[2n-1].
	*/
	int vecLis[2*NMAX] ;
	int n, nVec ;
	if ( argc != 2)
	{	fprintf(stderr,"usage: %s n\n",argv[0]) ;
		return 1;
	}

	/* get the command line argument */
	n = atoi(argv[1]) ;

	/* compute 2^n as the number of binary vectors of length n */
	nVec = 1<<n ;

	/* initialize the auxiliary lookup arrays with the binary vectors
	* and the array of their dot products
	*/
	binvec = initBin(n) ;
	coupl = initCoupl(n,nVec,binvec) ;

	/* call the main program which will recursively fill the first row,
	* the second column, the second row, the second column etc alternatively
	* to generate the normal matrices.
	*/
	have = fill(mat,n,0,nVec,coupl,binvec,vecLis,0) ;
	printf("n=%d %lld\n",n,have) ;

	free(binvec) ;
	free(coupl) ;
	return 0 ;
}