/*
 *  File:    GetOrthogonal.cpp
 *  Author:  S Harry White
 *  Created: 2020-08-02
 *  Updated: 2020-09-27
 *    Fix open_input when file name includes .txt.
 *  Updated: 2020-10-02
 *    Fix reading numbers > 20.
 *  Updated: 2021-05-12
 *    Fix isLS for orders > 31.
 *  Updated: 2021-05-16
 *    Improve storage allocation.
 *    Include extra code for x64 build, commented out for Win32.
 *  Updated: 2021-05-17
 *    For the first square with the maximum number of orthogonal mates, (if >= 10), copy
 *    the square and its orthogonal mates to a separate file.
 *  Updated: 2021-05-19
 *    Provide a choice of two formats:
 *      1) Get the counts of orthogonal pairs for each square and write the first square
 *         with the maximum number of pairs and its mates to a file.
 *      2) Prompt for input of a square number to write with its mates to a file.
 *  Updated: 2021-06-07
 *    Change opNumber to unsigned to accomodate a bigger number before overflow.
 *  Updated: 2023-01-28
 *    Change to compile with g++ for macOS.
 */

/*
 *  Finds orthogonal pairs in a file of Latin squares.
 *  There is a choice of two formats:
 *    1) Get the counts of orthogonal pairs for each square and a list of mates. 
 *       For the first square with the maximum number of pairs, write it and its mates to a file.
 *    2) Prompt for input of a square number to write with its mates to a file.
 */

//#include <ctype.h>
#include <curses.h>
#include <errno.h>
#include <fcntl.h>
#include <limits.h>
#include <math.h>
#include <stdbool.h>
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <sys/stat.h>
//#include <sys/types.h>
#include <time.h>
#include <unistd.h>

typedef unsigned char LSint;
typedef unsigned int   Uint;
typedef unsigned long Uint64;

void **blocks=NULL; LSint **xLS=NULL;
const int maxN=UCHAR_MAX, startLS=50000;
const int startNumBlocks=1,
          pointerSize=sizeof(void *),
          bufSize=1024, outSize=bufSize+50;
const Uint maxSquares=UINT_MAX/pointerSize;

Uint64 opNumber; 
int N, M, O, NN, maxLS, sqNumber, lsNumber, smallRead, bigRead, lsWithMaxOps, 
    allocatedBlocks, maxBlocksIncr, LSperBlock, LSsize, allocatedLS;
typedef struct { int ls; Uint num; } t_maxOps; t_maxOps maxOps; 

const bool F=false, T=true; bool **QR=NULL, zeroBase, reportAlloc;

void initGlobals() {
//   -----------
  M=N-1; O=N+1; NN=N*N; maxLS=0; sqNumber=0; lsNumber=0; maxOps.ls=0; maxOps.num=0;
  lsWithMaxOps=0; allocatedBlocks=0; allocatedLS=0; smallRead=INT_MAX; bigRead=INT_MIN;
  LSsize=NN; reportAlloc=F; zeroBase=T;

  const int OneHundredMillion=100000000, FiftyThousand=50000;
  const Uint oneB=1000000000;
  LSperBlock=OneHundredMillion/LSsize;
  if (LSperBlock>FiftyThousand) LSperBlock=FiftyThousand;
  else if (LSperBlock<1) LSperBlock=1;
  maxBlocksIncr=oneB/(LSperBlock*LSsize);
} // initGlobals
//========================================= store ===============================================

const char *storeAllocFail="Storage allocation failed"; bool bell=T;
bool reportError(const char *msg) {
//   -----------
  printf("%sError: %s.\n", bell?"\a\n":"",  msg);
  if (bell) bell=F; return F;
} // reportError

void freeLS(int allocated) {
//   ------
  if (blocks!=NULL) {
    if (allocated!=0) {
      /* printf(".. freeing LS\n"); */ for (int i=0; i<allocated; i++) free(blocks[i]);
      if (xLS!=NULL) { free(xLS); xLS=NULL; allocatedLS=0; }
    }
    free(blocks); blocks=NULL;
  }
  allocatedBlocks=0;
} // freeLS

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

