/*
 *  File:    Franklin.cpp
 *  Author:  S Harry White
 *  Created: 2012-01-01
 *  Updated: 2020-05-03
 *    Added howMany.
 *  Updated 2022-01-02
 *    Tidy code.
 *  Updated: 2023-01-07
 *    Change to compile with g++ for macOS.
 *  Updated: 2023-03-28
 *    Fix store allocation. Call freeStore() for each allocation.
 */

/*
 *  Makes Franklin normal magic squares of order 8k.
 */

#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 int Uint;
const int startSize=64; int **xSquare=NULL;
const bool F=false, T=true;
//================================================= store =================================================

void freeStore(const int size) {
//   ---------
  if (xSquare!=NULL) { for (int i=0; i<size; i++) free(xSquare[i]); free(xSquare); xSquare=NULL; }
} // freeStore

bool allocateStore(const int size) {
//   -------------
  bool ok; xSquare=(int**) malloc(size*sizeof(int*)); ok=(xSquare!=NULL);
  if (ok) {
    int numAllocated=size;
    for(int i=0; i<size; i++) {
      int *p=(int*) malloc(size*sizeof(int)); xSquare[i]=p;
      if (p==NULL) { numAllocated=i; ok=F; break; }
    }
    if (!ok) freeStore(numAllocated);
  }
  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');
} // getY

