#include <stdio.h>
#include <stdlib.h>
#include <string.h>

extern "C" {
   long random();
           };

class link;

class node 
{
  private:
   link ** links;
   int capacity;

  public:
   int nclusters;
   int state;
   int nlinks;
   double thcy;
   char * name;   
   void readName (FILE * in){
      char buf [1000];
      if (1 != fscanf (in, "%s", buf)){
         exit (-1);
      }
      name = strdup (buf);
   }      
   double toState (int s);     
   node (int ncl = 2){
      nclusters = ncl;
      capacity = 1000;
      links = new link * [capacity];
      nlinks = 0;
   }
   double getErr ();
   void push (link * l){
      nlinks ++;
      if (nlinks <= capacity){
         links [nlinks-1] = l;
      } else {
         capacity += 1000;
         link ** tmp = links;
         links = new link * [capacity];
         links [nlinks-1] = l;
         for (int i = 0; i < nlinks - 1; i++) links [i] = tmp [i];
         delete tmp;
      }
   }     
};

class link
{
  private:
   int capacity, nnodes;
   node ** nodes;
  public:
   int * cnt;
   int own, nclusters;
   long stmp;
   link () {
      nnodes = 0;
      cnt = NULL;
      capacity = 100;
      nodes = new node * [capacity];
   };
   void print (FILE * out){
      for (int k = 0; k < nnodes; k++) {
         fprintf (out, ";%s", nodes [k] ->name);
      }
      return;
   }  
   int sameState (){
      int bad=0;
      int st = nodes [0] ->state;
      for (int k = 1; k < nnodes; k++) {
         if (nodes [k] ->state != st){
            bad = 1;
            break;
         }
      }
      return !bad;
   }
   int allVal (int st){
      int bad = 0;
      for (int k = 0; k < nnodes; k++) {
         if (nodes [k] ->state != st){
            bad = 1;
            break;
         }
      }
      return !bad;
   }
   double change (int s0, int s1){
      if (own > nclusters || own < 0 || s0 < 0 || s0 > nclusters || s1 < 0 || s1 > nclusters)
         fprintf (stderr, "bad link: %d,%d,%d\n", own, s0, s1);
      
      double err = 0;
      if (cnt [own] < cnt [s1]+1 || own == s0){
         double tmp = getErr ();
         cnt [s0] --;
         cnt [s1] ++;
         own = 0; int ct = cnt [0];
         for (int j = 1; j < nclusters; j++){
            if (ct < cnt [j]){
               ct = cnt [j];
               own = j;
            }
         }
         err = getErr () - tmp;
      } else {
         cnt [s0] --;
         cnt [s1] ++;
      }
      return err;
   }
   double getErr (){
      double err = 0;
      for (int j = 0; j < nclusters; j++){
         if (j != own){
            err += cnt [j];
         } else {
            err -= cnt [j];
         }
      }
      return err;
   }
      
   void setNClusters (int ncl){
      nclusters = ncl;
      if (cnt) delete cnt;
      cnt = new int [ncl];
      for (int i = 0; i < ncl; i++) cnt [i] = 0;
   }
   void setOwn (){
      for (int j = 0; j < nnodes; j++)  {
         if (nodes [j] ->state >= nclusters  || nodes [j] ->state < 0)
            fprintf (stderr, "bad link: %d,%d\n", j, nodes [j] ->state);
         cnt [nodes [j] ->state] ++;
      }
   }
   void push (node * n){
      nnodes ++;
      if (nnodes <= capacity){
         nodes [nnodes-1] = n;
      } else {
         capacity += 100;
         node ** tmp = nodes;
         nodes = new node * [capacity];
         nodes [nnodes-1] = n;
         for (int i = 0; i < nnodes - 1; i++) nodes [i] = tmp [i];
         delete tmp;
      }
   }
   
};

double node ::toState (int s){
   double err = 0;
   for (int i = 0; i < nlinks; i++)
      err += links [i] ->change (state, s);
   if (s > nclusters || s < 0) 
      fprintf (stderr, "bad state %d\n", s);   
   state = s;
   return err;
}

double node ::getErr (){
   double err = 0;
   for (int i = 0; i < nlinks; i++)
      err += links [i] ->getErr ();
   return err;
}

double initializeNodes (double thcyFrom, double thcyTo,
           int nnodes,
           node * nodes)
{
   int selected [nnodes];
   double thcy=0;
   for (int i = 0; i < nnodes; i++) {
      selected[i] = 0;
      nodes [i] .state = 0;
   }
   while (thcy < thcyFrom){      
      int i = random () % nnodes;
      if (selected [i]) continue;
      if (thcy + nodes [i] .thcy > thcyTo) continue;
      nodes [i] .state = 1;
      thcy += nodes [i] .thcy;
   }
   return thcy;
}
double calcThcy (int nnodes,
           node * nodes)
{
   double thcy = 0;
   for (int i = 0; i < nnodes; i++)
      if (nodes [i] .state)
         thcy += nodes [i] .thcy;
   return thcy;
}

