/*
 *  File:    DoubleBorderedRectangles.cpp
 *  Author:  S Harry White
 *  Created: 2023-06-03
 *  Updated: 2023-06=19
 *    Changes to compile with g++.
 */

/*
 *  Adapted from a program written by Johan Öfverstedt that uses
 *  "Constraint-Based Local Search" and "Tabu Search" techniques:
 *
 *        "Water Retention on Magic Squares Solver",
 *        http://sourceforge.net/projects/wrmssolver/.
 *
 *  Makes magic rectangles with double borders.
 */

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

int Rend, Cend, RCend, // full rectangle
    R, C,      // of sub-rectangles
    N,         // sqrt(R*C)
    Rm1, Cm1,  // R-1, C-1
    Rm2, Cm2,  // R-2, C-2
    Rd2, Cd2,  // R/2, C/2
    Rp1,       // R+1
    RC,        // R*C
    RCm1,      // RC-1
    RCd2,      // RC/2
    S2,        // RC+1
    S1,        // S2/2 for odd R, C
    SR,        // R*(RC+1)/2
    SC,        // C*(RC+1)/2

    *xRectangle, *specific, *freeNum, *freeIx, *opp,
    *row, *col, *rowSum, *colSum, diagLsum, diagRsum,
    *tabu, tabuLen, Nfree, RCfree, RCfree_1, RCfree_2;
const bool F=false, T=true;
bool *diagL, *diagR, *numUsed, allZero, odd, fullSearch, RCswapped;

enum borderRectangleType { none, cyclic, opposed };
borderRectangleType borderType;
enum magicType { normalSemiMagic, normalMagic, Associative, Bordered };
const magicType firstType=normalSemiMagic, lastType=Bordered;
magicType rectangleType;

#define max(A,B) (A)>(B)?(A):(B)
#define min(A,B) (A)<(B)?(A):(B)

void seed_rand() { srand((unsigned int)time(NULL)); }
//   ---------
int random_(int x) { return rand()%x; }
//  -------

const int bufSize=1024; bool bell=T; char fmt[bufSize];
bool reportError(const char *msg) {
//   -----------
  printf("%sError: %s.\n", bell ? "\a\n" : "",  msg); return bell=F;
} // reportError
//============================================ store ==============================================

void freeStore() {
//   ---------
  if (xRectangle!=NULL) { free(xRectangle); xRectangle=NULL; }
  if (specific  !=NULL) { free(specific);   specific  =NULL; }
  if (freeNum   !=NULL) { free(freeNum);    freeNum   =NULL; }
  if (freeIx    !=NULL) { free(freeIx);     freeIx    =NULL; }
  if (opp       !=NULL) { free(opp);        opp       =NULL; }
  if (row       !=NULL) { free(row);        row       =NULL; }
  if (col       !=NULL) { free(col);        col       =NULL; }
  if (rowSum    !=NULL) { free(rowSum);     rowSum    =NULL; }
  if (colSum    !=NULL) { free(colSum);     colSum    =NULL; }
  if (tabu      !=NULL) { free(tabu);       tabu      =NULL; }
  if (numUsed   !=NULL) { free(numUsed);    numUsed   =NULL; }
  if (diagL     !=NULL) { free(diagL);      diagL     =NULL; }
  if (diagR     !=NULL) { free(diagR);      diagR     =NULL; }
}

bool allocateInts(int **x, const int size) { return (*x=(int*) malloc(size*sizeof(int)))!=NULL; }
//   ------------

bool allocateBools(bool **x, const int size) { return (*x=(bool*) malloc(size*sizeof(bool)))!=NULL; }
//   -------------