void freeStore() { freeLS(allocatedBlocks); freeBoolArray(&QR, N+1); }
//   ---------

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

bool allocateLS(int numNew) {
//   ----------
  bool ok=F; const int numOld=allocatedBlocks;
  union ovly { Uint64 i; LSint *p; };
  if ((numNew-numOld)>maxBlocksIncr) numNew=numOld+maxBlocksIncr;
  const Uint numNewLS=numNew*LSperBlock; // Uint for malloc(numNewLS*pointerSize)

  if (numNewLS<=maxSquares) {
    if (reportAlloc) {
      Uint nsq=numNewLS, msq=nsq/1000000, tsq=nsq%1000000/1000, rsq=nsq%1000;
      printf("    .. increasing LS store to ");
      if (msq==0)
        if (tsq==0) printf("%11d\n", rsq);
        else printf("%7d,%03d\n", tsq, rsq);
      else printf("%3d,%03d,%03d\n", msq, tsq, rsq);
    } else reportAlloc=T;

    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(LSperBlock*LSsize); 
        if (p==NULL) { ok=F; allocated=i; break; } tmp[i]=p;
      }
	    freeLS(0); blocks=tmp;
	    if (ok) {
        LSint **tmp=(LSint **)malloc(numNewLS*pointerSize);
        if ((ok=(tmp!=NULL))) {
          for (int i=0; i<lsNumber; i++) tmp[i]=xLS[i];
          free(xLS); xLS=tmp; int j=lsNumber;
          for (int i=numOld; i<numNew; i++) {
            int k=0; LSint *p=(LSint *)(blocks[i]);
            while (k++<LSperBlock) {
              xLS[j++]=p; ovly ip; ip.p=p; ip.i+=LSsize; p=ip.p;           
            }
          }
	        allocatedBlocks=numNew; allocatedLS=numNewLS; maxLS=numNewLS-1;
        }
      } else freeLS(allocated);
    }
  }
  if (!ok) reportError(storeAllocFail); return ok;
} // allocateLS

bool allocateStore() {
//   -------------
  bool ok=T;
  if ((ok=allocateLS(startNumBlocks))) ok=allocateBoolArray(&QR, N+1, N+1);

  if (!ok) { reportError(storeAllocFail); freeStore(); } return ok;
} // allocateStore
//========================================= input ===============================================

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

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=-1;

  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

int getInt() { int n=0; scanf("%d", &n); clearLine(getchar()); return n; }
//  -------  

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; }
  *s='\0'; return T;
} // getFileName

void adjustName(char *buf) {
//   ----------
  char *s=buf; bool txt=F;
  if ((*s!='.')&(*s!='/')) { // Assume in current folder.
    char tmp[bufSize]; snprintf(tmp, bufSize, "./%s", buf); strcpy(buf, tmp);
  }
  if (access(buf, R_OK)!=0) {
    while (*s++!='\0')
      ; --s;  if ((strlen(buf)>6)&&(*s--=='t')&&(*s--=='x')&&(*s--=='t')&&(*s=='.')) txt=T;
    if (!txt) strcat(buf, ".txt"); // no .txt entered, add it
  }
} // adjustName

FILE *openInput(char buf[bufSize]) {
//    ---------
  char *rFname=NULL; FILE *rfp=NULL;
  do {
    printf("\nEnter the name of the squares file: ");
    if (getFileName(buf, bufSize-6)) { rFname=buf; break; }
    printf("\a\nCan't read the file name. Try again? y (yes) or n (no) "); if (!getY()) break;
  } while (T);
  if (rFname!=NULL) {
    adjustName(buf);
    if ((rfp=fopen(buf, "r"))==NULL) {
      const int msgSize=bufSize+50; char msg[msgSize];
      snprintf(msg, msgSize, "\a\nCan't open for read %s", buf); perror(msg);
    }
  }
  return rfp;
} // openInput

