#include <stdio.h>
#include "sparse.h"

#define ARRAY_SIZE 9

static int dumperrors=0;

static double tempval;
static int phase=0;

#define correctval(r,c) (tempval=CORRECTVAL(r,c),tempval)
#define negval(r,c) (tempval=NEGVAL(r,c),tempval)

static double CORRECTVAL(int r,int c) {
  return ((r+1) + (0.1*(c+1)));
}

static double NEGVAL(int r,int c) {
  return -((r+1) + (0.1*(c+1)));
}

static void SetupTridiag(SparseArray Tridiag) {
  int row;
  int col;
  int diff;

  /* Set up a tridiagonal using a doubly-nested loop */
  for (row=0; row<ARRAY_SIZE; row++) {
    for (col=0; col<ARRAY_SIZE; col++) {
      diff = row-col;
      if (diff >= -1 && diff <= 1) {
	Store(Tridiag,row,col,correctval(row,col));
      }
    }
  }
}


static void SetupTridiag2(SparseArray Tridiag) {
  int rowcol;

  /* Set up a tridiagonal using a singly-nested loop */
  Store(Tridiag,ARRAY_SIZE-1,ARRAY_SIZE-1,correctval(ARRAY_SIZE-1,ARRAY_SIZE-1));
  for (rowcol=0; rowcol<ARRAY_SIZE-1; rowcol++) {
    Store(Tridiag,rowcol+1,rowcol,correctval(rowcol+1,rowcol));
    Store(Tridiag,rowcol,rowcol+1,correctval(rowcol,rowcol+1));
    Store(Tridiag,rowcol,rowcol,correctval(rowcol,rowcol));
  }
}


void Dump(SparseArray X,char *subname) {
  char name[80];
  FILE* outfile;

  if (phase==0) {
    sprintf(name,"dump.normal.%s.out",subname);
  } else {
    sprintf(name,"dump.unusual.%s.out",subname);
  }
  outfile = fopen(name,"w");
  WriteDense(X,outfile);
  fclose(outfile);
}


void TestRead(SparseArray X) {
  int row;
  int col;
  int diff;
  double val;
  int problem;

  problem = 0;

  /* Make sure that values in the tridiagonal are correct using a
     reverse doubly-nested loop */
  for (col=ARRAY_SIZE-1; col>=0; col--) {
    for (row=ARRAY_SIZE-1; row>=0; row--) {
      val = Read(X,row,col);
      diff = row-col;
      if (diff >= -1 && diff <= 1) {
	problem |= (val != correctval(row,col));
      } else {
	problem |= (val != 0.0);
      }
    }
  }

  if (problem) {
    fprintf(stderr,"ERROR: TestRead() failed\n");
    if (dumperrors) {
      Dump(X,"Read");
    }
  } else {
    fprintf(stderr,"       TestRead() passed\n");
  }
}


void TestStore(SparseArray X) {
  int rowcol;
  int row;
  int col;
  int diff;
  int problem;
  double val;

  problem = 0;

  /* Negate the main diagonal using a reverse single loop */
  for (rowcol=ARRAY_SIZE-1; rowcol>=0; rowcol--) {
    Store(X,rowcol,rowcol,negval(rowcol,rowcol));
  }

  /* Check that the main diagonal is properly negated using a
     doubly-nested loop */
  for (row=0; row<ARRAY_SIZE; row++) {
    for (col=0; col<ARRAY_SIZE; col++) {
      val = Read(X,row,col);
      diff = row-col;
      if (diff == -1 || diff == 1) {
	problem |= (val != correctval(row,col));
      } else if (diff == 0) {
	problem |= (val != negval(row,col));
      } else {
	problem |= (val != 0.0);
      }
    }
  }

  if (problem) {
    if (dumperrors) {
      Dump(X,"Store-halfway");
    }
  }

  /* Fix the main diagonal using a forwards loop */
  for (rowcol=0; rowcol<ARRAY_SIZE; rowcol++) {
    Store(X,rowcol,rowcol,correctval(rowcol,rowcol));
  }

  /* Check that it's correct using a reverse doubly-nested loop */
  for (col=ARRAY_SIZE-1; col>=0; col--) {
    for (row=ARRAY_SIZE-1; row>=0; row--) {
      val = Read(X,row,col);
      diff = row-col;
      if (diff >= -1 && diff <= 1) {
	problem |= (val != correctval(row,col));
      } else {
	problem |= (val != 0.0);
      }
    }
  }

  if (problem) {
    fprintf(stderr,"ERROR: TestStore() failed\n");
    if (dumperrors) {
      Dump(X,"Store");
    }
  } else {
    fprintf(stderr,"       TestStore() passed\n");
  }
}