bool allocateStore() {
//   -------------
  bool ok=T; const int R=Rend, C=Cend, RC=R*C;
  if ((ok=allocateInts(&xRectangle, RC)))
    if ((ok=allocateInts(&specific, RC)))
      if ((ok=allocateInts(&freeNum, RC)))
        if ((ok=allocateInts(&freeIx, RC)))
          if ((ok=allocateInts(&opp, RC))) 
            if ((ok=allocateInts(&row, RC)))
              if ((ok=allocateInts(&col, RC)))
                if ((ok=allocateInts(&rowSum, R)))
                  if ((ok=allocateInts(&colSum, C)))
                    if ((ok=allocateInts(&tabu, RC)))
                      if ((ok=allocateBools(&numUsed, RC+1)))
                        if ((ok=allocateBools(&diagL, RC)))
                          ok=allocateBools(&diagR, RC);
  if (!ok) { freeStore(); reportError("Storage allocation failed"); }
  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

int getInt(const char *prompt) {
//  -------
  printf("%s", prompt); int n=0; scanf("%d", &n); clearLine(getchar());
  if (n<=0) reportError("Invalid input"); return n;
} // getInt

bool getInts(int *p, int *q, int c) {
//   -------
  bool ok=F; *p=0; *q=0; if (c<0) do { c=getchar(); } while ((c==' ')|(c=='\t'));
  if ( ('1'<=c)&(c<='9') ) {
    int i=c-'0'; while ( ('0'<=(c=getchar()))&&(c<='9') ) i=i*10+c-'0'; *p=i;
    if ((c==' ')|(c=='\t')) {
      do { c=getchar(); } while ((c==' ')|(c=='\t'));
      if ( ('1'<=c)&(c<='9') ) {
        int i=c-'0'; while ( ('0'<=(c=getchar()))&&(c<='9') ) i=i*10+c-'0'; *q=i;
      }
    }
    if (*q==0) *q=*p; ok=T;
  }  
  clearLine(c); if (!ok) return reportError("Invalid input"); return T;
} // getInts

bool getYorInts(int *p, int *q) {
//   ----------
  bool ok=F; int c; *p=-1; do { c=getchar(); } while ((c==' ')|(c=='\t')|(c=='\n') );
  if ( (c=='Y')|(c=='y') ) ok=T; else if ( (c!='N')&(c!='n') ) return getInts(p, q, c);  
  clearLine(c); return ok;
} // getYorInts

char *getLine(char *s, const char nl) {
//    -------
  int c, i=0; char *sin=s; do { c=getchar(); } while ((c==' ')|(c=='\t')|(c==nl) );
  *s=c; if (c=='\n') return s; // nothing
  while (i++<bufSize) if ( (*++s=getchar())=='\n') break;
  if (*s!='\n') {
    printf("\nType list too long.\n"); clearLine(*s); return NULL;
  }
  while ( (*--s==' ')||(*s=='\t') ); // strip any trailing whitespace
  *++s='\0'; s=sin; const int length=strlen(s);
  for (int i=0; i<length; ++i) s[i]=tolower(s[i]); return sin;
} // getLine
//============================================ output ========================================

bool openDir() {
//   -------
  char buf[bufSize], msg[bufSize+50]; const char *baseName="./DoubleBorderedRectangles";
  strcpy(buf, baseName); int sub=0;
  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


int fieldWidth() {
//  ----------
  int i=Rend*Cend, 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; int FW0, FW;

bool filePrintCR(int *x, FILE *wfp) {
//   -----------
  for (int i=0; i<C; i++) {
    if (!fprintFW[FW0](wfp, x[i])) return F;
    for (int j=1; j<R; j++) if (!fprintFW[FW](wfp, x[j*C+i])) return F;
    if (fputc('\n', wfp)==EOF) return F;
  }
  return fputc('\n', wfp)!=EOF;
} // filePrintCR

bool filePrint(int *x, FILE *wfp) {
//   ---------
  if (RCswapped) return filePrintCR(x,wfp);
  for (int i=0; i<R; i++) {
    if (!fprintFW[FW0](wfp, x[i*C])) return F;
    for (int j=1; j<C; j++) if (!fprintFW[FW](wfp, x[i*C+j])) return F;
    if (fputc('\n', wfp)==EOF) return F;
  }
  return fputc('\n', wfp)!=EOF;
} // filePrint

bool consolePrint(int *x) { return filePrint(x, stdout); }
//   ------------

char OutputFileName[bufSize];
FILE *openOutput() {
//    ----------
  char buf[bufSize], defaultName[bufSize]; int sub=0; FILE *wfp=NULL;
  if (RCswapped) snprintf(defaultName, bufSize, "./%dx%d", Cend, Rend);
  else           snprintf(defaultName, bufSize, "./%dx%d", Rend, Cend);
  snprintf(buf, bufSize, "%s.txt", defaultName);
  
  do {
    if (access(buf, F_OK)!=0) break;
    snprintf(buf, bufSize, "%s_%d.txt", defaultName, ++sub);
  } while (T);
  if ((wfp=fopen(buf, "w"))==NULL) {
    char msg[bufSize+50]; snprintf(msg, bufSize+50, "\a\nCan't open for write %s", buf); perror(msg);
  } else {
    printf("\nOutput file is %s\n", buf); strcpy(OutputFileName, buf);
  }
  return wfp;
} // openOutput
//=========================================== make ============================================

bool getSpecific() {
//   -----------
  bool ok=T; for (int i=0; i<RC; ++i) specific[i]=0;
  if (!allZero) {
    int r=R-4, c=C-4, rc=r*c, s2=rc+1, add=(S2-s2)/2, incr=2*C+2;
    for (int a=0; a<r; ++a, incr+=4) {
      for (int b=0; b<c; ++b) {
        int i=a*c+b; specific[i+incr]=xRectangle[i]+add;
      }
    }
  }
  return ok;
} // getSpecific

bool getFree() {
//   -------
  int next=0;
  if (allZero) {
    if (odd&(rectangleType==Associative)) {
      for (int i=0; i<RCd2; ++i) freeIx[next++]=i; for (int i=RCd2+1; i<RC; ++i) freeIx[next++]=i;
    } else { for (int i=0; i<RC; ++i) freeIx[next++]=i; }
    RCfree=next; next=0; for (int i=1; i<=RC; ++i) if (!numUsed[i]) freeNum[next++]=i;
  } else {
    for (int i=0; i<RC; ++i) if (specific[i]==0) freeIx[next++]=i; else numUsed[specific[i]]=T;
    RCfree=next; next=0; for (int i=1; i<=RC; ++i) if (!numUsed[i]) freeNum[next++]=i;
  }
  Nfree=(int)sqrt((double)RCfree+1.); RCfree_1=RCfree-1;
  RCfree_2=RCfree-2; tabuLen=((N==4)&allZero) ? 5 : (Nfree+Nfree)/3; return T;
} // getFree

bool initializeCommon() {
//   ----------------
  Rm1=R-1; Rm2=R-2; Rp1=R+1; Cm1=C-1; Cm2=C-2; RC=R*C; Rd2=R/2; Cd2=C/2; RCd2=RC/2; RCm1=RC-1;
  S2=RC+1; S1=S2/2; SR=C*S2/2; SC=R*S2/2;
  N=(R==C) ? R : (int)sqrt((double)RC);

  for (int i=0; i<RC; ++i) {
    const int idC=i/C, imC=i%C; row[i]=idC; col[i]=imC;
    if ((rectangleType==normalMagic)|((R==C)&(rectangleType==Bordered)&(borderType!=none)))
      { diagL[i]=(idC==imC); diagR[i]=((idC+imC)==Cm1); }
  }
  for (int i=0; i<=RC; ++i) numUsed[i]=F;
  return getSpecific();
} // initializeCommon

bool initializeStart() {
//   ---------------
  RCswapped=(odd&&(Cend>Rend)) | ((borderType==opposed)&&(Rend<Cend));
  if (RCswapped) { int t=Cend; Cend=Rend; Rend=t; }
  if (!allocateStore()) return F; RCend=Rend*Cend;
  FW0=fieldWidth(); FW=FW0+1; if (FW>maxFieldWidth) return reportError("Output format overflow");
  return T;
} // initializeStart

bool initializeSemi() { return initializeCommon()&&getFree(); }
//   --------------
  
bool initializeMagic() { return initializeCommon()&&getFree(); }
//   ---------------

void initializeOpp() {
//   -------------
  if (rectangleType==Associative) {
    for (int i=0; i<RC; ++i) opp[i]=RCm1-i;
  } else { // Bordered
    for (int i=0; i<RC; ++i) {
      const int r=row[i], c=col[i]; int t=0;
      if (borderType==none) {
        if ((c>1)&(c<Cm2)) { if (r==0) t=i+C*Rm1; else if (r==1) t=i+C*(R-3); }
        else if ((r>1)&(r<Rm2)) { if (c==0) t=i+Cm1; else if (c==1) t=i+C-3; }
        if (r==0) { if (c==0) t=RCm1; else if (c==1) t=i+Cm1; else if (c==Cm2) t=i+C+1; else if (c==Cm1) t=i+C*(R-2)+1; }
        else if (r==1) { if (c==1) t=i+C*(R-2)-3; else if (c==Cm2) t=i+C*(R-4)+3; }
        else if (r==Rm2) { if (c==0) t=i+C+1; else if (c==Cm1) t=i+C-1; }
      } else {
        if ((c>1)&(c<Cm2)) { if ((r==0)|(r==Rm2)) t=i+C; } else if ((r>1)&(r<Rm2)) { if ((c==0)|(c==Cm2)) t=i+1; }
        if (borderType==cyclic) {
          if (r==0) { if (c<=1) t=i+C; else if (c==Cm2) t=i+1; } else if (r==1) { if (c==Cm2) t=i+1; }
          else if (r==Rm2) { if (c==0) t=i+1; else if (c>=Cm2) t=i+C; } else if (r==Rm1) { if (c==0) t=i+1; }
        } else { // opposed
          if ((r==0)|(r==Rm2)) { if ((c==0)|(c==1)|(c==Cm2)|(c==Cm1)) t=i+C; }
        }
      }
      if (t!=0) { opp[i]=t; opp[t]=i; }
    }
  }
} // initializeOpp

bool initializeSym() {
//   -------------
  bool ok=initializeCommon(); initializeOpp();
  if (ok&allZero&odd) { specific[RCd2]=S1; numUsed[S1]=T; } 
  return ok ? getFree() : F;
} // initializeSym

int lastViolation;
int updateViolationSemi(const int i1, const int i2) {
//  -------------------
  const int r1=row[i1], c1=col[i1], r2=row[i2], c2=col[i2],
	          x1=xRectangle[i1], x2=xRectangle[i2], d=x1-x2; // swapped

  int oldv=0, newv=0;
  if (r1!=r2) {
    oldv+=abs(rowSum[r1]-SR)+abs(rowSum[r2]-SR); rowSum[r1]+=d; rowSum[r2]-=d;
    newv+=abs(rowSum[r1]-SR)+abs(rowSum[r2]-SR);
  }
  if (c1!=c2) {
    oldv+=abs(colSum[c1]-SC)+abs(colSum[c2]-SC); colSum[c1]+=d; colSum[c2]-=d;
    newv+=abs(colSum[c1]-SC)+abs(colSum[c2]-SC);
  }
 	return lastViolation+=(newv-oldv);
} // updateViolationSemi

int updateViolationMagic(const int i1, const int i2) { // R==C
//  --------------------
  const int r1=row[i1], c1=col[i1], r2=row[i2], c2=col[i2],
	          x1=xRectangle[i1], x2=xRectangle[i2], d=x1-x2; // swapped

  int oldv=0, newv=0;
  if (r1!=r2) {
    oldv+=abs(rowSum[r1]-SR)+abs(rowSum[r2]-SR); rowSum[r1]+=d; rowSum[r2]-=d;
    newv+=abs(rowSum[r1]-SR)+abs(rowSum[r2]-SR);
  }
  if (c1!=c2) {
    oldv+=abs(colSum[c1]-SR)+abs(colSum[c2]-SR); colSum[c1]+=d; colSum[c2]-=d;
    newv+=abs(colSum[c1]-SR)+abs(colSum[c2]-SR);
  }
  if (r1==c1) {
    if (r2!=c2) { // First is in diag 1 but second is not.
      oldv+=abs(diagLsum-SR); diagLsum+=d; newv+=abs(diagLsum-SR);
    }
  } else if (r2==c2) { // First is not in diag 1 but second is.
    oldv+=abs(diagLsum-SR); diagLsum-=d; newv+=abs(diagLsum-SR);
  }

  if (r1==(Rm1-c1)) {
    if (r2!=(Rm1-c2)) { // First is in diag 2 but second is not.
      oldv+=abs(diagRsum-SR); diagRsum+=d; newv+=abs(diagRsum-SR);
    }
  } else if (r2==(Rm1-c2)) { // First is not in diag 2 but second is.
    oldv+=abs(diagRsum-SR); diagRsum-=d; newv+=abs(diagRsum-SR);
	 }
 	return lastViolation+=(newv-oldv);
} // updateViolationMagic

int updateViolationSym(int i1, int i2) {
//  ------------------
  int r1=row[i1], c1=col[i1], r2=row[i2], c2=col[i2];
	const int x1=xRectangle[i1], x2=xRectangle[i2], d=x1-x2; // swapped

  int oldv=0, newv=0;
  if (r1!=r2) {
    oldv+=abs(rowSum[r1]-SR)+abs(rowSum[r2]-SR); rowSum[r1]+=d; rowSum[r2]-=d;
    newv+=abs(rowSum[r1]-SR)+abs(rowSum[r2]-SR);
  }
  if (c1!=c2) {
    oldv+=abs(colSum[c1]-SC)+abs(colSum[c2]-SC); colSum[c1]+=d; colSum[c2]-=d;
    newv+=abs(colSum[c1]-SC)+abs(colSum[c2]-SC);
  }

  i1=opp[i1];
  if (i1!=i2) { i2=opp[i2];
    r1=row[i1]; c1=col[i1]; r2=row[i2]; c2=col[i2];
    if (r1!=r2) {
      oldv+=abs(rowSum[r1]-SR)+abs(rowSum[r2]-SR); rowSum[r1]-=d; rowSum[r2]+=d;
      newv+=abs(rowSum[r1]-SR)+abs(rowSum[r2]-SR);
    }
    if (c1!=c2) {
      oldv+=abs(colSum[c1]-SC)+abs(colSum[c2]-SC); colSum[c1]-=d; colSum[c2]+=d;
      newv+=abs(colSum[c1]-SC)+abs(colSum[c2]-SC);
    }
  }
 	return lastViolation+=(newv-oldv);
} // updateViolationSym

int updateViolationBord2(int i1, int i2) {
//  --------------------
  int r1=row[i1], c1=col[i1], r2=row[i2], c2=col[i2];
	const int x1=xRectangle[i1], x2=xRectangle[i2], d=x1-x2; // swapped

  int oldv=0, newv=0;
  if (r1!=r2) {
    if ((r1<2)|(r1>=Rm2)) { oldv+=abs(rowSum[r1]-SR); rowSum[r1]+=d; newv+=abs(rowSum[r1]-SR); }
    if ((r2<2)|(r2>=Rm2)) { oldv+=abs(rowSum[r2]-SR); rowSum[r2]-=d; newv+=abs(rowSum[r2]-SR); }
  }
  if (c1!=c2) {
    if ((c1<2)|(c1>=Cm2)) { oldv+=abs(colSum[c1]-SC); colSum[c1]+=d; newv+=abs(colSum[c1]-SC); }
    if ((c2<2)|(c2>=Cm2)) { oldv+=abs(colSum[c2]-SC); colSum[c2]-=d; newv+=abs(colSum[c2]-SC); }
  }
  if ((R==C)&(borderType!=none)) {
    if (r1==c1) {
      if (r2!=c2) { // First is in diag 1 but second is not.
        oldv+=abs(diagLsum-SR); diagLsum+=d; newv+=abs(diagLsum-SR);
      }
    } else if (r2==c2) { // First is not in diag 1 but second is.
      oldv+=abs(diagLsum-SR); diagLsum-=d; newv+=abs(diagLsum-SR);
    }
    if (r1==(Rm1-c1)) {
      if (r2!=(Rm1-c2)) { // First is in diag 2 but second is not.
        oldv+=abs(diagRsum-SR); diagRsum+=d; newv+=abs(diagRsum-SR);
      }
    } else if (r2==(Rm1-c2)) { // First is not in diag 2 but second is.
      oldv+=abs(diagRsum-SR); diagRsum-=d; newv+=abs(diagRsum-SR);
	  }
  }

  i1=opp[i1];
  if (i1!=i2) { i2=opp[i2];
    r1=row[i1]; c1=col[i1]; r2=row[i2]; c2=col[i2];
    if (r1!=r2) {
      if ((r1<2)|(r1>=Rm2)) { oldv+=abs(rowSum[r1]-SR); rowSum[r1]-=d; newv+=abs(rowSum[r1]-SR); }
      if ((r2<2)|(r2>=Rm2)) { oldv+=abs(rowSum[r2]-SR); rowSum[r2]+=d; newv+=abs(rowSum[r2]-SR); }
    }
    if (c1!=c2) {
      if ((c1<2)|(c1>=Cm2)) { oldv+=abs(colSum[c1]-SC); colSum[c1]-=d; newv+=abs(colSum[c1]-SC); }
      if ((c2<2)|(c2>=Cm2)) { oldv+=abs(colSum[c2]-SC); colSum[c2]+=d; newv+=abs(colSum[c2]-SC); }
    }
    if ((R==C)&(borderType!=none)) {
      if (r1==c1) {
        if (r2!=c2) { // First is in diag 1 but second is not.
          oldv+=abs(diagLsum-SR); diagLsum-=d; newv+=abs(diagLsum-SR);
        }
      } else if (r2==c2) { // First is not in diag 1 but second is.
        oldv+=abs(diagLsum-SR); diagLsum+=d; newv+=abs(diagLsum-SR);
      }
      if (r1==(Rm1-c1)) {
        if (r2!=(Rm1-c2)) { // First is in diag 2 but second is not.
          oldv+=abs(diagRsum-SR); diagRsum-=d; newv+=abs(diagRsum-SR);
        }
      } else if (r2==(Rm1-c2)) { // First is not in diag 2 but second is.
        oldv+=abs(diagRsum-SR); diagRsum+=d; newv+=abs(diagRsum-SR);
	    }
    }
  }
  return lastViolation+=(newv-oldv);
} // updateViolationBord2

int violationSemi() {
//  -------------
	int violation=0, *pr=xRectangle;

  // Row, column constraints
	for (int i=0; i<R; ++i) {
	 	int rSum=0; for (int j=0; j<C; ++j) rSum+=*pr++;rowSum[i]=rSum; violation+=abs(rSum-SR);
  }
  for (int i=0; i<C; ++i) {
	 	 int cSum=0, *pc=xRectangle+i;
	 	 for (int j=0; j<R; ++j) { cSum+=*pc; pc+=C; } colSum[i]=cSum; violation+=abs(cSum-SC);
  }
 	return lastViolation=violation;
} // violationSemi

int violationBord2() {
//  --------------
	int violation=0, *pr=xRectangle;

  // Row, column constraints
	for (int i=2; i<Rm2; ++i) rowSum[i]=SR;
  for (int i=0; i<2; ++i) {
	 	int rSum=0; for (int j=0; j<C; ++j) rSum+=*pr++; rowSum[i]=rSum; violation+=abs(rSum-SR);
  }
  pr+=(R-4)*C;
	for (int i=Rm2; i<R; ++i) {
	 	int rSum=0; for (int j=0; j<C; ++j) rSum+=*pr++; rowSum[i]=rSum; violation+=abs(rSum-SR);
  }
  for (int i=2; i<Cm2; ++i) colSum[i]=SC;
  for (int i=0; i<2; ++i) {
    int cSum=0, *pc=xRectangle+i;
	  for (int j=0; j<R; ++j) { cSum+=*pc; pc+=C; } colSum[i]=cSum; violation+=abs(cSum-SC);
  }
  for (int i=Cm2; i<C; ++i) {
    int cSum=0, *pc=xRectangle+i;
	  for (int j=0; j<R; ++j) { cSum+=*pc; pc+=C; } colSum[i]=cSum; violation+=abs(cSum-SC);
  }
  if ((R==C)&(borderType!=none)) { // diagonal constraints
    int Lsum=0, Rsum=0, *pd=xRectangle;
	  for (int i=0; i<R; ++i) { Lsum+=pd[i]; Rsum+=pd[Rm1-i]; pd+=R; }
    diagLsum=Lsum; diagRsum=Rsum; violation+=abs(Lsum-SR)+abs(Rsum-SR);
  }
 	return lastViolation=violation;
} // violationBord2

int violationMagic() { // R==C
//  --------------
	int violation=0, *pr=xRectangle, Lsum=0, Rsum=0;

	// Row, column, diagonal constraints
	for (int i=0; i<R; ++i) {
	  int rSum=0, cSum=0, *pc=xRectangle+i, *pd=pr;
	  for (int j=0; j<R; ++j) { rSum+=*pr++; cSum+=*pc; pc+=R; } rowSum[i]=rSum; colSum[i]=cSum;
    violation+=abs(rSum-SR)+abs(cSum-SR);	Lsum+=pd[i]; Rsum+=pd[Rm1-i];
  }
  diagLsum=Lsum; diagRsum=Rsum; violation+=abs(Lsum-SR)+abs(Rsum-SR);
 	return lastViolation=violation;
} // violationMagic

int randomRestartMagic(int *X) {
//  ------------------
  for (int i=0; i<(RCfree-1); ++i) {
    const int ri=i+random_(RCfree-i), tmp=freeNum[i]; freeNum[i]=freeNum[ri]; freeNum[ri]=tmp;
  }
  int next=0;
  for (int i=0; i<RC; ++i) if (specific[i]==0) X[i]=freeNum[next++]; else X[i]=specific[i];
  return rectangleType==normalMagic ? violationMagic() : violationSemi();
} // randomRestartMagic

int randomRestartSym(int *X) {
//  ----------------
  for (int i=0; i<RC; ++i) X[i]=specific[i]; const int halfRCfree=RCfree/2;
  for (int i=0; i<(halfRCfree-1); ++i) {
    const int ri=i+random_(halfRCfree-i), tmp=freeNum[i]; freeNum[i]=freeNum[ri]; freeNum[ri]=tmp;
  }
  int next=0;
  if (odd) X[RCd2]=S1;
  for (int i=0; i<RC; ++i) if (i<opp[i]) {
    int v=freeNum[next++]; if (random_(2)) v=S2-v; X[i]=v; X[opp[i]]=S2-v;
  }
  return violationSemi();
} // randomRestartSym

int randomRestartBord2(int *X) {
//  ------------------
  for (int i=0; i<RC; ++i) X[i]=specific[i]; const int halfRCfree=RCfree/2;
  for (int i=0; i<(halfRCfree-1); ++i) {
    const int ri=i+random_(halfRCfree-i), tmp=freeNum[i]; freeNum[i]=freeNum[ri]; freeNum[ri]=tmp;
  }
  int next=0;
  for (int i=0; i<RC; ++i) {
    const int r=row[i], c=col[i]; int v=0;
    if ((c>1)&(c<Cm2)) { if ((r==0)|(r==Rm2)) v=1; } else if ((r>1)&(r<Rm2)) { if ((c==0)|(c==Cm2)) v=1; }
    if (borderType==none) {
      if (r==0) { if ((c<2)|(c>=Cm2)) v=1; } else if (r==1) { if ((c==1)|(c==Cm2)) v=1; }
      else if (r==Rm1) { if ((c==1)|(c==Cm2)) v=1; }
    } else if (borderType==cyclic) {
      if (r==0) { if (c<=1) v=1; else if (c==Cm2) v=1; } else if (r==1) { if (c==Cm2) v=1; }
      else if (r==Rm2) { if (c==0) v=1; else if (c>=Cm2) v=1; } else if (r==Rm1) { if (c==0) v=1; }
    } else { // opposed
      if ((r==0)|(r==Rm2)) { if ((c==0)|(c==1)|(c==Cm2)|(c==Cm1)) v=1; }
    }
    if (v!=0) { v=freeNum[next++]; if (random_(2)) v=S2-v; X[i]=v; X[opp[i]]=S2-v; }
  }
  //filePrint(X,stdout);
  return violationBord2(); 
} // randomRestartBord2

void doSwapMagic(int *X, const int i1, const int i2) {
//   -----------
	 const int tmp=xRectangle[i1]; xRectangle[i1]=xRectangle[i2]; xRectangle[i2]=tmp;
} // doSwapMagic

void doSwapSym(int *X, int i1, int i2) {
//   ---------
	int tmp=xRectangle[i1]; xRectangle[i1]=xRectangle[i2]; xRectangle[i2]=tmp; i1=opp[i1];
  if (i1!=i2) { i2=opp[i2]; tmp=xRectangle[i1]; xRectangle[i1]=xRectangle[i2]; xRectangle[i2]=tmp; }
} // doSwapSym

int randomSwapAndCalcMagic(int *X, const int swaps) {
//   --------------------
  for (int swap=0; swap<swaps; ++swap) {
    const int a=freeIx[random_(RCfree)], b=freeIx[random_(RCfree)];
    doSwapMagic(X,a,b); updateViolationMagic(a, b);
  }
  return lastViolation;
} // randomSwapAndCalcMagic

int swapDeltaSemi(int *X, const int i1, const int i2) {
//  -------------
	 const int r1=row[i1], c1=col[i1], r2=row[i2], c2=col[i2],
	           x1=xRectangle[i1], x2=xRectangle[i2], y=x2-x1; int x, delta=0;

	 if (r1!=r2) { // Different rows
    x=rowSum[r1]-SR; delta+=abs(x+y)-abs(x); x=rowSum[r2]-SR; delta+=abs(x-y)-abs(x);
	 }
	 if (c1!=c2) { // Different columns
    x=colSum[c1]-SC; delta+=abs(x+y)-abs(x); x=colSum[c2]-SC; delta+=abs(x-y)-abs(x);
 	}
 return delta;
} // swapDeltaSemi

int swapDeltaMagic(int *X, const int i1, const int i2) { // R==C
//  --------------
  const int r1=row[i1], c1=col[i1], r2=row[i2], c2=col[i2],
	          x1=xRectangle[i1], x2=xRectangle[i2], y=x2-x1; int x, delta=0;

	if (r1!=r2) { // Different rows
    x=rowSum[r1]-SR; delta+=abs(x+y)-abs(x); x=rowSum[r2]-SR; delta+=abs(x-y)-abs(x);
  }
	if (c1!=c2) { // Different columns
    x=colSum[c1]-SR; delta+=abs(x+y)-abs(x); x=colSum[c2]-SR; delta+=abs(x-y)-abs(x);
 	}

  if (r1==c1) {
    if (r2!=c2) { /* First is in diag 1 but second is not. */ x=diagLsum-SR; delta+=abs(x+y)-abs(x); }
  } else if (r2==c2) { /* First is not in diag 1 but second is. */ x=diagLsum-SR; delta+=abs(x-y)-abs(x); }

  if (r1==(Rm1-c1)) {
    if (r2!=(Rm1-c2)) { /* First is in diag 2 but second is not. */ x=diagRsum-SR; delta+=abs(x+y)-abs(x); }
  } else if (r2==(Rm1-c2)) { /* First is not in diag 2 but second is. */ x=diagRsum-SR; delta+=abs(x-y)-abs(x); }
	return delta;
} // swapDeltaMagic

int swapDeltaSym(int *X, const int i1, const int i2) {
//  ------------
	const int r1=row[i1], c1=col[i1], r2=row[i2], c2=col[i2],
	          x1=xRectangle[i1], x2=xRectangle[i2], y=x2-x1;
  int x, d, delta=0; const bool cPair=opp[i1]==i2;

  if (r1!=r2) { // Different rows
    x=rowSum[r1]-SR; d=abs(x+y)-abs(x); delta+=d;
    if (cPair) delta+=d; else { x=rowSum[r2]-SR; delta+=abs(x-y)-abs(x); }
	 }
	 if (c1!=c2) { // Different columns
    x=colSum[c1]-SC; d=abs(x+y)-abs(x); delta+=d;
    if (cPair) delta+=d; else { x=colSum[c2]-SC; delta+=abs(x-y)-abs(x); }
 	}
  return delta;
} // swapDeltaSym

int swapDeltaBord2(int *X, const int i1, const int i2) {
//  --------------
	const int r1=row[i1], c1=col[i1], r2=row[i2], c2=col[i2],
	          x1=xRectangle[i1], x2=xRectangle[i2], y=x2-x1;
  int x, d, delta=0; const bool cPair=opp[i1]==i2;

  if (r1!=r2) { // Different rows
    if ((r1<2)|(r1>=Rm2)) { x=rowSum[r1]-SR; d=abs(x+y)-abs(x); delta+=d; }
    if (cPair) delta+=d; else if ((r2<2)|(r2>=Rm2)) { x=rowSum[r2]-SR; delta+=abs(x-y)-abs(x); }
	}
	if (c1!=c2) { // Different columns
    if ((c1<2)|(c1>=Cm2)) { x=colSum[c1]-SC; d=abs(x+y)-abs(x); delta+=d; }
    if (cPair) delta+=d; else if ((c2<2)|(c2>=Cm2)) { x=colSum[c2]-SC; delta+=abs(x-y)-abs(x); }
 	}
  if ((R==C)&(borderType!=none)) {
    if (r1==c1) {
      if (r2!=c2) { /* First is in diag 1 but second is not. */ x=diagLsum-SR; delta+=abs(x+y)-abs(x); }
    } else if (r2==c2) { /* First is not in diag 1 but second is. */ x=diagLsum-SR; delta+=abs(x-y)-abs(x); }

    if (r1==(Rm1-c1)) {
      if (r2!=(Rm1-c2)) { /* First is in diag 2 but second is not. */ x=diagRsum-SR; delta+=abs(x+y)-abs(x); }
    } else if (r2==(Rm1-c2)) { /* First is not in diag 2 but second is. */ x=diagRsum-SR; delta+=abs(x-y)-abs(x); }
  }
  return delta;
} // swapDeltaBord2

void tabuSumsSemi(int *X, const int iterations, int *violation) {
//   ------------
  int it=0, viol=*violation;
  if (viol) {
    for (int i=RCm1; i>=0; --i) tabu[i]=0;
    while (it<iterations) {
		  int i1=-1, i2=-1, bestDelta=0, minR=-1, maxR=-1, minC=-1, maxC=-1,
          minRSum=INT_MAX, maxRSum=0, minCSum=INT_MAX, maxCSum=0;
      if (!fullSearch) {
        for (int i=0; i<R; ++i) {
          int q=rowSum[i]; if (q<minRSum) { minRSum=q; minR=i; } if (q>maxRSum) { maxRSum=q; maxR=i; }
        }
        for (int i=0; i<C; ++i) {
          int q=colSum[i]; if (q<minCSum) { minCSum=q; minC=i; } if (q>maxCSum) { maxCSum=q; maxC=i; }
        }
      }
      for (int i=RCfree_2; i>=0; --i) if (tabu[i]<=it) {
        int fi=freeIx[i]; bool ido=fullSearch;
        if (!ido) { int iR=row[fi], iC=col[fi]; ido=(iR==minR)||(iR==maxR)||(iC==minC)||(iC==maxC); }
        if (ido) for (int j=i+1; j<RCfree; ++j) if (tabu[j]<=it) {
          int delta=swapDeltaSemi(X, fi, freeIx[j]);
          if ((i1<0)|(bestDelta>delta)) { i1=i; i2=j; bestDelta=delta; }
				  else if ((delta==bestDelta)&&random_(2)) { i1=i; i2=j; }
			  } // for (int j ...
		  } // for (int i ...
      ++it;
      if (i1>=0) {
		    doSwapMagic(X, freeIx[i1], freeIx[i2]); tabu[i1]=it+tabuLen;
        if (Nfree>5) tabu[i2]=it+tabuLen;
        if ((viol=updateViolationSemi(freeIx[i1],freeIx[i2]))==0) break;
      }
	   } // while (it ...
  } // if (viol ...
  *violation=viol;
} // tabuSumsSemi

void tabuSumsMagic(int *X, const int iterations, int *violation) { // R==C
//   -------------
  int it=0, viol=*violation;
  if (viol) {
    int zeroBD=0; const int noItGainX=Nfree<10 ? 1 : Nfree<20 ? 2 :
        Nfree<30 ? 3 : Nfree<40 ? 4 : Nfree<50 ? 5 : 6;
    for (int i=RCm1; i>=0; --i) tabu[i]=0;
    while (it<iterations) {
		  int i1=-1, i2=-1, bestDelta=0, minR=-1, maxR=-1, minC=-1, maxC=-1,
          minRSum=INT_MAX, maxRSum=0, minCSum=INT_MAX, maxCSum=0,
          minSum=INT_MAX, maxSum=0; bool doL=F, doR=F;
      if (!fullSearch) {
        for (int i=0; i<R; ++i) {
          int q=rowSum[i]; if (q<minRSum) { minRSum=q; minR=i; } if (q>maxRSum) { maxRSum=q; maxR=i; }
          q=colSum[i]; if (q<minCSum) { minCSum=q; minC=i; } if (q>maxCSum) { maxCSum=q; maxC=i; }
        }
        minSum=min(minRSum, minCSum); maxSum=max(maxRSum, maxCSum);
        doL=(diagLsum<=minSum)|(diagLsum>=maxSum); doR=(diagRsum<=minSum)|(diagRsum>=maxSum);
      }
      for (int i=RCfree_2; i>=0; --i) if (tabu[i]<=it) {
        const int fi=freeIx[i]; bool ido=fullSearch;
        if (!ido) {
          const int iR=row[fi], iC=col[fi];
          ido=(iR==minR)||(iR==maxR)||(iC==minC)||(iC==maxC)||(doL&&diagL[fi])||(doR&&diagR[fi]);
        }
        if (ido) for (int j=i+1; j<RCfree; ++j) if (tabu[j]<=it) {
          const int delta=swapDeltaMagic(X, fi, freeIx[j]);

          if ((i1<0)|(bestDelta>delta)) { i1=i; i2=j; bestDelta=delta; }
				  else if ((delta==bestDelta)&&random_(2)) { i1=i; i2=j; }
			  } // for (int j ...
		  } // for (int i ...
      ++it;
      if (i1>=0) {
        if (bestDelta==0) {
          if (++zeroBD>(2*tabuLen)) {
            viol=randomSwapAndCalcMagic(X,noItGainX); if (viol==0) break; zeroBD=0;
          } else {
            const int a=freeIx[i1], b=freeIx[i2]; doSwapMagic(X,a,b); viol=updateViolationMagic(a,b);
          }
        } else {
          const int a=freeIx[i1], b=freeIx[i2];
          doSwapMagic(X,a,b); viol=updateViolationMagic(a,b); if (viol==0) break; zeroBD=0;
        }   
        tabu[i1]=it+tabuLen; if (Nfree>5) tabu[i2]=it+tabuLen;
      }
	   } // while (it ...
  } // if (viol ...
  *violation=viol;
} // tabuSumsMagic

void tabuSumsSym(int *X, const int iterations, int *violation) {
//   -----------
  int it=0, viol=*violation;
  if (viol) {
    for (int i=RCm1; i>=0; --i) tabu[i]=0;
    while (it<iterations) {
		  int i1=-1, i2=-1, bestDelta=0, minR=-1, maxR=-1, minC=-1, maxC=-1,
          minRSum=INT_MAX, maxRSum=0, minCSum=INT_MAX, maxCSum=0;
      if (!fullSearch) {
        for (int i=0; i<R; ++i) {
          const int q=rowSum[i]; if (q<minRSum) { minRSum=q; minR=i; } if (q>maxRSum) { maxRSum=q; maxR=i; }
        }
        for (int i=0; i<C; ++i) {
          const int q=colSum[i]; if (q<minCSum) { minCSum=q; minC=i; } if (q>maxCSum) { maxCSum=q; maxC=i; }
        }
      }
      for (int i=RCfree_2; i>=0; --i) if (tabu[i]<=it) {
        const int fi=freeIx[i]; bool ido=fullSearch;
        if (!ido) {
          const int iR=row[fi], iC=col[fi]; ido=(iR==minR)||(iR==maxR)||(iC==minC)||(iC==maxC);
        }
        if (ido) for (int j=i+1; j<RCfree; ++j) if (tabu[j]<=it) {
          const int delta=swapDeltaSym(X, fi, freeIx[j]);
          if ((i1<0)|(bestDelta>delta)) { i1=i; i2=j; bestDelta=delta; }
				  else if ((delta==bestDelta)&&random_(2)) { i1=i; i2=j; }
			  } // for (int j ...
		  } // for (int i ...
      ++it;
      if (i1>=0) {
		    doSwapSym(X, freeIx[i1], freeIx[i2]); tabu[i1]=it+tabuLen; tabu[RCfree_1-i2]=it+tabuLen;
        if (Nfree>5) { tabu[i2]=it+tabuLen; tabu[RCfree_1-i1]=it+tabuLen; }
        if ((viol=updateViolationSym(freeIx[i1],freeIx[i2]))==0) break;
      }
	   } // while (it ...
  } // if (viol ...
  *violation=viol;
} // tabuSumsSym

void tabuSumsBord2(int *X, const int iterations, int *violation) {
//   -------------
  int it=0, viol=*violation;
  if (viol) {
    for (int i=RCm1; i>=0; --i) tabu[i]=0;
    while (it<iterations) {
		  int i1=-1, i2=-1, bestDelta=0, minR=-1, maxR=-1, minC=-1, maxC=-1, minLSum=INT_MAX, maxLSum=0, 
          minRSum=INT_MAX, maxRSum=0, minCSum=INT_MAX, maxCSum=0, minSum=INT_MAX, maxSum=0; bool doL=F, doR=F;
      if (!fullSearch) {
        for (int i=Rm2; i<R; ++i) {
          const int q=rowSum[i]; if (q<minRSum) { minRSum=q; minR=i; } if (q>maxRSum) { maxRSum=q; maxR=i; }
        }
        for (int i=Cm2; i<C; ++i) {
          const int q=colSum[i]; if (q<minCSum) { minCSum=q; minC=i; } if (q>maxCSum) { maxCSum=q; maxC=i; }
        }
        minSum=min(minRSum, minCSum); maxSum=max(maxRSum, maxCSum);
        if ((R==C)&(borderType!=none)) {
          doL=(diagLsum<=minSum)|(diagLsum>=maxSum); doR=(diagRsum<=minSum)|(diagRsum>=maxSum);
        }
      }
      for (int i=RCfree_2; i>=0; --i) if (tabu[i]<=it) {
        const int fi=freeIx[i]; bool ido=fullSearch;
        if (!ido) {
          const int iR=row[fi], iC=col[fi];
          ido=(iR==minR)||(iR==maxR)||(iC==minC)||(iC==maxC)||(doL&&diagL[fi])||(doR&&diagR[fi]);
        }
        if (ido) for (int j=i+1; j<RCfree; ++j) if (tabu[j]<=it) {
          const int delta=swapDeltaBord2(X, fi, freeIx[j]);
          if ((i1<0)|(bestDelta>delta)) { i1=i; i2=j; bestDelta=delta; }
				  else if ((delta==bestDelta)&&random_(2)) { i1=i; i2=j; }
			  } // for (int j ...
		  } // for (int i ...
      ++it;
      if (i1>=0) {
		    doSwapSym(X, freeIx[i1], freeIx[i2]); tabu[i1]=it+tabuLen; tabu[RCfree_1-i2]=it+tabuLen;
        if (Nfree>5) { tabu[i2]=it+tabuLen; tabu[RCfree_1-i1]=it+tabuLen; }
        if ((viol=updateViolationBord2(freeIx[i1],freeIx[i2]))==0) break;
      }
	   } // while (it ...
  } // if (viol ...
  *violation=viol;
} // tabuSumsBord2
//==================================================== main ================================================

typedef bool (*t_bNoParm)(); typedef int  (*t_iP)(int *); typedef void  (*t_iPcIiP)(int *, const int, int *);

struct t_type { const char *name; t_bNoParm initialize; t_iP randomRestart; t_iPcIiP  tabuSums; };
const int numTypes=lastType-firstType+1;
t_type rectangleTypes[numTypes]=
{
  { "Semimagic",   initializeSemi,  randomRestartMagic, tabuSumsSemi  },
  { "Magic",       initializeMagic, randomRestartMagic, tabuSumsMagic },
  { "Associative", initializeSym,   randomRestartSym,   tabuSumsSym   },
  { "Bordered",    initializeSym,   randomRestartBord2, tabuSumsBord2 }
};

bool checkRC(const int nr, const int nc) {
//   -------
  rectangleType=Bordered; int lo=odd ? (nr==nc) ? 5 : 7 : (nr==nc) ? 8 : 6;
  if (nr==nc) {
    if (nr==2) { printf("\aThere is no magic square of order %d.\n", nr); return F; }
    else if (nr<lo) { printf("\aThere is no double bordered magic square of order %d.\n",  nr); return F; }}
  else if ((nr&1)!=(nc&1)) { printf("\aThere is no odd, even magic rectangle.\n"); return F; }
  else if ((nr<lo)|(nc<lo)) {
    printf("\aThere is no double bordered magic rectangle of order %d %d.\n", nr, nc); return F; }
  return T;
}// checkRC

t_type *getStartType() {
//      ------------
  allZero=T; int tmp=(Rend<=Cend)?Rend:Cend; if (odd) ++tmp; bool deven=((tmp&3)==0);
  if (Rend==Cend) {
    R=odd ? deven ? 3 : (borderType==opposed) ? 5 : 1 : deven ? 4 : 6; C=R;
    rectangleType=(R==6)?normalMagic:random_(2)?normalMagic:Associative;
  } else {
    int t=odd ? deven ? 3 : 5 : deven ? 4 : 2;
    if (Rend<=Cend) { R=t; C=Cend-Rend+t; } else { R=Rend-Cend+t; C=t; }
    rectangleType=normalSemiMagic;
  }
  return &rectangleTypes[rectangleType];
} // getStartType

int makeRectangle(const int iterations) {
//  -------------
  int violation, loop=100;
  while (loop--) {
    t_type *sqType=getStartType();
    if (sqType->initialize()) {
      violation=sqType->randomRestart(xRectangle);
      sqType->tabuSums(xRectangle, iterations, &violation);
    }
    if (violation==0) {
      rectangleType=Bordered;
      sqType=&rectangleTypes[rectangleType]; allZero=F;
      for (int r=R+4; r<=Rend; r+=4) {
        R+=4; C+=4;
        if (sqType->initialize()) {
          int loop2=50;
          while (loop2--) {
            violation=sqType->randomRestart(xRectangle);
            sqType->tabuSums(xRectangle, iterations, &violation);
            if (violation==0) break;
          }
        }
        if (violation!=0) break;
      }
      if (violation==0) break;
    }
  }
  return violation;
} //makeRectangle

int getTypeAndHowMany(borderRectangleType *borderType) {
//  -----------------
  int tmp, max;
  if (Rend==5) { tmp=getInt("Border rectangles, (1 none, 2 cyclic)? "); max=2; }
  else { tmp=getInt("Border rectangles, (1 none, 2 cyclic, 3 opposed)? ");  max=3; }
  if (tmp<=0) tmp=1; else if (tmp>max) tmp=max; *borderType=(borderRectangleType)(tmp-1);
  int howMany=getInt("How many? "); if (howMany<=0) howMany=1; return howMany;
} // getTypeAndHowMany

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

int main() {
//  ----
  bool another=T, inputOrder=T, ok=F; seed_rand();
  if (openDir()) {
    do { 
      if (inputOrder) { printf("\nOrder? "); if (!getInts(&Rend, &Cend, -1)) break; } odd=Rend&1;
      const int howMany=getTypeAndHowMany(&borderType);
      if (checkRC(Rend, Cend)&&initializeStart()) {
        const int iterations=(Rend==Cend)?Rend*500:10000;
        FILE *wfp=openOutput();
        if (wfp!=NULL) {
          int good=0, runs=(howMany>1000)?howMany:1000; time_t startTime=time(NULL); fullSearch=F;
	        for (int i=0; i<runs; ++i) {
            int violation=makeRectangle(iterations);
            if (violation==0) {
              ok=T; printf("%d\n", ++good);
              if (!filePrint(xRectangle, wfp)) { reportError("Output failed"); break; } fflush(wfp);
              if (good==howMany) break;
            } // if (violation==0)
	        } // for (int i ...
          printf("Number of %ss: %d\n", (Rend==Cend) ? "square" : "rectangle", good);
          printElapsedTime(stdout, startTime);
          fclose(wfp); if (good==0) remove(OutputFileName);
        } // if (wfp ...
        freeStore();
      } // if (checkRC ....
      printf("\nAnother order? input y (yes), n (no) or the order: ");
      if (getYorInts(&Rend, &Cend)) inputOrder=(Rend<0); else break;
    } while (another);
  }
  printf("\nHit return to close the console.");
  while (T) if (getchar()=='\n') break; return ok ? EXIT_SUCCESS : EXIT_FAILURE;
} // main