/*
 *  File:    DLScanonizer.cpp
 *  Author:  S Harry White
 *  Created: 2020-08-01
 *  Updated: 2020-09-27
 *    Fix open_input when file name includes .txt.
 *  Updated: 2020-10-02
 *    Fix reading of numbers > 20.
 *  Updated: 2021-05-11
 *    Fix allocation of CF store.
 *  Updated: 2021-12-29
 *    Fix potential problem with use of | vs ||.
 * Updated: 2023-01-27
 *   Change outFile from bufSize to outSize.
 */

/*
 *  Converts DLS to canonical form. 
 *  Optionally the format is 1: natural order first row
 *                        or 2: natural order \diagonal
 *
 * All 8 aspects of all valid row/column swapping transformations of the DLS are searched.
 *
 * For format 1, the form is the smallest.
 *
 * For format 2, the form is the smallest that has the smallest /diagonal.
 *
 *             (n-2)/2
 *  There are 2       (n/2)! transformations for order N,
 *  where n=N for even orders and N-1 for odd orders:
 *
 *    order   4 5    6 7    8 9    10 11    12 13     14 15
 *           -----------------------------------------------
 *   xforms    4      24    192    1,920    23,040   322,560 
 *
 *  If more than one DLS canonical form results, they are sorted in
 *  ascending numerical order removing duplicates.
 */

#include "stdafx.h"
#include <conio.h>
#include <errno.h>
#include <fcntl.h>
#include <io.h>
#include <share.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <sys\stat.h>
#include <time.h>
#include <Windows.h>

typedef unsigned char CFint;
typedef unsigned long Uint;
//typedef unsigned long long Uint64; // for x64

const CFint maxCFint=UCHAR_MAX;
const int minN=4, maxN=15,            // supported orders
          nRot=8,                     // aspects
          CFintBytes=1,
          startNumBlocks=1,
          CFsPerBlock=100000,
          maxBlocksIncr=64,
          pointerSize=sizeof(void *);
const Uint maxCFs=MAXUINT/pointerSize;

int  N, NN, G, M, O, // order, N*N, N/2, N-1, N+1
  nXforms,           // number of xforms
  CFsSize, CFsIndex, // number of canonical forms, index    
  Xindex,            // xforms index
  ROindex,           // row orders index
  Qindex, QueueSize, // Queue of row/column numbers to swap
  **nfr_1_N=NULL, n1_nfr_1_N, n2_nfr_1_N, // for format 1
    // the numbers at 1 and N in a xform, followed by the indexes in X of all xforms
    // that have the same numbers at 1 and N
  **nfr_2_N=NULL, n1_nfr_2_N, n2_nfr_2_N, // for format 1
    // the numbers at 2 and N in a xform, followed by the indexes in X of all xforms
    // that have the same numbers at 2 and N
  **nbd_M_O=NULL, n1_nbd_M_O, n2_nbd_M_O, // for format 2
    // the numbers at M and O in a xform, followed by the indexes in X of all xforms
    // that have the same numbers at M and O
  smallRead, bigRead, dlsNumber,          // input
  numBlocksAllocated, CFsize,             // store allocated
  dlsBytes, nDLSs, writeSize, outDLSs;      // output
  

struct t_xform { CFint w, x; }; t_xform *Queue=NULL;
void **blocks=NULL;
CFint **CFs=NULL,          // canonical forms store
      **A=NULL,            // the 8 aspects of the input DLS
      **rot=NULL,          // the 7 other aspects of 0 .. NN-1
      **aDLS=NULL,         // 0 .. NN-1, changed to enter each xform into X
      **rowOrders=NULL,    // used in getting xforms
       *rowOrderTmp=NULL,  // used in getting xforms
       *swapCounts=NULL,   // used in getting xforms
       *sDLS=NULL,         // used in finding smallest form
       *iDLS=NULL,         // input DLS
       *X=NULL;            // the xforms of 0 .. NN-1

char *outBuffer=NULL, *outPointer=NULL;

const bool F=false, T=true;
bool valueFree[maxN], // used in getting xforms
     Odd,             // odd order
     zeroBase;        // input DLS range is 0 to N-1

//void outLocalTime() {
////   ------------
//  time_t startTime=time(NULL);
//  struct tm local; localtime_s(&local, &startTime); char dateTime[100];
//  size_t slen=strftime(dateTime, 100, "%A %Y-%m-%d %X %Z\n\0", &local);
//  printf(dateTime);
//} // outLocalTime

void printElapsedTime(time_t startTime, const char *c) {
//   ---------------- 
  int et=(int)difftime(time(NULL), startTime);
  if (et>0) printf("%selapsed time %d:%02d:%02d\n", c, et/3600, et%3600/60, et%60);
} // printElapsedTime
//====================================== allocate store ============================================

char *storeAllocFail="Storage allocation failed"; FILE *wfpE=NULL;
void reportError(char *msg) { fprintf(wfpE,"\a\nError: %s.\n", msg); }
//   -----------

void freeIntArray(int ***a, const int n1) {
//   ------------
  if (*a!=NULL) { for (int i=0; i<n1; i++) free((*a)[i]); free(*a); *a=NULL; }
} // freeIntArray