void TestAccess(SparseArray X) {
  int rowcol;
  int row;
  int col;
  int diff;
  int problem;
  double val;
  double* valptr;

  problem = 0;

  /* Negate the main diagonal using a reverse single loop */
  for (rowcol=ARRAY_SIZE-1; rowcol>=0; rowcol--) {
    valptr = Access(X,rowcol,rowcol);
    if (valptr == NULL) {
      problem = 1;
    } else {
      problem |= (*valptr != correctval(rowcol,rowcol));
      *valptr = negval(rowcol,rowcol);
    }
  }

  if (problem) {
    fprintf(stderr,"ERROR: simple TestAccess() failed\n");
    if (dumperrors) {
      Dump(X,"Access-simple");
    }
  } else {
    fprintf(stderr,"       simple TestAccess() passed\n");
  }

  /* Check that the main diagonal is properly negated using a
     doubly-nested loop */
  for (row=0; row<ARRAY_SIZE; row++) {
    for (col=0; col<ARRAY_SIZE; col++) {
      val = Read(X,row,col);
      valptr = Access(X,row,col);
      diff = row-col;
      if (diff == -1 || diff == 1) {
	problem |= (val != correctval(row,col));
	problem |= (valptr == NULL);
      } else if (diff == 0) {
	problem |= (val != negval(row,col));
	problem |= (valptr == NULL);
      } else {
	problem |= (val != 0.0);
	problem |= (valptr != NULL);
      }
    }
  }

  if (problem) {
    if (dumperrors) {
      Dump(X,"Access-halfway");
    }
  }

  /* Fix the main diagonal using a forwards loop */
  for (rowcol=0; rowcol<ARRAY_SIZE; rowcol++) {
    valptr = Access(X,rowcol,rowcol);
    if (valptr == NULL) {
      problem = 1;
    } else {
      problem |= (*valptr != negval(rowcol,rowcol));
      *valptr = correctval(rowcol,rowcol);
    }
  }

  /* Check that it's correct using a reverse doubly-nested loop */
  for (col=ARRAY_SIZE-1; col>=0; col--) {
    for (row=ARRAY_SIZE-1; row>=0; row--) {
      val = Read(X,row,col);
      valptr = Access(X,row,col);
      diff = row-col;
      if (diff >= -1 && diff <= 1) {
	problem |= (val != correctval(row,col));
	problem |= (valptr == NULL);
      } else {
	problem |= (val != 0.0);
	problem |= (valptr != NULL);
      }
    }
  }

  if (problem) {
    fprintf(stderr,"ERROR: TestAccess() failed\n");
    if (dumperrors) {
      Dump(X,"Access");
    }
  } else {
    fprintf(stderr,"       TestAccess() passed\n");
  }
}


void TestURV(SparseArray X) {
  int row;
  int col;
  int diff;  
  double val;
  int problem;

  problem = 0;

  val = GetURV(X);
  problem |= (val != 0.0);
  SetURV(X,-1.0);
  val = GetURV(X);
  problem |= (val != -1.0);
  
  /* Check that it's correct using a reverse doubly-nested loop */
  for (col=ARRAY_SIZE-1; col>=0; col--) {
    for (row=ARRAY_SIZE-1; row>=0; row--) {
      val = Read(X,row,col);
      diff = row-col;
      if (diff >= -1 && diff <= 1) {
	problem |= (val != correctval(row,col));
      } else {
	problem |= (val != -1.0);
      }
    }
  }

  if (problem) {
    if (dumperrors) {
      Dump(X,"URV-halfway");
    }
  }

  SetURV(X,0.0);
  val = GetURV(X);
  problem |= (val != 0.0);

  /* Check that it's correct using a reverse doubly-nested loop */
  for (col=0; col<ARRAY_SIZE; col++) {
    for (row=0; row<ARRAY_SIZE; row++) {
      val = Read(X,row,col);
      diff = row-col;
      if (diff >= -1 && diff <= 1) {
	problem |= (val != correctval(row,col));
      } else {
	problem |= (val != 0.0);
      }
    }
  }

  if (problem) {
    fprintf(stderr,"ERROR: TestURV() failed\n");
    if (dumperrors) {
      Dump(X,"URV");
    }
  } else {
    fprintf(stderr,"       TestURV() passed\n");
  }
}