bool isLS(LSint *x) { // a Latin square?
//   ----
  if (!(((smallRead==0)&(bigRead==M))|((smallRead==1)&(bigRead==N)))) {
    printf("Error square %d: Numbers must be in the range 0..%d or 1..%d.\n", sqNumber, M, N); return F;
  }
  bool u[maxN+1];
  for (int i=0; i<NN; i+=N) {
    for (int a=0; a<=N; ++a) u[a]=F;
   	for (int j=0; j<N; j++) {
      if (u[x[i+j]])
        { printf("Error square %d: Duplicate %d in row %d.\n", sqNumber, x[i+j], i/N+1); return F; } u[x[i+j]]=T;
	  	}
	 }
 	for (int i=0; i<N; i++) {
    for (int a=0; a<=N; ++a) u[a]=F;
	  for (int j=0; j<NN; j+=N) {
	  	if (u[x[i+j]])
        { printf("Error square %d: Duplicate %d in column %d.\n", sqNumber, x[i+j], i+1); return F; } u[x[i+j]]=T;
	  }
 	}
 	return T;
} // isLS

bool errRead() { printf("\nError reading square %d.\n", ++sqNumber); return T; }
//   -------

bool readLS(FILE *rfp) {
//   ------
  if ((lsNumber>maxLS)&&!allocateLS(allocatedBlocks+allocatedBlocks)) return F;
  int smallest=INT_MAX, biggest=INT_MIN; LSint *x=xLS[lsNumber];
  int c=fgetc(rfp), d; if (c==EOF) return F;

  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 T; }
    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; x[i]=c;
    if (ceqd) c=d; else c=fgetc(rfp);
  }
  smallRead=smallest; bigRead=biggest; ++sqNumber;
  if (isLS(x)) { if (bigRead==N) zeroBase=F; ++lsNumber; } return T;
} // readDLS
//=================================================== output ================================================

void stripName(char *inFname, char *obuf) {
//   ---------
  char *s=inFname; while (*s++!='\0');
  while (--s!=inFname)
    if (*s=='.') *s='\0';             // Remove .txt
    else if (*s=='/') { ++s; break; } // Remove any directory path
  strcpy(obuf, s); 
} // stripName

char outFile1[outSize], outFile2[outSize], outFile3[outSize];
FILE *openOutput(char inFname[bufSize], const int which, const int target) {
//    ----------
  int sub=0; FILE *wfp=NULL; const int baseSize=bufSize+25;
  char baseName[baseSize], buf[outSize]; stripName(inFname, buf);
  if (which==1) snprintf(baseName, baseSize, "./%s-orthCounts", buf);
  else if (which==2) snprintf(baseName, baseSize, "./%s-orthNos", buf);
  else snprintf(baseName, baseSize, "./%s-%dorths", buf, target);
  snprintf(buf, outSize, "%s.txt", baseName);
  do {
    if (access(buf, F_OK)!=0) break;
    snprintf(buf, outSize, "%s_%d.txt", baseName, ++sub);
  } while (T);
  if ((wfp=fopen(buf, "w"))==NULL) {
    char msg[outSize+50]; snprintf(msg, outSize+50, "\a\nCan't open for write %s", buf); perror(msg);
  } else {
    printf("..output file %s\n", buf);
    if (which==1) snprintf(outFile1, outSize, "%s", buf);
    else if (which==2) snprintf(outFile2, outSize, "%s", buf);
    else snprintf(outFile3, outSize, "%s", buf);
  }
  return wfp;
} // openOutput

bool printLS3(LSint *x, FILE *wfp) {
//   --------
  for (int r=0, i=0; r<N; r++) {
    if (fprintf(wfp, "3%d", x[i++])<0) return F;
    for (int c=1; c<N; c++) { if (fprintf(wfp, " %3d", x[i++])<0) return F; }
    if (fputc('\n', wfp)==EOF) return F;
  }
  if (fputc('\n', wfp)==EOF) return F;
  return T;
} // printLS3

bool printLS2(LSint *x, FILE *wfp) {
//   --------
  for (int r=0, i=0; r<N; r++) {
    if (fprintf(wfp, "%2d", x[i++])<0) return F;
    for (int c=1; c<N; c++) { if (fprintf(wfp, " %2d", x[i++])<0) return F; }
    if (fputc('\n', wfp)==EOF) return F;
  }
  if (fputc('\n', wfp)==EOF) return F;
  return T;
} // printLS2