void freeCFs(const int allocated) {
//   -------
  if (blocks!=NULL) {
    if ((allocated!=0)) {
      for (int i=0; i<allocated; i++) free(blocks[i]);
      if (CFs!=NULL) { free(CFs); CFs=NULL; CFsSize=0; }
    }
    free(blocks); blocks=NULL;
  }
  numBlocksAllocated=0;
} // freeCFs

void freeCFints(CFint **a) { if (*a!=NULL) { free(*a); *a=NULL; }}
//   ----------

void freeCFintArray(CFint ***a, const int n1) {
//   --------------
  if (*a!=NULL) { for (int i=0; i<n1; i++) free((*a)[i]); free(*a); *a=NULL; }
} // freeCFintArray

void freeChars(char **a) { if (*a!=NULL) { free(*a); *a=NULL; }}
//   ---------

void freeQueue() { if (Queue!=NULL) { free(Queue); Queue=NULL; }}
//   ---------

void freeStore() {
//   ---------
  freeQueue();
  freeIntArray(&nfr_1_N, n1_nfr_1_N); freeIntArray(&nfr_2_N, n1_nfr_2_N);
  freeIntArray(&nbd_M_O, n1_nbd_M_O); freeCFs(numBlocksAllocated);
  freeCFintArray(&rowOrders, nXforms); freeCFintArray(&aDLS, N);
  freeCFintArray(&A, nRot); freeCFintArray(&rot, nRot-1);
  freeChars(&outBuffer); freeCFints(&X); freeCFints(&swapCounts);
  freeCFints(&sDLS); freeCFints(&iDLS); freeCFints(&rowOrderTmp);
} // freeStore

bool allocateCFints(CFint **l, const int n) {
//   --------------
  *l=(CFint*)malloc(n*sizeof(CFint)); return *l!=NULL;
} // allocateCFints

bool allocateChars(char **l, const int n) {
//   -------------
  *l=(char*)malloc(n*sizeof(char)); return *l!=NULL;
} // allocateChars

bool allocateIntArray(int ***a, const int n1, const int n2) {
//   ----------------
  bool ok=((*a=(int**) malloc(n1*sizeof(int*)))!=NULL);
  if (ok) {
    int numAllocated=n1;
    for (int i=0; i<n1; i++) {
      int *p=(int*) malloc(n2*sizeof(int));    
      if (p==NULL) { numAllocated=i; ok=F; break; } (*a)[i]=p;
    }
    if (!ok) freeIntArray(a, numAllocated);
  }
  return ok;
} // allocateIntArray

bool allocateCFintArray(CFint ***a, const int n1, const int n2) {
//   ------------------
  bool ok=((*a=(CFint**) malloc(n1*sizeof(CFint*)))!=NULL);
  if (ok) {
    int numAllocated=n1;
    for (int i=0; i<n1; i++) {
      CFint *p=(CFint*) malloc(n2*sizeof(CFint));   
      if (p==NULL) { numAllocated=i; ok=F; break; } (*a)[i]=p;
    }
    if (!ok) freeCFintArray(a, numAllocated);
  }
  return ok;
} // allocateCFintArray

bool allocateCFs(int numNew) {
//   -----------
  bool ok=F; union ovly { Uint i; void *p; }; // for Win32
  //union ovly { Uint64 i; void *p; }; // for x64
  const int numOld=numBlocksAllocated;
  if ((numNew-numOld)>maxBlocksIncr) numNew=numOld+maxBlocksIncr;
  const Uint numNewCFs=numNew*CFsPerBlock; // Uint for malloc(numNewCFs*pointerSize)
  if (numNew!=startNumBlocks) printf("increase CFs store to %u\n", numNewCFs);
  if (numNewCFs<=maxCFs) {
    void **tmp=(void **)malloc(numNew*pointerSize);
    if (ok=(tmp!=NULL)) {
      int allocated=numNew;
	    for (int i=0; i<numOld; i++) tmp[i]=blocks[i];
      for (int i=numOld; i<numNew; i++) {
        void *p=malloc(CFsPerBlock*CFsize);
        if (p==NULL) { ok=F; allocated=i; break; } tmp[i]=p;
      }
	    freeCFs(0); blocks=tmp;
	    if (ok) {
        tmp=(void **)malloc(numNewCFs*pointerSize);
        if (ok=(tmp!=NULL)) {
          for (int i=0; i<CFsIndex; i++) tmp[i]=CFs[i]; free(CFs); CFs=(CFint**)tmp; int j=CFsIndex;
          for (int i=numOld; i<numNew; i++) {
            int k=0; void *p=blocks[i];
            while (k++<CFsPerBlock) {
              CFs[j++]=(CFint*)p; ovly ip; ip.p=p; ip.i+=NN*CFintBytes; p=ip.p;           
            }
          }
	        numBlocksAllocated=numNew; CFsSize=numNewCFs;
        }
	    } else freeCFs(allocated);
    }
  }
  if (!ok) reportError(storeAllocFail); return ok;
} // allocateCFs

bool allocateQueue() {
//   -------------
  Queue=(struct t_xform*)malloc(QueueSize*sizeof(struct t_xform)); return Queue!=NULL;
} // allocateQueue