void TestDenseIO(SparseArray X,char *name) {
  int row;
  int col;
  double val1;
  double val2;
  int problem;
  
  problem = 0;
  {
    FILE *outfile;

    outfile = fopen(name,"w");
    WriteDense(X,outfile);
    fclose(outfile);
  }
  {
    FILE *infile;
    SparseArray Y;
    
    infile = fopen(name,"r");

    Y = ReadDense(infile,ARRAY_SIZE,ARRAY_SIZE,0.0);
    fclose(infile);

    for (col=ARRAY_SIZE-1; col>=0; col--) {
      for (row=ARRAY_SIZE-1; row>=0; row--) {
	val1 = Read(X,row,col);
	val2 = Read(Y,row,col);
	problem |= (val1 != val2);
      }
    }
    if (problem) {
      fprintf(stderr,"ERROR: TestDenseIO() failed\n");
      if (dumperrors) {
	Dump(Y,"Dense");
      }
    } else {
      fprintf(stderr,"       TestDenseIO() passed\n");
    }
  }
}


void TestSparseIO(SparseArray X,char *name) {
  int row;
  int col;
  double val1;
  double val2;
  int problem;

  problem=0;
  {
    FILE* outfile;

    outfile = fopen(name,"w");
    WriteSparse(X,outfile);
    fclose(outfile);
  }
  {
    FILE* infile;
    SparseArray Y;

    infile = fopen(name,"r");
    Y = ReadSparse(infile);
    fclose(infile);

    for (col=ARRAY_SIZE-1; col>=0; col--) {
      for (row=ARRAY_SIZE-1; row>=0; row--) {
	val1 = Read(X,row,col);
	val2 = Read(Y,row,col);
	problem |= (val1 != val2);
      }
    }
    if (problem) {
      fprintf(stderr,"ERROR: TestSparseIO() failed\n");
      if (dumperrors) {
	Dump(Y,"Sparse");
      }
    } else {
      fprintf(stderr,"       TestSparseIO() passed\n");
    }
  }
}


void TestIterateRight(SparseArray X) {
  int row;
  int col;
  int problem;
  int i;
  int lo;
  int hi;

  problem = 0;

  for (row=ARRAY_SIZE-1; row>=0; row--) {
    col = SetCursorRowStart(X,row);
    lo = -1;
    hi = 1;
    if (row == 0) {
      lo = 0;
    } else if (row == ARRAY_SIZE-1) {
      hi = 0;
    }
    for (i=lo;i<=hi;i++) {
      problem |= (col != row+i);
      col = NextInRow(X);
    }
    problem |= (col != -1);
  }

  if (problem) {
    fprintf(stderr,"ERROR: TestIterateRight() failed\n");
  } else {
    fprintf(stderr,"       TestIterateRight() passed\n");
  }
}


void TestIterateLeft(SparseArray X) {
  int row;
  int col;
  int problem;
  int i;
  int lo;
  int hi;

  problem = 0;

  for (row=ARRAY_SIZE-1; row>=0; row--) {
    col = SetCursorRowEnd(X,row);
    lo = -1;
    hi = 1;
    if (row == 0) {
      lo = 0;
    } else if (row == ARRAY_SIZE-1) {
      hi = 0;
    }
    for (i=hi;i>=lo;i--) {
      problem |= (col != row+i);
      col = PrevInRow(X);
    }
    problem |= (col != -1);
  }

  if (problem) {
    fprintf(stderr,"ERROR: TestIterateLeft() failed\n");
  } else {
    fprintf(stderr,"       TestIterateLeft() passed\n");
  }
}


void TestIterateDown(SparseArray X) {
  int row;
  int col;
  int problem;
  int i;
  int lo;
  int hi;

  problem = 0;

  for (col=ARRAY_SIZE-1; col>=0; col--) {
    row = SetCursorColStart(X,col);
    lo = -1;
    hi = 1;
    if (col == 0) {
      lo = 0;
    } else if (col == ARRAY_SIZE-1) {
      hi = 0;
    }
    for (i=lo;i<=hi;i++) {
      problem |= (row != col+i);
      row = NextInCol(X);
    }
    problem |= (row != -1);
  }

  if (problem) {
    fprintf(stderr,"ERROR: TestIterateDown() failed\n");
  } else {
    fprintf(stderr,"       TestIterateDown() passed\n");
  }
}