bool printLS1(LSint *x, FILE *wfp) {
//   --------
  for (int r=0, i=0; r<N; r++) {
    if (fprintf(wfp, "%d", x[i++])<0) return F;
    for (int c=1; c<N; c++) { if (fprintf(wfp, " %d", x[i++])<0) return F; }
    if (fputc('\n', wfp)==EOF) return F;
  }
  if (fputc('\n', wfp)==EOF) return F;
  return T;
} // printLS1

bool printLS(LSint *x, FILE *wfp) {
//   -------
  const int max1=zeroBase ? 10 : 9, max2=zeroBase ? 100 : 99;
  return (N<=max1) ? printLS1(x, wfp) : (N<=max2) ? printLS2(x, wfp) : printLS3(x, wfp);
} // printLS
//============================================ get pairs =============================================

bool isOLS(LSint *x, LSint *y, bool **b) {
//   -----
  for (int k=0; k<O; ++k) for (int l=0; l<O; ++l) b[k][l]=F;
  for (int k=0; k<NN; ++k) { const int u=x[k], v=y[k]; if (b[u][v]) return F; b[u][v]=T; } 
  return T;
} // isOLS

bool getOrts(LSint **X, const int i, FILE *wfpc, FILE *wfpn) {
//   -------
  bool first=T; LSint *x=X[i]; int n=0, m=0; Uint ops=0; const int nper=10;
  for (int j=0; j<lsNumber; ++j) if (j!=i) {
    LSint *y=X[j];
    if (isOLS(x, y, QR)) {
      if (first) { if (fprintf(wfpn, "%d: ", i+1)<0) return F; first=F; }
      if (fprintf(wfpn, " %d", j+1)<0) return F; ++n;
      if (++m==nper) { if (fputc('\n', wfpn)==EOF) return F; m=0; }
      ++ops;
    }
  }
  if (n>0) if (fputc('\n', wfpn)==EOF) return F;
  if (fprintf(wfpc, "%12d %12u%s\n", i+1, ops, (i+1)%10==0?"\n":"")<0) return F;
  if (ops>=maxOps.num) {
    if (ops==maxOps.num) ++lsWithMaxOps;
    else { lsWithMaxOps=1; maxOps.ls=i; maxOps.num=ops; }
  }
  opNumber+=ops; return T;
} // getOrts

bool getOrtsSageMath(LSint **X, const int i, FILE *wfpc, FILE *wfpn) {
//   ---------------
  bool first=T; LSint *x=X[i]; int n=0, m=0; Uint ops=0; const int nper=11;
  for (int j=0; j<lsNumber; ++j) if (j!=i) {
    LSint *y=X[j];
    if (isOLS(x, y, QR)) { ++m;
      if (first) { if (fprintf(wfpn, "%s%d: [%d", opNumber==0?"":"],\n",i+1,j+1)<0) return F; first=F; }
      else if (m==nper) { if (fprintf(wfpn, "],\n%d: [%d",i+1,j+1)<0) return F; m=1; }
      else if (fprintf(wfpn, ",%d", j+1)<0) return F; ++n;
      ++ops;
    }
  }
  if (i==(lsNumber-1)) if (fprintf(wfpn, "]\n")<0) return F;

  if (fprintf(wfpc, "%12d %12u%s\n", i+1, ops, (i+1)%10==0?"\n":"")<0) return F;
  if (ops>=maxOps.num) {
    if (ops==maxOps.num) ++lsWithMaxOps;
    else { lsWithMaxOps=1; maxOps.ls=i; maxOps.num=ops; }
  }
  opNumber+=ops; return T;
} // getOrtsSageMath