bool allocateStore() {
//   -------------
  bool ok=allocateQueue();

  if (ok) ok=allocateCFs(startNumBlocks);
    if (ok) ok=allocateCFintArray(&rowOrders, nXforms, G);
      if (ok) ok=allocateCFintArray(&aDLS, N, N);
        if (ok) ok=allocateCFintArray(&A, nRot, NN);
          if (ok) ok=allocateCFintArray(&rot, nRot-1, NN);

            if (ok) ok=allocateChars(&outBuffer, writeSize);
              if (ok) ok=allocateCFints(&X, nXforms*NN);
                if (ok) ok=allocateCFints(&swapCounts, nXforms);
                  if (ok) ok=allocateCFints(&sDLS, NN);
                    if (ok) ok=allocateCFints(&iDLS, NN);
                      ok=allocateCFints(&rowOrderTmp, G);
  if (!ok) { freeStore(); reportError(storeAllocFail); }
  return ok;
} // allocateStore
//====================================== input ============================================

void clearLine(int c) { while (c!='\n') c=getchar(); }
//   ---------

int getInt(char *q) {
//  -------
  int n=0; printf("%s? ", q); scanf_s("%d", &n); clearLine(getchar()); return n;
}

bool getY() {
//   ----
  int c; do { c=getchar(); } while ((c==' ')|(c=='\t')|(c=='\n'));
  clearLine(c); return (c=='Y')|(c=='y');
}

bool getYorOrder(int *n) {
//   -----------
  bool result=F; int c; *n=0; do { c=getchar(); } while ((c==' ')|(c=='\t')|(c=='\n') );
  if ( (c=='Y')|(c=='y') ) result=T;
  else if ( (c!='N')&(c!='n') )
    if ( ('1'<=c)&(c<='9') ) {
      int i=c-'0'; while ( ('0'<=(c=getchar()))&&(c<='9') ) i=i*10+c-'0'; *n=i; result=T;
    }   
  clearLine(c); return result;
} // getYorOrder

const int bufSize=1024, outSize=bufSize+50;
bool getFileName(char *buf, const int size) {
//   -----------
  int c, i=0; char *s=buf; do { c=getchar(); } while ((c==' ')|(c=='\t')|(c=='\n'));
  *s=c; while (i++<size) if ((*++s=getchar())=='\n') break;
  if (*s!='\n') { printf("\nFile name too long.\n"); clearLine(*s); return F; }
  while ((*--s==' ')||(*s=='\t'));  /* strip trailing whitespace */ *++s='\0'; return T;
} // getFileName

void check_txt(char *buf, const int size) {
//   ---------
  char *s=buf; bool txt=F; while (*s++!='\0');
  while (--s!=buf) if (*s=='.') { txt=(*++s=='t')&&(*++s=='x')&&(*++s=='t')&&(*++s=='\0'); break; }
  if (!txt) strcat_s(buf, size, ".txt");
} // check_txt

FILE *openInput() {
//    ---------
  char buf[bufSize]; printf("File name? "); FILE *rfp=NULL;
  if (getFileName(buf, bufSize)) {
    check_txt(buf, bufSize);
    if (fopen_s(&rfp, buf, "r")!=0) {
      const int msgSize=bufSize+50; char msg[msgSize];
      sprintf_s(msg, msgSize, "\a\nCan't open for read %s", buf); perror(msg);
    }
  }
  return rfp;
} // openInput

bool errRead() { fprintf(wfpE,"\nError reading DLS.\n"); return F; }
//   -------

bool readDLS(FILE *rfp) {
//   ----------
  int smallest=LONG_MAX, biggest=LONG_MIN; int c=fgetc(rfp), d; if (c==EOF) return F;

  while ((c==' ')|(c=='\t')|(c=='\r')|(c=='\n')) c=fgetc(rfp);
  if (c=='[') { while ((c!=']')&(c!=EOF)) c=fgetc(rfp); if (c==EOF) return F; c=fgetc(rfp); }

  for (int i=0; i<NN; i++) {
    bool ceqd=F;
    while ((c==' ')|(c=='\t')|(c=='\r')|(c=='\n')) c=fgetc(rfp);
    if (c==EOF) { if (i>0) errRead(); return F; }
    if ((c>='a')&(c<='f')) c+=10-'a'; else if ((c>='A')&(c<='F')) c+=10-'A';
    else {
      if ((c<'0')|(c>'9')) return errRead(); c-='0'; d=fgetc(rfp); 
      while (('0'<=d)&(d<='9')) { c=c*10+d-'0'; d=fgetc(rfp); } ceqd=T;
    }
    if (c<smallest) smallest=c; if (c>biggest) biggest=c; iDLS[i]=c;
    if (ceqd) c=d; else c=fgetc(rfp);
  }
  smallRead=smallest; bigRead=biggest; zeroBase=(smallRead==0);
  ++dlsNumber; return T;
} // readDLS