void TestIterateUp(SparseArray X) {
  int row;
  int col;
  int problem;
  int i;
  int lo;
  int hi;

  problem = 0;

  for (col=ARRAY_SIZE-1; col>=0; col--) {
    row = SetCursorColEnd(X,col);
    lo = -1;
    hi = 1;
    if (col == 0) {
      lo = 0;
    } else if (col == ARRAY_SIZE-1) {
      hi = 0;
    }
    for (i=hi;i>=lo;i--) {
      problem |= (row != col+i);
      row = PrevInCol(X);
    }
    problem |= (row != -1);
  }

  if (problem) {
    fprintf(stderr,"ERROR: TestIterateUp() failed\n");
  } else {
    fprintf(stderr,"       TestIterateUp() passed\n");
  }
}


void IterateAround(SparseArray X) {
  int rowcol;
  int problem;
  int row;
  int col;

  problem = 0;

  row = 0;
  col = SetCursorRowStart(X,row);
  for (rowcol=0; rowcol<ARRAY_SIZE-1; rowcol++) {
    problem |= (row != col);
    col = NextInRow(X);
    problem |= (col != row+1);
    row = NextInCol(X);
  }
  problem |= (row != ARRAY_SIZE-1);
  row = NextInCol(X);
  problem |= (row != -1);

  col = ARRAY_SIZE-1;
  row = SetCursorColEnd(X,col);
  for (rowcol=ARRAY_SIZE-1; rowcol>0; rowcol--) {
    problem |= (row != col);
    col = PrevInRow(X);
    problem |= (col != row-1);
    row = PrevInCol(X);
  }
  problem |= (col != 0);
  col = PrevInCol(X);
  problem |= (col != -1);

  if (problem) {
    fprintf(stderr,"ERROR: IterateAround() failed\n");
  } else {
    fprintf(stderr,"       IterateAround() passed\n");
  }
}


void TestDensity(SparseArray X) {
  double val;
  int problem;

  problem = 0;

  val = Density(X);
  if (problem) {
    fprintf(stderr,"ERROR: simple TestDensity() failed: %f\n",val);
  } else {
    fprintf(stderr,"       simple TestDensity() passed\n");
  }
  problem |= (val != (double)(2*2 + 3*(ARRAY_SIZE-2))/(ARRAY_SIZE*ARRAY_SIZE));
  Store(X,0,0,-1.0);
  Store(X,0,ARRAY_SIZE-1,-2.0);
  Store(X,ARRAY_SIZE-1,0,-3.0);
  Store(X,ARRAY_SIZE-1,ARRAY_SIZE-1,-4.0);
  val = Density(X);
  problem |= (val != (double)(3*ARRAY_SIZE)/(ARRAY_SIZE*ARRAY_SIZE));

  if (problem) {
    fprintf(stderr,"ERROR: TestDensity() failed: %f\n",val);
    if (dumperrors) {
      Dump(X,"Density");
    }
  } else {
    fprintf(stderr,"       TestDensity() passed\n");
  }
}


void TestNumRepresented(SparseArray X) {
  int val;
  int problem;

  problem = 0;

  val = NumRepresented(X);
  if (problem) {
    fprintf(stderr,"ERROR: simple TestNumRepresented() failed: %d\n",val);
  } else {
    fprintf(stderr,"       simple TestNumRepresented() passed\n");
  }
  problem |= (val != 2*2 + 3*(ARRAY_SIZE-2));
  Store(X,0,0,-1.0);
  Store(X,0,ARRAY_SIZE-1,-2.0);
  Store(X,ARRAY_SIZE-1,0,-3.0);
  Store(X,ARRAY_SIZE-1,ARRAY_SIZE-1,-4.0);
  val = NumRepresented(X);
  problem |= (val != 3*ARRAY_SIZE);

  if (problem) {
    fprintf(stderr,"ERROR: TestNumRepresented() failed: %d\n",val);
    if (dumperrors) {
      Dump(X,"NumRepresented");
    }
  } else {
    fprintf(stderr,"       TestNumRepresented() passed\n");
  }
}