double err (int nlinks, link * links)
{
   double err = 0;
   for (int i = 0; i < nlinks; i++)
      err += links [i] .getErr ();
   return err;
}   

double err1 (int nnodes, node * nodes)
{
   double err = 0;
   for (int i = 0; i < nnodes; i++)
      err += nodes [i] .getErr ();
   return err;
}   


double change (int i, int j, node * nodes, double *thcy)
{
   if (nodes [i] .state == nodes [j] .state)
      return 0;
   
   int a = i, b = j;
   if (nodes [b] .state){
      a = j; b = i;
   }
   *thcy += nodes [b] .thcy - nodes [a] .thcy;
   
   int tmp = nodes [j] .state;
   double err = nodes [j] .toState (nodes [i] .state);
   err += nodes [i] .toState (tmp);
   if (tmp > 1 || tmp < 0) fprintf (stderr, "bad state %d,%d,%d\n", i, j, tmp);
   return err;
}

double change (int i, node * nodes, int state, double * thcy)
{
   if (nodes [i] .state == state)
      return 0;
   *thcy += state ? (nodes [i] .thcy):(-nodes [i] .thcy);
   if (state > 1 || state < 0) fprintf (stderr, "bad state %d,%d\n", i, state);
   return nodes [i] .toState (state);
}

void setLinks (int nclusters, int nlinks, link * links)
{
    // following needs rerun if the number of clusters changes 
   for (int i = 0; i < nlinks; i++){
       //set cluster size
      links [i] .setNClusters (nclusters);
       // fill in ownership counts
      links [i] .setOwn ();
   }

    //calculate link ownership
   for (int i = 0; i < nlinks; i++){
      links [i] .own = 0;
      int cnt = links [i] .cnt [0];
      for (int j = 1; j < nclusters; j++){
         if (cnt < links [i] .cnt [j]){
            cnt = links [i] .cnt [j];
            links [i] .own = j;
         }
      }
   }
}

double Min (int n, int * x)
{
   double res = x[0];
   for (int i = 1; i < n; i++) res = (res < x [i] ? res : x [i]);
   return res;
}

int setRange (int nclusters, int nnodes, int * clustersizes, int ranges)
{
   int res = int ((Min (nclusters, clustersizes) / nnodes * nclusters)*ranges);
   if (res >= ranges) res = ranges - 1;
   if (res < 0) res = 0;
   return res;
}   

int setRange (double thcy, double thcyFrom, double thcyTo, int ranges)
{
   int res = int((thcy-thcyFrom)/(thcyTo-thcyFrom)*ranges);
   if (res >= ranges) res = ranges - 1;
   if (res < 0) res = 0;
   return res;
}   