bool isLS(CFint *Q) { // a Latin square?
//   ----
  for (int i=0, b; i<NN; i+=N) { b=0;
	  	for (int j=0, t; j<N; j++) {
      if (b&(t=(1<<Q[i+j])))
        { fprintf(wfpE,"Error: Duplicate %d in row %d.\n", Q[i+j], i/N+1); return F; } b|=t;
	  	}
	 }
 	for (int i=0, b; i<N; i++) { b=0;
	  	for (int j=0, t; j<NN; j+=N) {
	  	 	if (b&(t=(1<<Q[i+j])))
          { fprintf(wfpE,"Error: Duplicate %d in column %d.\n", Q[i+j], i+1); return F; } b|=t;
	  	}
 	}
 	return T;
} // isLS

bool isDLS(CFint *Q) { // a diagonal Latin square?
//   -----
  if (!(((smallRead==0)&(bigRead==M))|((smallRead==1)&(bigRead==N)))) {
    fprintf(wfpE,"Error: Numbers must be in the range 0..%d or 1..%d.\n", M, N); return F;
  }
  if (!zeroBase) for(int i=0; i<NN; i++ ) --Q[i]; zeroBase=T; if (!isLS(Q)) return F;
	 for (int j=0, b=0, t; j<NN; j+=O) {
	  	if (b&(t=(1<<Q[j]))) {
        fprintf(wfpE,"Error: Duplicate %d in \\diagonal.\n", Q[j]); return F; }	b|=t;
	 }
	 for (int j=M, b=0, t; j<=(NN-N); j+=M) {
	  	if (b&(t=(1<<Q[j]))) {
        fprintf(wfpE,"Error: Duplicate %d in /diagonal.\n", Q[j]); return F; }	b|=t;
	 }
	 return T;
} // isDLS
//========================================== output ==============================================

char outFile[outSize];
int openOutput(char *which) {
//  ----------
  int wfd=-1; int sub=0; const int baseSize=bufSize+25; char baseName[baseSize], buf[outSize]; 
  sprintf_s(baseName, baseSize, "%s", which); sprintf_s(buf, outSize, "%s.txt", baseName);
  do {
    if (_access_s(buf, 00)==ENOENT) break; sprintf_s(buf, outSize, "%s_%d.txt", baseName, ++sub);
  } while (T);
  if (_sopen_s(&wfd, buf, _O_CREAT|_O_WRONLY, _SH_DENYNO, _S_IWRITE)==0) {
    fprintf(wfpE,".. writing DLS to file %s\n", buf); sprintf_s(outFile, outSize, "%s", buf);
  } else {
    char msg[outSize+50]; sprintf_s(msg, outSize+50, "\a\nCan't open for write %s", buf); perror(msg);
  }
  return wfd;
} // openOutput

bool printCF4_10(const int i, const int wfd) {
//   -----------
  char *s=outPointer; CFint *p=CFs[i];
  for (int r=0; r<N; ++r) {
    int t=*p; *s++='0'+t;  p++;
    for (int c=1; c<N; ++c) { t=*p; *s++=' '; *s++='0'+t; p++; } *s++='\n';
  }
  *s++='\n'; outPointer=s; ++outDLSs;
  if (outDLSs>=nDLSs) {
    if (_write(wfd, outBuffer, writeSize)!=writeSize) return F;
    outPointer=outBuffer; outDLSs=0;
  }
  return T;
} // printCF4_10

bool printCF11_15(const int i, const int wfd) {
//   ------------
  char *s=outPointer; CFint *p=CFs[i];
  for (int r=0; r<N; ++r) {
    int t=*p; if (t<10) { *s++=' '; *s++='0'+t; } else { *s++='1'; *s++='0'+t-10; } p++;
    for (int c=1; c<N; ++c) {
      t=*p; *s++=' '; if (t<10) { *s++=' '; *s++='0'+t; }else { *s++='1'; *s++='0'+t-10; } p++;
    }
    *s++='\n';
  }
  *s++='\n'; outPointer=s; ++outDLSs;
  if (outDLSs>=nDLSs) {
    if (_write(wfd, outBuffer, writeSize)!=writeSize) return F;
    outPointer=outBuffer; outDLSs=0;
  }
  return T;
} // printCF11_15

bool printCF(const int i, const int wfd) { return (N<11) ? printCF4_10(i,wfd) : printCF11_15(i,wfd); }
//   -------

bool outLast(const int wfd) {
//   -------
  if (outDLSs>0) {
    int outBytes=outDLSs*dlsBytes;
    if (_write(wfd, outBuffer, outBytes)!=outBytes) return F;
  }
  return T;
} // outLast

int cmpCFs(const void *pp, const void *pq) {
//  ------
  CFint *p=*(CFint **)pp, *q=*(CFint **)pq; int i=0;
  while(i++<NN) if (*p<*q) return -1; else if (*p++>*q++) return 1; return 0;
} // cmpCFs

void sortCFs() { qsort(CFs, CFsIndex, sizeof(void *), cmpCFs); }
//   -------

int sortAndwriteCFs(const int wfd) {
//  ---------------
  int dup=0; sortCFs(); outPointer=outBuffer; outDLSs=0; printCF(0, wfd);
  for (int i=1; i<CFsIndex; i++)
    if (cmpCFs(&CFs[i-1], &CFs[i])==0) ++dup; else printCF(i, wfd);
  if (!outLast(wfd)) return 0; return CFsIndex-dup;     
} // sortAndwriteCFs