bool getYorOrder(int *n) { // if user inputs size instead of 'y' or 'n', use it
//   -----------
  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

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

const int bufSize=1024;
bool openDir() {
//   -------
  char buf[bufSize], msg[bufSize+50]; int sub=0;
  const char *baseName="./FranklinNormalMagicSquares"; strcpy(buf, baseName);
  do {
    if (mkdir(buf, 0775)==0) break;
    if (errno!=EEXIST) {
      snprintf(msg, bufSize+50, "Can't make folder %s", buf); perror(msg); return F;
    }
    snprintf(buf, bufSize, "%s_%d", baseName, ++sub);
  } while (T); 
  if (chdir(buf)!=0) {
    snprintf(msg, bufSize+50, "Can't open folder %s", buf); perror(msg); return F;
  }
  printf("Output files are in folder %s\n", buf); return T;
} // openDir

FILE *open_File(const int n) {
//    ---------
  FILE *wfp=NULL; char buf[bufSize]; snprintf(buf, bufSize, "./Order%d.txt", n);
  if ((wfp=fopen(buf, "a"))==NULL) {
    char msg[bufSize+50]; snprintf(msg, bufSize+50, "\a\nCan't open for write %s", buf); perror(msg);
  }
  return wfp;
} // open_File

int fieldWidth(const int n) {
//  ----------
  if (n==1) return 1; int i=n*n, width=1; while ((i=i/10)!=0) ++width; return width;
} // fieldWidth

typedef bool (*t_fprintFW)(FILE *fp, const int i);

bool fprintFW1(FILE *fp, const int i) { return fprintf(fp, "%1d",  i)>0; }
bool fprintFW2(FILE *fp, const int i) { return fprintf(fp, "%2d",  i)>0; }
bool fprintFW3(FILE *fp, const int i) { return fprintf(fp, "%3d",  i)>0; }
bool fprintFW4(FILE *fp, const int i) { return fprintf(fp, "%4d",  i)>0; }
bool fprintFW5(FILE *fp, const int i) { return fprintf(fp, "%5d",  i)>0; }
bool fprintFW6(FILE *fp, const int i) { return fprintf(fp, "%6d",  i)>0; }
bool fprintFW7(FILE *fp, const int i) { return fprintf(fp, "%7d",  i)>0; }
bool fprintFW8(FILE *fp, const int i) { return fprintf(fp, "%8d",  i)>0; }
bool fprintFW9(FILE *fp, const int i) { return fprintf(fp, "%9d",  i)>0; }
bool fprintFWa(FILE *fp, const int i) { return fprintf(fp, "%10d", i)>0; }

static t_fprintFW fprintFW[]={ NULL,
  fprintFW1, fprintFW2, fprintFW3, fprintFW4, fprintFW5,
  fprintFW6, fprintFW7, fprintFW8, fprintFW9, fprintFWa
};
const int maxFieldWidth=10;

bool printSquare(FILE *wfp, const int n) {
//   -----------
  const int fw0=fieldWidth(n), fw=fw0+1; if (fw>maxFieldWidth) return F;
  for (int i=0; i<n; i++) {
    if (!fprintFW[fw0](wfp, xSquare[i][0])) return F;
    for (int j=1; j<n; j++) if (!fprintFW[fw](wfp, xSquare[i][j])) return F;
    if (fputc('\n', wfp)==EOF) return F;
  }
  return fputc('\n', wfp)!=EOF;
} // printSquare
//================================================= check =================================================

bool isHalfMagic(int **x, const int n, Uint halfSum) {
//   -----------
  Uint sumXY=0, sumYX=0; const int nd2=n/2;
  for (int i=0; i<n; ++i) {
    Uint rowH=0, colH=0; for (int j=0; j<nd2; ++j) { rowH+=x[i][j]; colH+=x[j][i]; }
    if ((rowH!=halfSum)|(colH!=halfSum)) return F; sumXY+=x[i][n-i-1]; sumYX+=x[i][i];
  }
  Uint magicSum=halfSum+halfSum; return (sumXY==magicSum)&(sumYX==magicSum);
} // isHalfMagic

bool isBentDiagonal(int **x, const int n, Uint magicSum) {
//   --------------
  const int nm1=n-1, nd2=n/2;

  for (int i=n-1; i>=0; --i) { // bent left ?
    Uint sum=0;
    for (int r=0, c=i; r<nd2; ++r, --c) { const int cModn=c<0 ? c+n : c; sum+=x[r][cModn]+x[nm1-r][cModn]; }
    if (sum!=magicSum) return F;
  }
  for (int i=0; i<n; ++i) {  // bent right ?
    int sum=0;
    for (int r=0, c=i; r<nd2; ++r, ++c) { const int cModn=c<n ? c : c-n; sum+=x[r][cModn]+x[nm1-r][cModn]; }
    if (sum!=magicSum) return F;
  }
  for (int i=n-1; i>=0; --i) { // bent up ?
    int sum=0;
    for (int c=0, r=i; c<nd2; ++c, --r) { const int rModn=r<0 ? r+n : r; sum+=x[rModn][c]+x[rModn][nm1-c]; }
    if (sum!=magicSum) return F;
  }
  for (int i=0; i<n; ++i) { // bent down ?
    int sum=0;
    for (int c=0, r=i; c<nd2; ++c, ++r) { const int rModn=r<n ? r : r-n; sum+=x[rModn][c]+x[rModn][nm1-c]; }
    if (sum!=magicSum) return F;
  }
  return T;
} // isBentDiagonal

bool isCompact(int **x, const int n, Uint S4) {
//   ---------
  const int nm1=n-1;

  for (int r=0; r<nm1; ++r) {
    for (int c=0; c<nm1; ++c) { if ((x[r][c]+x[r][c+1]+x[r+1][c]+x[r+1][c+1])!=S4) return F; }
    if ((x[r][nm1]+x[r][0]+x[r+1][nm1]+x[r+1][0])!=S4) return F;
    if ((x[nm1][r]+x[0][r]+x[nm1][r+1]+x[0][r+1])!=S4) return F;
  }
  return T;
} // isCompact

const char *notMagicFranklin =
  "\a\nERROR: Square is NOT magic Franklin. Please report by email\n";

bool checkMagicFranklin(const int n) {
//   ------------------
  // Calculation as magicSum/2 might give F result on overflow.
  const Uint halfSum=n/4*(n*n+1), S2=n*n+1;

  if (isHalfMagic(xSquare, n, halfSum)
        &&isBentDiagonal(xSquare, n, halfSum+halfSum) // Overflow is ok.
          &&isCompact(xSquare, n, S2+S2)) return T;
  else { printf("%s",notMagicFranklin); return F; }
} // checkMagicFranklin
//=============================================== make ===============================================

void makeActual(const int n) {
//   ----------
  const int pPlus=n*n/2, mPlus=pPlus+1;
  for (int i=0; i<n; i++) for (int j=0; j<n; j++) xSquare[i][j]+=xSquare[i][j]>0 ? pPlus : mPlus;
} // makeActual

void seed_rand() { srand((Uint)time(NULL)); }
//   ---------

int random_(const int x) { return rand()%x; }
//  -------

void complement(const int n) {
//   ----------
  if (random_(2)) { for (int r=0; r<n; ++r) for (int c=0; c<n; ++c) xSquare[r][c]=-xSquare[r][c]; }
} // complement

void swapRowsCols(const int n) {
//   ------------
  int **x=xSquare; const int nd2=n/2;
  // Randomly opt to swap top and bottom halves.
  if (random_(2)) { for (int r=0; r<nd2; ++r) { int *tmp=x[r]; x[r]=x[r+nd2]; x[r+nd2]=tmp; } }
  // Randomly choose first 2 rows/columns from those that start in the first quarter.
  int k=random_(nd2), a=0, b=k%2==1 ? k-1 : k;
  if (a!=b) {
    int loop=2; 
    while (loop--!=0) {
      int *tmp=x[a]; x[a]=x[b]; x[b]=tmp;
      for (int r=0; r<n; ++r) { const int tmp=x[r][a]; x[r][a]=x[r][b]; x[r][b]=tmp; } ++a; ++b;
    }
  }
  // Randomly swap pairs of rows/columns that start in the first quarter.
  k=random_(nd2); a=k%2==1 ? k-1 : k; k=random_(nd2); b=k%2==1 ? k-1 : k;
  if (a!=b) {
    int loop=2;
    while (loop--!=0) {
      int *tmp=x[a]; x[a]=x[b]; x[b]=tmp;
      for (int r=0; r<n; ++r) { const int tmp=x[r][a]; x[r][a]=x[r][b]; x[r][b]=tmp; } ++a; ++b;
    }
  }
} // swapRowsCols

/*
 * Each index represents 8 bones values, as follows:
 *
 *    0:  +-( 0.5, 1.5, 2.5, 3.5 )
 *    1:  +-( 4.5, 5.5, 6.5, 7.5  )
 *       ...
 *
 *  Insert each of indexes 0 to n*n/8-1 into 8 appropriately selected cells.
 *  Convert indexes to bones values: 4 times index plus factor 0.5 to 3.5.
 *  Change sign of alternating pairs of cells to minus.
 */
void makeSquare8k(const int n, FILE *wfp) {
//   ------------
  int **x=xSquare; const int nm1=n-1, nd2=n/2, nd2m2=nd2-2; int u=0, v=n*n/8-1;  

  // insert indexes 0 .. n*n/16-1
  for (int r=0; r<nd2; r+=3) {
    int  a, b, d, q=nm1-r;
    for (int c=0; c<nd2; ++c) {
      a=nd2m2-c; b=nm1-a; d=nm1-c; x[r][c]=u; x[q][c]=u; x[r][d]=u; x[q][d]=u;
      x[r+2][a]=u; x[q-2][a]=u; x[r+2][b]=u; x[q-2][b]=u++;
      c+=3; a=nd2-c; b=nm1-a; d=nm1-c; x[r][c]=u; x[q][c]=u; x[r][d]=u; x[q][d]=u;
      x[r+2][a]=u; x[q-2][a]=u; x[r+2][b]=u; x[q-2][b]=u++;
    }
    ++r; q=nm1-r;
    for (int c=nd2-2; c>=0;  c-=3) {
      a=nd2m2-c; b=nm1-a; d=nm1-c; x[r][c]=u; x[q][c]=u; x[r][d]=u; x[q][d]=u;
      x[r+2][a]=u; x[q-2][a]=u; x[r+2][b]=u; x[q-2][b]=u++;
      c--; a=nd2-c; b=nm1-a; d=nm1-c; x[r][c]=u; x[q][c]=u; x[r][d]=u; x[q][d]=u;
      x[r+2][a]=u; x[q-2][a]=u; x[r+2][b]=u; x[q-2][b]=u++;
    }
  }
  // insert indexes n*n/8-1 .. n*n/16
  for (int r=0; r<nd2; r+=3) {
    int a, b, d, q=nm1-r;
    for (int c=1; c<nd2; c+=3) {
      a=nd2-c; b=nm1-a; d=nm1-c; x[r][c]=v; x[q][c]=v; x[r][d]=v; x[q][d]=v;
      x[r+2][a]=v; x[q-2][a]=v; x[r+2][b]=v; x[q-2][b]=v--;
      c++; a=nd2m2-c; b=nm1-a; d=nm1-c; x[r][c]=v; x[q][c]=v; x[r][d]=v; x[q][d]=v;
      x[r+2][a]=v; x[q-2][a]=v; x[r+2][b]=v; x[q-2][b]=v--;
    }
    ++r; q=nm1-r;
    for (int c=nd2-1; c>=0; --c) {
      a=nd2-c; b=nm1-a; d=nm1-c; x[r][c]=v; x[q][c]=v; x[r][d]=v; x[q][d]=v;
      x[r+2][a]=v; x[q-2][a]=v; x[r+2][b]=v; x[q-2][b]=v--;
      c-=3; a=nd2m2-c; b=nm1-a; d=nm1-c; x[r][c]=v; x[q][c]=v; x[r][d]=v; x[q][d]=v;
      x[r+2][a]=v; x[q-2][a]=v; x[r+2][b]=v; x[q-2][b]=v--;
    }
  }

  const int a[2][4]={ {  1,    4,   3,   2 }, {   3,   2,   1,   4 } };
                  //{ { 0.5, 3.5, 2.5, 1.5 }, { 2.5, 1.5, 0.5, 3.5 } };
  const int k[2]={ 1, 0 };  int i=1, j=0; bool incr2=F;
  // convert indexes to values
  for (int r=0; r<nd2; ++r) {
    const int q=nm1-r; i=k[i];
    for (int c=0; c<n; ++c) {
      int e=a[i][j], f=a[k[i]][j]; x[r][c]=4*x[r][c]+e; x[r+2][c]=4*x[r+2][c]+f;
      x[q][c]=4*x[q][c]+f; x[q-2][c]=4*x[q-2][c]+e; j=j==3 ? 0 : j+1;
    }
    if (incr2) r+=2; incr2=!incr2;
  }

  // sign
  for (int r=0; r<n; r+=2) for (int c=0; c<n; c+=4) {
    x[r][c]*=-1; x[r][c+1]*=-1; x[r+1][c+2]*=-1; x[r+1][c+3]*=-1;
  }
  swapRowsCols(n); complement(n); makeActual(n);
}  // makeSquare8k

bool makeSquare(const int n, FILE *wfp) {
//   ----------
  makeSquare8k(n, wfp); checkMagicFranklin(n); return !printSquare(wfp, n);
} // makeSquare
//================================================== main ============================================

bool validOrder(const int n) {
//   ----------
  if ((n<=0)|(n%8!=0)) {
    printf("\aERROR: Order %d is not a positive integral multiple of 8.\n", n); return F;  
  }
  return T;
} // validOrder

int main() {
//  ----
  bool another=T, inputSize=T, writeError=F, ok=F; seed_rand();
  if (openDir()) {
    int n=0; 
    do {     
      if (inputSize) { printf("\nInput the order of the square: "); n=getSize(); }
      if (validOrder(n)) {
      	if (allocateStore(n)) {
	        FILE *wfp=NULL;
          if ( (wfp=open_File(n))!=NULL) {
            printf("\nHow many? "); int howMany=getSize();
            while (howMany--) { writeError=makeSquare(n, wfp); if (writeError) break; }
            if (!writeError) ok=T; fclose(wfp);
	        }
          freeStore(n); 
        } else printf("\a\nERROR: Storage allocation failed.\n");
      }
      if (writeError) { perror("\a\nError writing file"); another=F; }
      else {
        printf("\nMake another square? input y (yes) or n (no) or the order of the square: ");
        if (getYorOrder(&n)) inputSize=(n==0); else another=F; 
      }
    } while (another);
  } // if (openDir()
  printf("\nHit return to close the console.");
  while (T) if (getchar()=='\n') break; return ok ? EXIT_SUCCESS : EXIT_FAILURE;
} // main