int main (int argc, char * argv [])
{
   char fname [250];
   sprintf (fname, "%s.dsc", argv[1]);
   FILE * in = fopen (fname, "r"), * in1;   
   int nnodes, nlinks;
   if (2 != fscanf (in, "%d %d", &nnodes, &nlinks))
      exit (-1);      
   fclose (in);
   
   int nclusters = 2;
   double thcyFrom = atof (argv[2]);
   double thcyTo = atof (argv[3]);

   node * nodes = new node [nnodes];
   sprintf (fname, "%s.cds", argv[1]);
   in = fopen (fname, "r");
   for (int i = 0; i < nnodes; i++){
      int j;
      double thcy;
      if (2 != fscanf (in, "%d %le", &j, &thcy))
         exit (-1); 
      nodes [j] .thcy = thcy;
      nodes [j] .readName (in);
   }
   fclose (in);
   double thcy = initializeNodes (thcyFrom, thcyTo, nnodes, nodes);
   
   link * links = new link [nlinks];
    // read links
   sprintf (fname, "%s.asc", argv[1]);
   in = fopen (fname, "r");
   int i = 0;
   while (!feof(in)){
      int n,  from;
      long stmp;
      if (2 != fscanf (in, "%d %ld", &n, &stmp)){
         break;
      }
      links [i] .stmp = stmp;
      for (int j = 0; j < n; j++){
         if (1 != fscanf (in, "%d", &from))
            break;
         if (from >= nnodes){
            fprintf (stderr, "number of nodes=%d is less than read=%d\n", nnodes, from);
            exit (-1);
         }   
         links [i] .push (nodes + from);
         nodes [from] .push (links + i);
      }
      i++;
   }
   fclose (in);
   if (i != nlinks){
      fprintf (stderr, "number of links=%d ne specified=%d\n", i, nlinks);
      exit (-1);
   }
   setLinks (nclusters, nlinks, links);

    // optimize
   int ranges = 10;
   double best [ranges], curr [ranges], touched [ranges];
   int state [ranges][nnodes];
   for (int j = 0; j < ranges; j++) touched [j] = 0;

   int range = setRange (thcy, thcyFrom, thcyTo, ranges);
   curr [range] = best [range] = err (nlinks, links);
   touched [range] ++;

   for (int j = 0; j < nnodes; j++)
      state [range][j] = nodes [j] .state; 
   
   long niter = (long) 1e7;
    if (argc > 4)
       niter = (long) atof (argv[4]);
   long nrefresh = 0;
   if (argc > 5)
      nrefresh = niter / ((long) atof (argv[5]));

   for (long k = 0; k < niter; k++){
      if (nrefresh && k%nrefresh == nrefresh-1){
         int doRest = 1;
         if (int (random () % 2)){
             // set up best solution from existing ranges
             range = int (random () % ranges);
             if (touched [range]){
                curr [range] = best [range] = err (nlinks, links);
                for (int j = 0; j < nnodes; j++){
                   nodes [j] .state = state [range][j];
                }
                thcy = calcThcy (nnodes, nodes);
                setLinks (nclusters, nlinks, links);
                doRest = 0;
             }
         }
         if (doRest){
             // initialize a random solution
            thcy = initializeNodes (thcyFrom, thcyTo, nnodes, nodes);
            range = setRange (thcy, thcyFrom, thcyTo, ranges);
            setLinks (nclusters, nlinks, links);
            if (!touched [range])
               for (int j = 0; j < nnodes; j++)
                  state [range][j] = nodes [j] .state;
            touched [range] ++;
      }
         continue;            
      }
      double ch = 0;
      int a = int (random () % nnodes);
      int b = int(random () % nnodes);
      if (a == b) continue;
      if (random () % 3 > 0){//probability of not changing N
         if (nodes [a] .state == nodes [b] .state) continue;
         
         ch = change (a, b, nodes, &thcy);
         if (thcy < thcyFrom || thcy > thcyTo || // if thcy is out of range or
             (ch > 0 && random () % 3 > 0)){     //probability of not decreaseing criteria
            change (b, a, nodes, &thcy);
            continue;
         }

         range = setRange (thcy, thcyFrom, thcyTo, ranges);
         if (!touched [range]){
            curr [range] = best [range] = err (nlinks, links);
            for (int j = 0; j < nnodes; j++){
               state [range][j] = nodes [j] .state;
            }
         }
         touched [range]++;
      }else{
          // do add/subtract based on size?
         int c = random () % nclusters;
         if (c != nodes [a] .state){
            int prevState = nodes [a] .state;
            ch = change (a,  nodes, c, &thcy);
            if (thcy < thcyFrom || thcy > thcyTo || // if thcy is out of range or
                (ch > 0 && random () % 2 > 0)){     //probability of not decreaseing criteria
               change (a, nodes, prevState, &thcy);               
               continue;
            }            
            range = setRange (thcy, thcyFrom, thcyTo, ranges);
            if (!touched [range]){
               curr [range] = best [range] = err (nlinks, links);
               for (int j = 0; j < nnodes; j++){
                  state [range][j] = nodes [j] .state;
               }
            }
            touched [range]++;
         }else
            continue;
      }      

      curr [range] += ch;
      if (best [range] > curr [range]) {
         best [range] = curr [range];
         for (int j = 0; j < nnodes; j++){
            state [range][j] = nodes [j] .state;
         }
      }
   }

    // result
   sprintf (fname, "%s.res", argv[1]);
   in = fopen (fname, "w");   
   sprintf (fname, "%s.freq", argv[1]);
   in1 = fopen (fname, "w");
   
   for (int i = 0; i < ranges; i++){
      if (touched [i]){
         int sum = 0;
         for (int j = 0; j < nnodes; j++) {
            fprintf (in, "%d %d %s\n", i, state [i][j], nodes [j] .name);
            change (j, nodes, state [i][j], &thcy);
            sum += state [i][j] == 1;
         }
         int sumBad = 0, sumGood = 0, sumTouch = 0;         
         for (int j = 0; j < nlinks; j++) {
            int bad = ! links [j] .sameState();
            int inside = links [j] .allVal (1);
             //fprintf (stderr, "%d;%d;%d", i, j, bad);            
             //links [j] .print (stderr);
             //fprintf (stderr, "\n");            
            sumBad += bad ? 1 : 0;
            sumGood += bad ? 0 : 1;
            sumTouch += (bad||inside) ? 1 : 0;
         }
          // cluster, nodes in, mrs crossing, remaining mrs, mrs in clstr 0?, mrs in clstr 1?, niter
          //fprintf (in1, "%d %d %d %d %d %d %.1lf\n", i, sum, sumBad, sumGood, sumClstr [0], sumClstr [1], touched [i]);
          // cluster, niter, nodes in, ..., mrs crossing, remaining mrs
         fprintf (in1, "%d %.1lf %d %.1lf %d %d %.4lf\n", i, touched [i], sum, err (nlinks, links), sumBad, sumTouch, (sumBad+0.0)/(sumTouch));
      }
   }
   fclose (in);
   fclose (in1);
}