//void printDLS(CFint *a) {
////   --------
//  for (int r=0; r<N; ++r) { for (int c=0; c<N; ++c) printf(" %d", a[r*N+c]); putchar('\n'); }
//} // printDLS
//====================================== initialize ============================================

bool putRowOrder(CFint *ro) {
//   -----------
  //if ((ROindex)>=nXforms) { printf("\aROindex O\\F %d\n", ROindex); return F; }
  CFint *p=rowOrders[ROindex++]; for (int i=0; i<G; ++i) p[i]=ro[i]; return T;
} // putRowOrder

bool getRowOrders() {
//   ------------
  int count=0; CFint *p=rowOrderTmp; ROindex=0;
  int Mid=Odd ? N/2 : -1;  // -1 indicates no Mid for even dls.
  int incr=nXforms/G, iMax;
  if (incr==0) { iMax=nXforms-1; incr=1; } else iMax=G;
  for (int i=0; i<iMax; ++i) {
    //Set up the first combo that starts with i.
    bool *Free=valueFree; for (int v=0; v<N; ++v) Free[v]=T;
    if (Mid>=0) Free[Mid]=F; p[0]=i; Free[i]=F; Free[M-i]=F;
    for (int y=1; y<G; ++y) { int x=0; while (!Free[x]) ++x; p[y]=x; Free[x]=F; Free[M-x]=F; }
    int stop=i==(iMax-1) ? nXforms-1 : (i+1)*incr, j=G-1;
    while(T) {
      if (!putRowOrder(p)) return F; if (count++>=stop) break;
      Free[p[j]]=T; Free[M-p[j]]=T; int end=j<1;
      while (j>=1) {
        int next=F; Free[p[j]]=T; Free[M-p[j]]=T;
        while (T) {
          int x=++p[j]; if (x>=N) { if (j==1) end=T; --j; next=F; break; }        
          if (Free[x]) {
            next=T; Free[x]=F; Free[M-x]=F;
            for (int l=j+1; l<G; ++l) {
              for (int m=0; m<N; ++m) {
                if (!Free[m]) continue; Free[m]=F; Free[M-m]=F;
                if (p[l]<N) { Free[p[l]]=T; Free[M-p[l]]=T; } p[l]=m; break;
              }
            }
            j=G-1; break;
          } // if (Free[x]
        } // while (T
        if (next) break;
      } // while (j>=1
      if (end) break;
    } // while(T)
  } // for (int i=0;
  //printf("ROindex %d ", ROindex);
  return T;
} // getRowOrders

bool pushQueue(const CFint i, const CFint j) {
//   ---------
  //if ((Qindex)>=QueueSize) { printf("\aQueue O\\F %d\n", Qindex); return F;  }
  t_xform *p=&Queue[Qindex++]; p->w=i; p->x=j; return T;
} // pushQueue

int getPos(int v, bool *isOpp) {
//  ------
  for (int i=0; i<G; ++i) {
    int x=rowOrderTmp[i]; if ((x==v)|(x==(M-v))) { *isOpp=(x!=v); return i; }
  }
  return 0;
} // getPos

bool queueSwaps() {
//   ----------
  CFint *p=rowOrderTmp; Qindex=0, Xindex=0; for (int i=0; i<G; ++i) p[i]=i;
  for (int i=1; i<ROindex; ++i) { // 0 is the original DLS already written
    CFint *q=rowOrders[i]; int swaps=0;
    for (CFint j=0; j<G; ++j) {
      if (q[j]!=p[j]) {
        if (q[j]==(M-p[j])) {
          if (!pushQueue(j, maxCFint)) return F; ++swaps; p[j]=M-p[j];
        } else {
          bool isOpp; CFint k=getPos(q[j], &isOpp);
          if (isOpp) { pushQueue(k, maxCFint); ++swaps; p[k]=M-p[k]; }
          if (!pushQueue(j, k)) return F; ++swaps; CFint x=p[j]; p[j]=p[k]; p[k]=x;
        }
      }
    }
    swapCounts[Xindex++]=swaps;
  }
  //printf("Xindex %d Qindex %d \n", Xindex, Qindex);
  return T;
} // queueSwaps

void swap1RandC(const int i) {
//   ----------
  CFint w=Queue[i].w, n=N, y=M-w, **a=aDLS, *t=a[w]; a[w]=a[y]; a[y]=t;
  for (int r=0; r<n; ++r) { CFint t=a[r][w]; a[r][w]=a[r][y]; a[r][y]=t; }
} // swap1RandC

void swap2RandC(const int i) {
//   ----------
  t_xform *p=&Queue[i]; CFint w=p->w, x=p->x, n=N, y=M-w, z=M-x;
  CFint **a=aDLS, *t=a[w]; a[w]=a[x]; a[x]=t; t=a[y]; a[y]=a[z]; a[z]=t;
  for (int r=0; r<n; ++r) {
    CFint t=a[r][w]; a[r][w]=a[r][x]; a[r][x]=t; t=a[r][y]; a[r][y]=a[r][z]; a[r][z]=t;
  }
} // swap2RandC