bool getTargetOrts(LSint **X, const int i, Uint *count, FILE *wfpt) {
//   -------------
  Uint num=0; LSint *x=X[i]; if (!printLS(x,wfpt)) return F;
  for (int j=0; j<lsNumber; ++j) if (j!=i) {
    LSint *y=X[j]; if (isOLS(x, y, QR)) { ++num; if (!printLS(y,wfpt)) return F; }
  } *count=num; return T;
} // getTargetOrts
//============================================== main ===============================================

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

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

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

bool startCountsFile(FILE *wfpc) {
//   ---------------
  if (fprintf(wfpc, "%27s", "orthogonal\n")<0) return F;
  if (fprintf(wfpc, "%12s %12s", "square", "pairs\n")<0) return F;
  if (fprintf(wfpc, "%12s %15s", "------", "----------\n\n")<0) return F;
  return T;
} // startCountsFile

int main() {
//  ----
  //Uint64 y=UINT_MAX; printf("y %lu y/2 %lu sizeof(Uint64) %lu\n", y, y/2u, sizeof(Uint64) );
  bool inputOrder=T, ok=T; char buf[bufSize]; //outputLocalTime();
  do {
    if (inputOrder) { printf("\nOrder? "); N=getInt(); } initGlobals();
    if (validOrder()&&allocateStore()) {
      FILE *rfp=openInput(buf); time_t startTime=time(NULL);
      if (rfp!=NULL) {
        int target=0;
        printf("Choose 1 - get counts and maximum pairs, or 2 - get pairs for one square: ");
        int mode=getInt(); if (mode<1) mode=1; else if (mode>2) mode=2;
        FILE *wfpc=NULL, *wfpn=NULL; if (mode==1) { wfpc=openOutput(buf, 1, 0); wfpn=openOutput(buf, 2, 0); }
        if ((mode==2)|((wfpc!=NULL)&(wfpn!=NULL))) { lsNumber=0;
          while (readLS(rfp))
            ; printf("squares %d ", lsNumber);
          if (mode==1) {
            startCountsFile(wfpc);
            for (int i=0; i<lsNumber; ++i) { if (!getOrts(xLS, i, wfpc, wfpn)) { ok=F; break; }}
            //for (int i=0; i<lsNumber; ++i) { if (!getOrtsSageMath(xLS, i, wfpc, wfpn)) { ok=F; break; }}
            fclose(wfpc); fclose(wfpn);
            printf("total orthogonal pairs %lu\n", opNumber/2u);
            if (ok) {
              if ((lsNumber==0)|(opNumber==0)) { remove(outFile1); remove(outFile2); }
              else {
                if (opNumber>0) { target=maxOps.ls+1;
                  printf("Maximum pairs for square %d: %u\n", target, maxOps.num);
                  if (lsWithMaxOps>1)
                    printf("There %s %d other square%s with this maximum number of pairs.\n",
                            lsWithMaxOps==2 ? "is" : "are", lsWithMaxOps-1, lsWithMaxOps==2 ? "" : "s");
                    else printf("This is the only square with this maximum number of pairs.\n");
                }
              }
            } else { remove(outFile1); remove(outFile2); perror("\n\aError writing file"); }
          } else {
            printf("\nGet pairs for square number, (1 .. %d)? ", lsNumber); target=getInt();
            if ((target<1)|(target>lsNumber)) { printf("\aBad number %d\n", target); ok=F; }
          }
          if (ok&(target>0)) {
            FILE *wfpp=openOutput(buf, 3, target);
            if (wfpp!=NULL) {
              Uint count=0;
              ok=getTargetOrts(xLS, target-1, &count, wfpp); fclose(wfpp);
              printf("Pairs for square %d: %u\n", target, count);
              if (!ok) perror("\n\aError writing file"); if (!ok|(count==0)) remove(outFile3);
            }
          }
        }
        fclose(rfp);
      } // rfp!=NULL
      printElapsedTime(startTime); freeStore();
    } // allocateStore
    printf("\nAnother order? input y (yes), n (no) or the order: ");
    if (getYorOrder(&N)) inputOrder=(N<0); else break;
  } while (T);
  printf("\nHit return to close the console.");
  while (T) if (getchar()=='\n') break; return ok ? EXIT_SUCCESS : EXIT_FAILURE;
} // main