int main(int argc,char *argv[]) {
  SparseArray T1A;
  SparseArray T1B;
  SparseArray T2A;
  SparseArray T2B;
  SparseArray T3A;
  SparseArray T3B;
  SparseArray T4A;
  SparseArray T4B;
  SparseArray T5A;
  SparseArray T5B;
  SparseArray T6A;
  SparseArray T6B;
  SparseArray T7A;
  SparseArray T7B;
  SparseArray T8A;
  SparseArray T8B;
  SparseArray T9A;
  SparseArray T9B;
  SparseArray TaA;
  SparseArray TaB;
  SparseArray TbA;
  SparseArray TbB;
  SparseArray TcA;
  SparseArray TcB;
  SparseArray TdA;
  SparseArray TdB;

  if (argc > 1) {
    dumperrors = 1;
  }

  T1A = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T1B = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T2A = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T2B = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T3A = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T3B = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T4A = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T4B = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T5A = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T5B = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T6A = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T6B = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T7A = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T7B = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T8A = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T8B = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T9A = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  T9B = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  TaA = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  TaB = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  TbA = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  TbB = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  TcA = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  TcB = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  TdA = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);
  TdB = NewSparseArray(ARRAY_SIZE,ARRAY_SIZE,0.0);

  fprintf(stderr,"Setting up all normal arrays...");
  fflush(stderr);
  SetupTridiag(T1A);
  SetupTridiag(T2A);
  SetupTridiag(T3A);
  SetupTridiag(T4A);
  SetupTridiag(T5A);
  SetupTridiag(T6A);
  SetupTridiag(T7A);
  SetupTridiag(T8A);
  SetupTridiag(T9A);
  SetupTridiag(TaA);
  SetupTridiag(TbA);
  SetupTridiag(TcA);
  SetupTridiag(TdA);
  fprintf(stderr,"done.\n");

  fprintf(stderr,"Testing all operations on normal arrays...\n");
  TestRead(T1A);
  TestStore(T2A);
  TestAccess(T3A);
  TestURV(T4A);
  TestIterateRight(T7A);
  TestIterateLeft(T8A);
  TestIterateDown(T9A);
  TestIterateUp(TaA);
  IterateAround(TbA);
  TestDensity(TcA);
  TestNumRepresented(TdA);
  TestDenseIO(T5A,"dense1.out");
  TestSparseIO(T6A,"sparse1.out");


  fprintf(stderr,"Deleting all normal arrays...");
  fflush(stderr);
  DeleteSparseArray(T1A);
  DeleteSparseArray(T2A);
  DeleteSparseArray(T3A);
  DeleteSparseArray(T4A);
  DeleteSparseArray(T5A);
  DeleteSparseArray(T6A);
  DeleteSparseArray(T7A);
  DeleteSparseArray(T8A);
  DeleteSparseArray(T9A);
  DeleteSparseArray(TaA);
  DeleteSparseArray(TbA);
  DeleteSparseArray(TcA);
  DeleteSparseArray(TdA);
  fprintf(stderr,"done.\n");

  
  phase=1;
  fprintf(stderr,"Setting up all unusual arrays...");
  fflush(stderr);
  SetupTridiag(T1B);
  SetupTridiag(T2B);
  SetupTridiag(T3B);
  SetupTridiag(T4B);
  SetupTridiag(T5B);
  SetupTridiag(T6B);
  SetupTridiag(T7B);
  SetupTridiag(T8B);
  SetupTridiag(T9B);
  SetupTridiag(TaB);
  SetupTridiag(TbB);
  SetupTridiag(TcB);
  SetupTridiag(TdB);
  fprintf(stderr,"done.\n");

  fprintf(stderr,"Testing all operations on unusual arrays...\n");
  TestRead(T1B);
  TestStore(T2B);
  TestAccess(T3B);
  TestURV(T4B);
  TestIterateRight(T7B);
  TestIterateLeft(T8B);
  TestIterateDown(T9B);
  TestIterateUp(TaB);
  IterateAround(TbB);
  TestDensity(TcB);
  TestNumRepresented(TdB);
  TestDenseIO(T5B,"dense2.out");
  TestSparseIO(T6B,"sparse2.out");

  fprintf(stderr,"Deleting all unusual arrays...");
  fflush(stderr);
  DeleteSparseArray(T1B);
  DeleteSparseArray(T2B);
  DeleteSparseArray(T3B);
  DeleteSparseArray(T4B);
  DeleteSparseArray(T5B);
  DeleteSparseArray(T6B);
  DeleteSparseArray(T7B);
  DeleteSparseArray(T8B);
  DeleteSparseArray(T9B);
  DeleteSparseArray(TaB);
  DeleteSparseArray(TbB);
  DeleteSparseArray(TcB);
  DeleteSparseArray(TdB);
  fprintf(stderr,"done.\n");
}