bool pushX(const int i) {
//   -----
  //if (i>((nXforms-1)*NN)) { printf("\a nXforms O\\F %d\n", i); return F; }
  CFint *t=X+i; for (int r=0, j=0; r<N; ++r) for (int c=0; c<N; ++c) *t++=aDLS[r][c];
  return T;
} // pushX

/*
 * All transforms are made of aDLS, the square of natural order numbers 0 .. NN-1,
 * and entered into the string of tranforms in X.
 * Transforms of iDLS, the input DLS, are made as required by indexing into it with
 * numbers from X.
 */
bool getXforms() {
//   ---------
  int j=0, iX=0;
  for (CFint r=0, i=0; r<N; ++r) for (CFint c=0; c<N; ++c) aDLS[r][c]=i++;
  if (!pushX(iX)) return F; iX+=NN;
  for (int i=0; i<Xindex; ++i) {
    int n=swapCounts[i];
    for (int i=0; i<n; ++i) if (Queue[j].x==maxCFint) swap1RandC(j++); else swap2RandC(j++);
    if (!pushX(iX)) return F; iX+=NN;
  }
  return T;
} // getXforms

/*
 *  The 7 other aspects of the natural order, 0 .. N-1, square..
 */
void fillRotate() {
//   ----------
  CFint *r=rot[0];

  for (int i=0; i<N; i++) for (int j=0; j<N; j++) r[j*N+M-i]    =i*N+j; r=rot[1];
  for (int i=0; i<N; i++) for (int j=0; j<N; j++) r[(M-i)*N+M-j]=i*N+j; r=rot[2];
  for (int i=0; i<N; i++) for (int j=0; j<N; j++) r[(M-j)*N+i]  =i*N+j; r=rot[3];
  for (int i=0; i<N; i++) for (int j=0; j<N; j++) r[i*N+M-j]    =i*N+j; r=rot[4];
  for (int i=0; i<N; i++) for (int j=0; j<N; j++) r[(M-j)*N+M-i]=i*N+j; r=rot[5];
  for (int i=0; i<N; i++) for (int j=0; j<N; j++) r[(M-i)*N+j]  =i*N+j; r=rot[6];
  for (int i=0; i<N; i++) for (int j=0; j<N; j++) r[j*N+i]      =i*N+j;
} // fillRotate

bool initialize() {
//   ----------
  const int n=(N&1) ? N-1 : N, nm2d2=(n-2)/2; int i=2, f=1;
  
  NN=N*N; G=N/2; M=N-1; O=N+1; Odd=N&1;

  for (int j=1; j<nm2d2; ++j) i*=2; for (int j=1; j<=G; ++j) f*=j; nXforms=i*f;
  QueueSize=(int)(double)nXforms*1.65;
  n1_nfr_1_N=(n/2-1)*n; n1_nfr_2_N=(n==4)?4:(n-4)*n1_nfr_1_N; n1_nbd_M_O=n1_nfr_1_N;
  n2_nfr_1_N=(nXforms)/n1_nfr_1_N+2; n2_nfr_2_N=(nXforms)/n1_nfr_2_N+2;
  n2_nbd_M_O=(nXforms)/n1_nbd_M_O+2;

  dlsBytes=NN+NN+1; nDLSs=100000; if (N>10) { dlsBytes+=NN; nDLSs=500000/N; }
  writeSize=nDLSs*dlsBytes;

  CFsize=NN*CFintBytes; CFsIndex=0; numBlocksAllocated=0;

  if (allocateStore()) {
    if (!getRowOrders()) return F; if (!queueSwaps()) return F;
    if (!getXforms()) return F; fillRotate(); return T;
  }
  return F;
} // initialize
//================================== get canonical form ========================================

void getSmallNFR(CFint *a, CFint *t) { // natural order first row zeroBase
//   -----------
  CFint n[maxN], *b=sDLS+N; bool smaller=F; for (int c=0; c<N; ++c) n[a[t[c]]]=c;
  for (int i=N; i<NN; ++i) { int x=n[a[t[i]]], y=*b++; if (x!=y) { smaller=(x<y); break; }}
  if (smaller) { b=sDLS+N; for (int i=N; i<NN; ++i) *b++=n[a[t[i]]]; }
} // getSmallNFR

