#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 ; }