#include "TSP.h"

void Perturb4Opt(int);
void TwoOpt(int, int, int, int);
void LocalSearch();
void FDD_Diversification();
bool AcceptanceTest();

int main(int argc,char** argv) {
  if (PreProcess(argc,argv) == 1)
    return 1;

  MAX_ITERATIONS = n*n; // n^2
  ReadConfigurationFile(ConfigurationFile);
  ReadInstanceFile(InstanceFile);
  LOGGER = new VizLogger("TSP",true,InstanceFile,ConfigurationFile,RunID,n,"ILS","ILS for TSP",MAX_ITERATIONS);

  CreateGreedySolution();

  LocalSearch();
  UpdateStatistics();

  // We stop if gap = 0, we know the BK values :), this is to save time for our experiments...
  for (Cur_Iter=1; Cur_Iter<=MAX_ITERATIONS && gap > 0; Cur_Iter++) { // O(MAX_ITERATIONS * k * n)
    Perturb4Opt(PERTURBATION_STRENGTH); // O(n)
    LocalSearch(); // O(k*n)

    LOGGER->LogCurrentSolution(solution,1,Cur_OV); // record solution after perturb+LS phase
    char str_Perturb[255],str_Accept[255];
    sprintf(str_Perturb,"%d",PERTURBATION_STRENGTH);
    sprintf(str_Accept,"%d",0);
    LOGGER->LogAlgorithmSpecific(str_Perturb,str_Accept,"");

    bool Accepted = AcceptanceTest(); // O(n), record solution after acceptance phase...
    LOGGER->LogCurrentSolution(solution,1,Cur_OV);
    sprintf(str_Accept,"%d",Accepted?2:0);
    LOGGER->LogAlgorithmSpecific(str_Perturb,str_Accept,"");

    UpdateStatistics(); // O(n)
  }

  PrintStatistics();

  return 0;
}

void Perturb4Opt(int Strength) {
  // Purpose    : Perturb the current TSP tour using 4-Opt double bridge move
  //              This 4-Opt move can't be reversed easily by 2-Opt/3-Opt.
  // Complexity : O(n)

  int pos1,pos2,pos3,newSolution[MAX_PROBLEM_SIZE];

  Prev_OV = Cur_OV;
  for (i=0; i<n; i++) { // backup...
    tempSolution[i] = solution[i];
    previousSolutionMap[i] = solutionMap[i];
  }

  for (l=0; l<Strength; l++) {
    // 4-Opt double bridge move... pick 3 split points randomly...
    pos1 = 1 + rand() % (n / 4);
    pos2 = pos1 + 1 + rand() % (n / 4);
    pos3 = pos2 + 1 + rand() % (n / 4);

    // Perturb from current solution, in ILS, current solution is the LO...
    for (i=j=0 ; i<pos1; i++,j++)
      newSolution[j] = solution[i]; // Part A
    for (i=pos3; i<n   ; i++,j++) 
      newSolution[j] = solution[i]; // Part D
    for (i=pos2; i<pos3; i++,j++) 
      newSolution[j] = solution[i]; // Part C
    for (i=pos1; i<pos2; i++,j++) 
      newSolution[j] = solution[i]; // Part B

    for (i=0; i<n; i++) {
      solution[i] = newSolution[i]; // put back
      solutionMap[solution[i]] = i;
    }
  }

  assert(IsValidPermutation());
  Cur_OV = Evaluate();
}

void TwoOpt(int t1,int t2,int t3,int t4) {
  // Purpose    : Given the two edges t1-t2, t3-t4, change the permutation so that the edges become t1-t3, t2-t4
  // Complexity : O(n)

  int i = solutionMap[t2];
  int j = solutionMap[t3];

  if (i>j) { // the other way around..
    i = solutionMap[t4];
    j = solutionMap[t1];
  }

  while (i<j) { // flip all along the route
    solution[i] ^= solution[j];
    solution[j] ^= solution[i];
    solution[i] ^= solution[j];
    solutionMap[solution[i]] ^= solutionMap[solution[j]];
    solutionMap[solution[j]] ^= solutionMap[solution[i]];
    solutionMap[solution[i]] ^= solutionMap[solution[j]];
    i++;
    j--;
  }
}