void getXformedDLSs1() {
//   ---------------
  bool hit=F;

  // Try to find a Xform aspect that has a[N]==a[1], which will become value 1.
  for (int i=0; i<n1_nfr_1_N; ++i) { // CF[1][0]==1
    const int J=nXforms/n1_nfr_1_N; int *p=nfr_1_N[i], x=*p++, y=*p++, *q=p; CFint *a;
    a=A[0]; if (a[x]==a[y]) { hit=T;
      for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); getSmallNFR(A[7],t); }}
    a=A[1]; if (a[x]==a[y]) { hit=T; p=q;
      for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); getSmallNFR(A[6],t); }}
    a=A[2]; if (a[x]==a[y]) { hit=T; p=q;
      for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); getSmallNFR(A[5],t); }}
    a=A[3]; if (a[x]==a[y]) { hit=T; p=q;
      for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); getSmallNFR(A[4],t); }}
  }

  // Try to find a Xform aspect that has a[N]==a[2], which will become value 2.
  if (!hit) {
    for (int i=0; i<n1_nfr_2_N; ++i) { // CF[1][0]==2
      const int J=nXforms/n1_nfr_2_N; int *p=nfr_2_N[i], x=*p++, y=*p++, *q=p; CFint *a;
      a=A[0]; if (a[x]==a[y]) { hit=T; p=q; for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); }}
      a=A[1]; if (a[x]==a[y]) { hit=T; p=q; for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); }}
      a=A[2]; if (a[x]==a[y]) { hit=T; p=q; for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); }}
      a=A[3]; if (a[x]==a[y]) { hit=T; p=q; for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); }}
      a=A[4]; if (a[x]==a[y]) { hit=T; p=q; for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); }}
      a=A[5]; if (a[x]==a[y]) { hit=T; p=q; for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); }}
      a=A[6]; if (a[x]==a[y]) { hit=T; p=q; for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); }}
      a=A[7]; if (a[x]==a[y]) { hit=T; p=q; for (int j=0; j<J; ++j) { CFint *t=X+*p++; getSmallNFR(a,t); }}
    }
  }

  // Find any Xform aspect that makes the smallest result.
  if (!hit) { CFint *t=X;
    for (int i=0; i<nXforms; ++i) {
      getSmallNFR(A[0],t); getSmallNFR(A[1],t); getSmallNFR(A[2],t); getSmallNFR(A[3],t);
      getSmallNFR(A[4],t); getSmallNFR(A[5],t); getSmallNFR(A[6],t); getSmallNFR(A[7],t); t+=NN;
    }
  }
} // getXformedDLSs1

void getSmallNBD(CFint *a, CFint *t) { // natural order first diagonal zeroBase
//   -----------
  CFint n[maxN], *b=sDLS; bool smaller=F, equal=T; for (int c=0, i=0; c<NN; c+=O) n[a[t[c]]]=i++;
  for (int c=M; c<(NN-M); c+=M) { int x=n[a[t[c]]], y=b[c]; if (x!=y) { smaller=(x<y); equal=F; break; }}
  if (equal) for (int i=0; i<NN; ++i) { int x=n[a[t[i]]], y=*b++; if (x!=y) { smaller=(x<y); break; }}
  if (smaller) { b=sDLS; for (int i=0; i<NN; ++i) *b++=n[a[t[i]]];
  }
} // getSmallNBD

void getXformedDLSs2() {
//   ---------------
  // Try to find a Xform aspect that has a[M]==a[O], which will become value 1.
  bool hit=F;
  for (int i=0; i<n1_nbd_M_O; ++i) { // CF[0][M]==1
    const int J=nXforms/n1_nbd_M_O; int *p=nbd_M_O[i], x=*p++, y=*p++, *q=p; CFint *a;
    a=A[0]; if (a[x]==a[y]) { hit=T; for (int j=0; j<J; ++j) getSmallNBD(a,X+*p++); } p=q;
    a=A[1]; if (a[x]==a[y]) { hit=T; for (int j=0; j<J; ++j) getSmallNBD(a,X+*p++); } p=q;
    a=A[2]; if (a[x]==a[y]) { hit=T; for (int j=0; j<J; ++j) getSmallNBD(a,X+*p++); } p=q;
    a=A[3]; if (a[x]==a[y]) { hit=T; for (int j=0; j<J; ++j) getSmallNBD(a,X+*p++); } p=q;
    a=A[4]; if (a[x]==a[y]) { hit=T; for (int j=0; j<J; ++j) getSmallNBD(a,X+*p++); } p=q;
    a=A[5]; if (a[x]==a[y]) { hit=T; for (int j=0; j<J; ++j) getSmallNBD(a,X+*p++); } p=q;
    a=A[6]; if (a[x]==a[y]) { hit=T; for (int j=0; j<J; ++j) getSmallNBD(a,X+*p++); } p=q;
    a=A[7]; if (a[x]==a[y]) { hit=T; for (int j=0; j<J; ++j) getSmallNBD(a,X+*p++); }
  }

  // Find any Xform aspect that makes the smallest result, (with the smallest /diagonal)..
  if (!hit) { CFint *t=X;
    for (int i=0; i<nXforms; ++i) {
      getSmallNBD(A[0],t); getSmallNBD(A[1],t); getSmallNBD(A[2],t); getSmallNBD(A[3],t);
      getSmallNBD(A[4],t); getSmallNBD(A[5],t); getSmallNBD(A[6],t); getSmallNBD(A[7],t); t+=NN;
    }
  }
} // getXformedDLSs2

bool pushCF(CFint *CF) {
//   ------
  if ((CFsIndex>=CFsSize)&&!allocateCFs(numBlocksAllocated*2)) return F;
  CFint *p=CFs[CFsIndex++]; for (int i=0; i<NN; ++i) *p++=*CF++; return T;
} // pushCF

/*
 *  The 8 aspects of iDLS, the input DLS.
 */
void getAspects() {
//   ----------
  CFint *a=A[0], *r=rot[0];

  for (int i=0; i<NN; ++i) a[i]=iDLS[i];    a=A[1];
  for (int i=0; i<NN; ++i) a[i]=iDLS[r[i]]; a=A[2]; r=rot[1];
  for (int i=0; i<NN; ++i) a[i]=iDLS[r[i]]; a=A[3]; r=rot[2];
  for (int i=0; i<NN; ++i) a[i]=iDLS[r[i]]; a=A[4]; r=rot[3];
  for (int i=0; i<NN; ++i) a[i]=iDLS[r[i]]; a=A[5]; r=rot[4];
  for (int i=0; i<NN; ++i) a[i]=iDLS[r[i]]; a=A[6]; r=rot[5];
  for (int i=0; i<NN; ++i) a[i]=iDLS[r[i]]; a=A[7]; r=rot[6];
  for (int i=0; i<NN; ++i) a[i]=iDLS[r[i]];
} // getAspects

bool getCF(const int format) { 
//   -----
  bool ok=T;
  if (isDLS(iDLS)) {
    getAspects();
    if (format==1) { sDLS[N]=N; getXformedDLSs1(); }
    else { sDLS[0]=N; sDLS[M]=N; getXformedDLSs2(); }
    ok=pushCF(sDLS);
  } else printf("Error: Square %d is not an acceptable DLS.\n\n", dlsNumber);
  return ok;
} // getCF

bool getTmatches1() {
//   ------------
  if (allocateIntArray(&nfr_1_N, n1_nfr_1_N, n2_nfr_1_N)
        && allocateIntArray(&nfr_2_N, n1_nfr_2_N, n2_nfr_2_N)) {
    CFint *t=X, a=maxCFint, b=maxCFint; int *p=NULL, j=-1;
    for (int i=0; i<nXforms; ++i) { // CF[1][0]==CF[0][1]
      CFint x=t[1], y=t[N];
      if ((x==a)&(y==b)) *p++=i*NN; else { p=nfr_1_N[++j]; *p++=x; *p++=y; *p++=i*NN; }
      a=x; b=y; t+=NN;
    }

    t=X; a=maxCFint; b=maxCFint; j=-1;
    for (int i=0; i<nXforms; ++i) { // CF[1][0]==CF[0][2]
      CFint x=t[2], y=t[N];
      if ((x==a)&(y==b)) *p++=i*NN; else { p=nfr_2_N[++j]; *p++=x; *p++=y; *p++=i*NN; }
      a=x; b=y; t+=NN;
    }
    return T;
  }
  return F;
} // getTmatches1

bool getTmatches2() {
//   ------------
  if (allocateIntArray(&nbd_M_O, n1_nbd_M_O, n2_nbd_M_O)) {
    CFint *t=X, a=maxCFint, b=maxCFint; int *p=NULL, j=-1;

    for (int i=0; i<nXforms; ++i) { // CF[0][M]==CF[1][1]
      CFint x=t[M], y=t[O];
      if ((x==a)&(y==b)) *p++=i*NN; else { p=nbd_M_O[++j]; *p++=x; *p++=y; *p++=i*NN; }
      a=x; b=y; t+=NN;
    }
    return T;
  }
  return F;
} // getTmatches2

bool setupFormat(int *format) {
//   -----------
  if (*format<1) *format=1; else if (*format>2) *format=2;

  freeQueue(); freeCFintArray(&rowOrders, nXforms); freeCFints(&rowOrderTmp);
  freeCFints(&swapCounts); freeCFintArray(&aDLS, N);
  
  if (*format==1) { for (int i=0; i<N; ++i) sDLS[i]=i; return getTmatches1(); }
  else return getTmatches2();
} // setupFormat
//============================================ main ===========================================

bool validOrder() {
//   ----------
  if ((N<minN)|(N>maxN)) {
    printf("\aERROR: Supported orders are %d to %d.\n", minN, maxN); return F; }
  return T;
} // validOrder

int main() {
//  ----
  //printf("sizeof(Uint) %d sizeof(void *) %d\n", sizeof(Uint), sizeof(void *)); //outLocalTime();
  bool another=T, inputSize=T, ok=F; wfpE=stdout;

  do {
    if (inputSize) N=getInt("\nOrder");
    if (validOrder()&&initialize()) {
      int format=getInt("Format, (1: first row or 2: \\diagonal)");
      if (setupFormat(&format)) {
        FILE *rfp=openInput();
        if (rfp!=NULL) {
          char buf[bufSize]; sprintf_s(buf, bufSize, "output%dCF%d", N, format);
          int wfd=openOutput(buf);
          if (wfd!=-1) {
            int count=0; ok=T; CFsIndex=0; dlsNumber=0;
            time_t startTime=time(NULL);
            while (readDLS(rfp)) { if (!getCF(format)) { ok=F; break; }}
            if (ok) {
              count=sortAndwriteCFs(wfd);
              printf("number of DLS %d CFs %d\n", CFsIndex, count);
            }
            printElapsedTime(startTime, "\n");
            _close(wfd); if (count==0) remove(outFile); 
          }
          fclose(rfp);
        } // rfp
      } // setupFormat
      freeStore();
      printf("\nContinue? y (yes) or n (no) or the order: ");
      if (getYorOrder(&N)) inputSize=(N==0); else another=F;
    }
  } while (another);
  
  printf("\nPress a key to close the console");
  while (!_kbhit()) Sleep(250);
  return ok ? EXIT_SUCCESS : EXIT_FAILURE;
} // main