void LocalSearch() {
  // Purpose    : Given the current solution, do steepest descent to TSP local optima by considering potential improving moves only.
  // Complexity : sub-quadratic, O(k*n)

  int x,t1,t2,t3,t4;
  int d_t1_t2,d_t2_t4;
  vector<int> bestT1,bestT2,bestT3,bestT4;

  while (1) { // O(K*n), just few steps since the perturbation is not that strong...
    Best_Neighbor_OV = 2147483647;

    for (t1=0; t1<n; t1++) { // O(k*n, sub-quadratic 2-Opt neighborhood)
      t2 = solution[(solutionMap[t1]+1)%n]; // O(1)
      d_t1_t2 = distance_matrix[t1][t2];

      // scan candidate list of city t2
      for (x=0; x<n && x<CANDIDATE_LIST_SIZE; x++) { // O(k)... not all edges will be checked
        t3 = candidate_list[t2][x];
        if (t3 == t1 || t3 == t2)
          continue; // not a good candidate

        t4 = solution[(solutionMap[t3]+1)%n]; // O(1)
        if (t4 == t1 && t4 == t2)
          continue; // not a good candidate too

        d_t2_t4 = distance_matrix[t2][t4];
        if (d_t2_t4 >= d_t1_t2) 
          break; // sub-quadratic achieved here..., this is time to stop...

        Neighbor_OV = EvaluateNeighbor(t1,t2,t3,t4); // O(1)

        if (Neighbor_OV < Best_Neighbor_OV) { // O(1)
          Best_Neighbor_OV = Neighbor_OV;
          bestT1.clear();
          bestT2.clear();
          bestT3.clear();
          bestT4.clear();
        }

        if (Neighbor_OV == Best_Neighbor_OV) {
          bestT1.push_back(t1);
          bestT2.push_back(t2);
          bestT3.push_back(t3);
          bestT4.push_back(t4);
        }
      }
    }

    if (Best_Neighbor_OV < Cur_OV) { // O(n), there exist improving move
      int randPos = (int) (rand()%bestT1.size()); // throw fair coin to break ties
      TwoOpt(bestT1[randPos],bestT2[randPos],bestT3[randPos],bestT4[randPos]); // O(n), go there...
      Cur_OV = Best_Neighbor_OV;
    }
    else // local optima, quit...
      break;
  }
}

void FDD_Diversification() {
  int Temporary_Solution[20][MAX_PROBLEM_SIZE];
  double Temporary_OV[20];
  int Sorted_Temporary_Solution[20][MAX_PROBLEM_SIZE];
  double Sorted_Temporary_OV[20];
  int Sorted_Temporary_Distance[20];
  int maxDistance = -1, maxID = -1;
  int a,b;
  int counter = 0;

  // enlarge from n/4 to n/3
  while (maxDistance < n/3 && counter < 10) {
    counter++;
    for (a=0; a<10; a++) {
      // Set solution = bestSolution, we want to perturb from here...
      for (b=0; b<n; b++)
        solution[b] = bestSolution[b];

      // Perturb randomly and optimize it...
      if (n < 200) {
        Perturb4Opt(1);
        LocalSearch();
      }
      else {
        Perturb4Opt(1);
        LocalSearch();
        Perturb4Opt(1); // 2 times... :O, to make it further away...
        LocalSearch();
      }

      // Save current solution and the objective value
      for (b=0; b<n; b++)
        Temporary_Solution[a][b] = solution[b];
      Temporary_OV[a] = Cur_OV;
    }

    for (a=0; a<5; a++) { // sort it, using selection sort
      double Min = 32767*32767;
      int MinID = -1;

      for (b=0; b<10; b++) // selection sort here
        if (Temporary_OV[b] != -1 && Temporary_OV[b] < Min) { // not flagged and better
          Min = Temporary_OV[b];
          MinID = b;
        }

      for (b=0; b<n; b++)
        Sorted_Temporary_Solution[a][b] = Temporary_Solution[MinID][b];
      Sorted_Temporary_OV[a] = Temporary_OV[MinID]; 

      Temporary_OV[MinID] = -1; // flag as used...
    }

    for (a=0; a<5; a++) { // measure the distance of the top 5
      for (b=0; b<n; b++) // restore
        solution[b] = Sorted_Temporary_Solution[a][b];
      // compute
      Sorted_Temporary_Distance[a] = BondDistanceToBF();
      if (Sorted_Temporary_Distance[a] > maxDistance) {
        maxDistance = Sorted_Temporary_Distance[a];
        maxID = a;
      }
    }
  }

  for (a=0; a<n; a++) { // set the new solution into this
    solution[a] = Sorted_Temporary_Solution[maxID][a];
    solutionMap[Sorted_Temporary_Solution[maxID][a]] = a;
  }

  Cur_OV = Evaluate();
}

bool AcceptanceTest() {
  // Purpose    : Accept or reject a new local optima
  // Complexity : Omega(1) - O(n)

  if (BETTER_ACCEPTANCE_CRITERIA && (Cur_OV >= Prev_OV)) { // If BETTER_ACCEPTANCE_CRITERIA is used but currently new local optima is poorer or equal...
    if (NonImprovingMoves > NON_IMPROVING_MOVES_TOLERANCE * n) { // If it exceeds the tolerance limit...
      // strong perturb not total restart!
      printf("FDD Diversification\n");
      FDD_Diversification();
      LocalSearch();
      NonImprovingMoves = 0; // must reset this...
      return true; // just move here...
    }

    // If worse and distance is "too far", restore
    // Note that, even for equal OV (usually returning to itself), I reject that move...
    Cur_OV = Prev_OV;
    for (i=0; i<n; i++) { // O(n)
      solution[i] = tempSolution[i];
      solutionMap[i] = previousSolutionMap[i];
    } 
    return false; // if this is executed, I am actually rejecting this...
  }
  else if (!BETTER_ACCEPTANCE_CRITERIA) {
    // if random walk, always accept...
    return true;
  }

  return true;
}