/******************************************************************
 *                                                                *
 *         Nest version 2 - 3D Wasp Nest Simulator                *
 *            Copyright (C) 1997 Sylvain GUERIN                   *
 *        LIASC, ENST Bretagne & Santa Fe Institute               *
 *                                                                *
 ******************************************************************
 
 E-mail: Sylvain.Guerin@enst-bretagne.fr
 or: Sylvain Guerin
 13 rue des Monts Lorgeaux, 51460 L'EPINE, FRANCE
 
 ******************************************************************
 
 This program is free software; you can redistribute it and/or
 modify it under the terms of the GNU General Public License
 as published by the Free Software Foundation; either version 2
 of the License, or (at your option) any later version.
 
 This program is distributed in the hope that it will be useful,
 but WITHOUT ANY WARRANTY; without even the implied warranty of
 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 GNU General Public License for more details.
 
 You should have received a copy of the GNU General Public License
 along with this program; if not, write to the Free Software
 Foundation, Inc., 59 Temple Place - Suite 330, 
 Boston, MA  02111-1307, USA.
 
 ******************************************************************
     
      Name          : graph.c
      Component     : fitness evaluation
      Fonctionality : contains methods for
                      class Graph and GraphNode
  
 ******************************************************************/

#include "graph.h"
#include "classify.h"
#include <math.h>

/****************************** class GraphNode ******************************/

GraphNode::GraphNode ()
{
  RandomGenerator randomGenerator (PSEUDO, 0, 7999);
  int i;
  
  ruleIndex = -1;
  nodeIndex = -1;
  for (i=0; i<26; i++) {
    neighbour[i] = -1;
  }
  for (i=0; i<26; i++) {
    inOut[i] = Input;
  }

  px = (double)randomGenerator.getNumber ();
  //px = 4000.0;
  py = (double)randomGenerator.getNumber ();

  numberOfInputEdges = 0;
  numberOfOutputEdges = 0;

  sequencesNumber = 0;
  for (i=0; i<defaultMaxCells; i++) {
    sequencesLength[i] = 0;
  }
  
  belongsToPattern = 0;
  isPatternInitialized = 0;
}

GraphNode::~GraphNode ()
{
  if (isPatternInitialized) {
    delete pattern;
    isPatternInitialized=0;
  }
}

void GraphNode::initGraph (Graph* aGraph)
{
  graph = aGraph;
}

void GraphNode::initSequences ()
{
  int i;

  sequencesNumber = 0;
  for (i=0; i<defaultMaxCells; i++) {
    sequencesLength[i] = 0;
  }
}

void GraphNode::searchSequences (GraphNode* firstNode, 
				 int searchedRuleId, 
				 int length)
{
  int i;
  int firstCall;
  int l;

  /* printf ("Search sequences for node %d - length: %d - currentNode: %d\n",
	  firstNode->getNodeId (), length, getNodeId ()); */


  /* test si debut */
  
  if (length == 0) {
    // printf ("On initialise\n");
    for (i=0; i<graph->getNodesNumber (); i++) {
      (graph->nodeAtIndex (i))->foundLength = -1;
      (graph->nodeAtIndex (i))->exploredAtLength = -1;
      }
    firstCall = 1;
  }
  else {
    firstCall = 0;
  }
  
  if ((exploredAtLength == -1) || 
      ((exploredAtLength != -1) && (length < exploredAtLength))) {
    
    exploredAtLength = length;
    
    if (!firstCall) {
      // printf ("On regarde un peu...\n");
      firstCall = 0;
      if (searchedRuleId == getRuleId ()) {
	// printf ("On a trouve\n");
	/* on a trouve */
	if (foundLength == -1) {
	  foundLength = length;
	}
	else if (length < foundLength) {
	  foundLength = length;
	}
	return;
      }
    }
    
    /* je crois que j'ai trouve: il faut mettre un marqueur pour
       dire qu'on est deja passe ici */
    
    if (foundLength != -1) {
      // printf ("On s'en retourne par la...\n");
      return;
    }  /* verifier ce truc, pas vraiment sur de cet algo */
    
    /* on continue a chercher */
    
    for (i=0; i<outputEdgesNumber (); i++) {
      /* printf ("On relance pour le noeud %d\n", 
	      (outputNeighbour(i))->getNodeId ()); */
      if (outputNeighbour(i) == NULL) {
	printf ("Probleme au niveau du noeud %d (pour l'index %d)\n",
		getNodeId (), i);
	return;
      }
      (outputNeighbour(i))->searchSequences (firstNode, 
					     searchedRuleId, 
					     length+1);
    }
    
    if (firstCall == 1) {
      // printf ("C'est fini\n");
      /* on met a jour toutes les longueurs trouvees */
      for (i=0; i<graph->getNodesNumber (); i++) {
	l = (graph->nodeAtIndex (i))->foundLength;
	if (l != -1) {
	  sequencesLength[sequencesNumber] = l;
	  sequencesNumber++;
	  /* printf ("Node %d (rule %d): adding %d\n", 
	     getNodeId (), getRuleId (), l); */
	}
      }
    }
    
    // printf ("C'est normal qu'on passe par la ?????\n");
    
    return;
    
  }
  
  else { 
    // printf ("Pas la peine, deja explore...\n");
    return;
  }

}

void GraphNode::debugSequences ()
{
  int i;

  printf ("Node %d (rule %d) - %d sequ. : ",
	  getNodeId (), getRuleId (), sequencesNumber);
  for (i=0; i<sequencesNumber; i++) {
    printf ("[%d]", sequencesLength[i]);
  }
  printf ("\n");
}

void GraphNode::freePattern ()
{
  if (isPatternInitialized) {
    delete pattern;
    isPatternInitialized=0; 
  }
}

void GraphNode::initPattern (int aMaxSize)
{
  int i;
  
  if (isPatternInitialized) {
    delete pattern;
    isPatternInitialized=0;
  }
  pattern = new Pattern (aMaxSize, getRuleId (), graph->getSimuType ());
  isPatternInitialized = 1;
  sequencesNumber = 0;
  for (i=0; i<defaultMaxCells; i++) {
    sequencesLength[i] = 0;
  }
  for (i=0; i<graph->getNodesNumber (); i++) {
    (graph->nodeAtIndex (i))->foundLength = -1;
    (graph->nodeAtIndex (i))->exploredAtLength = -1;
    (graph->nodeAtIndex (i))->belongsToPattern = 0;
  } 
}
 
int GraphNode::searchPattern (GraphNode* firstNode, int length)
     /* return 1 if something was found, else 0 */
{
  int i;
  int firstCall;
/* USELESS  int l; */
  int searchedRuleId;
  
  searchedRuleId = firstNode->getRuleId ();
  
  firstCall = 0;
  if (length == 0) {
    firstCall = 1;
  }
  
  if ((exploredAtLength == -1) || 
      ((exploredAtLength != -1) && (length < exploredAtLength))) {
    
    exploredAtLength = length;
    
    if (!firstCall) {
      if (searchedRuleId == getRuleId ()) {
	/* on a trouve */
	if (foundLength == -1) {
	  foundLength = length;
	}
	else if (length < foundLength) {
	  foundLength = length;
	}
	else {
	  return 0;
	}
	if (length <= maxLength) {
	  belongsToPattern = 1;
	  return 1;
	}
	else {
	  return 0;
	} 
      }
    }
    
    /* on continue a chercher */
    
    for (i=0; i<outputEdgesNumber (); i++) {
      if (outputNeighbour(i) == NULL) {
	printf ("Probleme au niveau du noeud %d (pour l'index %d)\n",
		getNodeId (), i);
	return 0;
      }
      if ((outputNeighbour(i))->searchPattern (firstNode,length+1) == 1) {
	belongsToPattern = 1;
      }
    }
    
    if (firstCall == 1) {
      if (belongsToPattern) {
	/* C'est fini, on met a jour le pattern trouve */
	for (i=0; i<graph->getNodesNumber (); i++) {
	  if ((graph->nodeAtIndex (i))->belongsToPattern) {
	    firstNode->pattern->set ((graph->nodeAtIndex (i))->x,
				     (graph->nodeAtIndex (i))->y,
				     (graph->nodeAtIndex (i))->z,
				     (graph->nodeAtIndex (i))->cellType,
				     (graph->nodeAtIndex (i))->getRuleId (),
				     (graph->nodeAtIndex (i))->getNodeId ());
	  }
	}
      }
      else {
	firstNode->freePattern ();
      }
    }
    
    if (belongsToPattern) {
      return 1;
    }
    else {
      return 0;
    }
    
  }
  
  else { 
    if (length <= maxLength) {
      belongsToPattern = 1;
      return 1;
    }
    else {
      return 0;
    } 
  }
  
}

void GraphNode::set (int aRuleIndex, int aNodeIndex)
{
  ruleIndex = aRuleIndex;
  nodeIndex = aNodeIndex;
}

int GraphNode::inputEdgesNumber ()
{
  return numberOfInputEdges;
}

int GraphNode::outputEdgesNumber ()
{
  return numberOfOutputEdges;
}

int GraphNode::indexOfNeighbour (int aPosCode)
/* return -1 when not available */
{
  if (aPosCode < 26) {
    return neighbour[aPosCode];
  }
  else {
    printf ("ERROR: invalid posCode\n");
    return -1;
  }
}

InOutType GraphNode::directionOfNeighbour (int aPosCode)
{
  if (aPosCode < 26) {
    return inOut[aPosCode];
  }
  else {
    printf ("ERROR: invalid posCode\n");
    return Input; /* par defaut, on aurait pu choisir n'importe quoi */
  }
}

int GraphNode::indexOfInputNeighbour (int anIndex)
/* return -1 when not available */
{
  int i,j;
  
  i=0;
  j=-1;
  while ((j<anIndex) && (i<26)) {
    if ((neighbour[i] != -1) && (inOut[i] == Input)) {
      j++;
    }
    i++;
  }
  
  if (i<27) {
    return neighbour[i-1];
  }
  else {
    printf ("ERROR in GraphNode::indexOfInputNeighbour\n");
    return -1;
  }
}

int GraphNode::indexOfOutputNeighbour (int anIndex)
/* return -1 when not available */
{
  int i,j;
  
  i=0;
  j=-1;
  while ((j<anIndex) && (i<26)) {
    if ((neighbour[i] != -1) && (inOut[i] == Output)) {
      j++;
    }
    i++;
  }
  
  // printf ("i=%d\n", i);

  if (i<27) {
    return neighbour[i-1];
  }
  else {
    printf ("ERROR in GraphNode::indexOfOutputNeighbour\n");
    return -1;
  }
}

GraphNode* GraphNode::inputNeighbour (int anIndex)
{
  if (indexOfInputNeighbour (anIndex) != -1) {
    return graph->nodeAtIndex (indexOfInputNeighbour (anIndex));
  }
  else {
    printf ("ERROR: node undefined\n");
    return NULL;
  }
}

GraphNode* GraphNode::outputNeighbour (int anIndex)
{
  if (indexOfOutputNeighbour (anIndex) != -1) {
    return graph->nodeAtIndex (indexOfOutputNeighbour (anIndex));
  }
  else {
    printf ("ERROR: node undefined\n");
    debug ();
    return NULL;
  }
}

void GraphNode::addInputEdge (int posCode, int nodeIndex)
{
  if (posCode<26) {
    if (neighbour[posCode] != -1) {
      printf ("ERROR: incoherent data in edge building\n");
    }
    neighbour[posCode] = nodeIndex;
    inOut[posCode] = Input;
    numberOfInputEdges++;
  }
  else {
    printf ("ERROR: invalid position code (posCode = %d)\n", posCode);
  }
}

void GraphNode::addOutputEdge (int posCode, int nodeIndex)
{
  if (posCode<26) {
    if (neighbour[posCode] != -1) {
      printf ("ERROR: incoherent data in edge building\n");
    }
    neighbour[posCode] = nodeIndex;
    inOut[posCode] = Output;
    numberOfOutputEdges++;
  }
  else {
    printf ("ERROR: invalid position code\n");
  }
}

double GraphNode::getPx ()
{
  return px;
}

double GraphNode::getPy ()
{
  return py;
}

void GraphNode::setPx (double newPx)
{
  px = newPx;
}

void GraphNode::setPy (double newPy)
{
  py = newPy;
}  

int GraphNode::getNodeId ()
{
  return nodeIndex;
}

int GraphNode::getRuleId ()
{
  return ruleIndex;
}

void GraphNode::debug()
{
  int i;
  
  printf ("NodeId: %d RuleId: %d Inputs: %d Outputs: %d [%f/%f]\n",
	  nodeIndex, ruleIndex, 
	  numberOfInputEdges, numberOfOutputEdges, 
	  px, py);

  for (i=0; i<26; i++) {
    printf ("[%.3d]", i);
  }
  printf ("\n");
  for (i=0; i<26; i++) {
    printf ("[%.3d]", neighbour[i]);
  }
  printf ("\n");
  for (i=0; i<26; i++) {
    printf ("[%.3d]", inOut[i]);
  }
  printf ("\n");
}

/****************************** class Graph ******************************/

Graph::Graph (SimulationType aSimuType, int aMaxCells, int anArchSize)
{
  int i, j;

  simuType = aSimuType;
  maxCells = aMaxCells;
  archSize = anArchSize;
  nodesNumber = 0;
  edgesNumber = 0;
  node = new GraphNode*[maxCells];
  /* for (i=0; i<maxCells; i++) {
     node[i]->initGraph (this);
     } */

  pCode[0][0][0] = 0;  pCode[1][0][0] = 1; pCode[2][0][0] = 2; 
  pCode[0][1][0] = 3;  pCode[1][1][0] = 4; pCode[2][1][0] = 5; 
  pCode[0][2][0] = 6;  pCode[1][2][0] = 7; pCode[2][2][0] = 8; 
  
  pCode[0][0][1] = 9;  pCode[1][0][1] = 10; pCode[2][0][1] = 11; 
  pCode[0][1][1] = 12; pCode[1][1][1] = -1; pCode[2][1][1] = 13; 
  pCode[0][2][1] = 14; pCode[1][2][1] = 15; pCode[2][2][1] = 16; 
  
  pCode[0][0][2] = 17; pCode[1][0][2] = 18; pCode[2][0][2] = 19; 
  pCode[0][1][2] = 20; pCode[1][1][2] = 21; pCode[2][1][2] = 22; 
  pCode[0][2][2] = 23; pCode[1][2][2] = 24; pCode[2][2][2] = 25;
  
  /* correspond a x%2 == 1 */
  
  p1Code[0][0][0] = 2;  p1Code[1][0][0] = 1; p1Code[2][0][0] = 6; 
  p1Code[0][1][0] = 3;  p1Code[1][1][0] = 0; p1Code[2][1][0] = 5; 
  p1Code[0][2][0] = -1; p1Code[1][2][0] = 4; p1Code[2][2][0] = -1; 
  
  p1Code[0][0][1] = 8;  p1Code[1][0][1] = 7;  p1Code[2][0][1] = 9; 
  p1Code[0][1][1] = 10; p1Code[1][1][1] = -1; p1Code[2][1][1] = 11; 
  p1Code[0][2][1] = -1; p1Code[1][2][1] = 12; p1Code[2][2][1] = -1; 
  
  p1Code[0][0][2] = 14; p1Code[1][0][2] = 15; p1Code[2][0][2] = 16; 
  p1Code[0][1][2] = 13; p1Code[1][1][2] = 19; p1Code[2][1][2] = 17; 
  p1Code[0][2][2] = -1; p1Code[1][2][2] = 18; p1Code[2][2][2] = -1;
  
  /* correspond a x%2 == 0 */
  
  p2Code[0][0][0] = -1; p2Code[1][0][0] = 1; p2Code[2][0][0] = -1; 
  p2Code[0][1][0] = 2;  p2Code[1][1][0] = 0; p2Code[2][1][0] = 6; 
  p2Code[0][2][0] = 3;  p2Code[1][2][0] = 4; p2Code[2][2][0] = 5; 
  
  p2Code[0][0][1] = -1; p2Code[1][0][1] = 7;  p2Code[2][0][1] = -1; 
  p2Code[0][1][1] = 8;  p2Code[1][1][1] = -1; p2Code[2][1][1] = 9; 
  p2Code[0][2][1] = 10; p2Code[1][2][1] = 12; p2Code[2][2][1] = 11; 
  
  p2Code[0][0][2] = -1; p2Code[1][0][2] = 15; p2Code[2][0][2] = -1; 
  p2Code[0][1][2] = 14; p2Code[1][1][2] = 19; p2Code[2][1][2] = 16; 
  p2Code[0][2][2] = 13; p2Code[1][2][2] = 18; p2Code[2][2][2] = 17;

  /* init de la matrice de co-occurences */
  
  for (i=0; i<maxRules; i++) {
    for (j=0; j<maxRules; j++)
      cooccurence[i][j] = 0;
    used[i] = 0;
  }
  maxRulesIndex = 0;
 
  /* init de kohonen */

  resetKohonen ();

  /* init des patterns */
  
  patternsNumber = 0;

  /* init des fitness */
  
  complexityFitness = 0.0;
  coherenceFitness = 0.0;
  isComplexityComputed = 0;
  isCoherenceComputed = 0;
  arePatternsComputed = 0;
  isPatternMatchingDone = 0;
}

Graph::~Graph ()
{
  int i;

  freeAllPatterns ();

  for (i=0; i<nodesNumber; i++) {
    delete node[i];
  }

  delete[] node;

}

GraphNode* Graph::nodeAtIndex (int aNodeId)
{
  return node[aNodeId];
}

void Graph::addNode (int aRuleIndex, typeOfCell aCellType, 
		     int ax, int ay, int az)
{
  if (nodesNumber >= maxCells) {
    printf ("ERROR: maxCells reached for graph building\n");
    return;
  }
  
  node[nodesNumber] = new GraphNode;
  node[nodesNumber]->initGraph (this);
  
  node[nodesNumber]->set (aRuleIndex, nodesNumber);
  node[nodesNumber]->cellType = aCellType;
  node[nodesNumber]->x = ax;
  node[nodesNumber]->y = ay;
  node[nodesNumber]->z = az;
  
  nodesNumber++;
  
  if (aRuleIndex >= maxRulesIndex) {
    maxRulesIndex = aRuleIndex;
  }

  if (used[aRuleIndex] == 0) {
    used[aRuleIndex] = 1;
  }

}

void Graph::addEdge (int startNodeId, int endNodeId, int posCode)
{
  int ruleId1, ruleId2;

  /* on met les connections bien comme il faut */

  if (simuType == cubic) {
    node[startNodeId]->addOutputEdge (25-posCode, endNodeId);
  }
  else if (simuType == hexa) {
    node[startNodeId]->addOutputEdge (19-posCode, endNodeId);
  }
  
  node[endNodeId]->addInputEdge (posCode, startNodeId);

  edgesNumber++;
  
  /* on met maintenant a jour la matrice de co-occurence */
  
  ruleId1 = node[startNodeId]->getRuleId ();
  ruleId2 = node[endNodeId]->getRuleId ();
  
  if ((ruleId1 >= 0) && (ruleId1 <= maxRules)
      && (ruleId2 >= 0) && (ruleId2 <= maxRules)) {

    (cooccurence[ruleId1][ruleId2])++;
  }
  
}

int Graph::posCodeCubic (int x, int y, int z)
{
    if (simuType == cubic)
	{
	    return pCode[x+1][y+1][z+1];
	}

/*    printf ("ERROR: incoherent architecture type\n"); */
    exit(1);
}

int Graph::posCodeHexa (int x, int y, int z, int decX)
{
    if (simuType == hexa)
	{
	    if (decX == 1) /* xxx */
		return p1Code[x+1][y+1][z+1];
	    else 
		return p2Code[x+1][y+1][z+1];
	}

/*    printf ("ERROR: incoherent architecture type\n"); */
    exit(1);
    
}

int Graph::getNodesNumber ()
{
  return nodesNumber;
}

int Graph::getEdgesNumber ()
{
  return edgesNumber;
}

SimulationType Graph::getSimuType ()
{
  return simuType;
}

void Graph::classifyRules (RuleArray* aRuleArray)
{
  int i, j;
  int nb;
  Rule* currentRule;
  int usedRuleNumber;
  RuleGraph* ruleGraph;
  int weight;

  nb = aRuleArray->getElementsNb ();
  usedRuleNumber = 0;

  for (i=0; i<nb; i++) {
    currentRule = aRuleArray->get (i);
    if (!used[i]) {
      currentRule->setRuleClass (unusedClass);
    }
    else {
      usedRuleNumber++;
      if (usedRuleNumber % 2 ==0) {
	currentRule->setRuleClass (classA);
      }
      else {
	currentRule->setRuleClass (classB);
      }
    }
  }

  /* on instancie le graphe */
  
  ruleGraph = new RuleGraph (usedRuleNumber+1);

  /* on cree les noeuds */

  for (i=0; i<nb; i++) {
    currentRule = aRuleArray->get (i);
    if (used[i]) {
      if (currentRule->getRuleClass () == classA) {
	ruleGraph->addNode (i, clusterA);
      }
      else {
	ruleGraph->addNode (i, clusterB);
      }
    }
  }

  /* on y colle des aretes */

  for (i=0; i<=maxRulesIndex; i++) {
    if (used[i]) {
      for (j=i+1; j<=maxRulesIndex; j++) {
	if (used[j]) {
	  weight = cooccurence[i][j]+cooccurence[j][i];
	  if (weight != 0) {
	    ruleGraph->addEdge (i, j, weight);
	  }
	}
      }
    }
  }
  
  /* algo de fidducia, maintenant */
  
  ruleGraph->partitionning ((int)((float)usedRuleNumber/3));

  /* on recupere les classes */
  
  for (i=0; i<nb; i++) {
    currentRule = aRuleArray->get (i);
    if (used[i]) {
      if ((ruleGraph->nodeWhichIdIs (i))->getCluster () == clusterA) {
	currentRule->setRuleClass (classA);
      }
      else {
	currentRule->setRuleClass (classB);
      }
    }
  }

  delete ruleGraph;
  
}

void Graph::displayCooccurenceMatrix (RuleArray* aRuleArray)
{
  Rule* currentRule1;
  Rule* currentRule2;
  int i, j;
  int c1, c2;
  
  printf ("Co-occurences matrix (after classification)\n");
  printf ("    ");
  
  printf ("|");
  for (c1=0; c1<2; c1++) {
    for (i=0; i<=maxRulesIndex; i++) {
      currentRule1 = aRuleArray->get (i);
      if ((int)currentRule1->getRuleClass () == c1) {
	printf ("| %.2d ", i+1);
      }
    }
    printf ("|");
  }
  printf ("|\n");
  
  for (i=0; i<=maxRulesIndex; i++) {
    if (aRuleArray->get(i)->getRuleClass () != 2) {
    printf ("=====");
    }
  }
  printf ("========\n");
  
  for (c1=0; c1<2; c1++) {
    for (i=0; i<=maxRulesIndex; i++) {
      currentRule1 = aRuleArray->get (i);
      if ((int)currentRule1->getRuleClass () == c1) {
	printf (" %.2d |", i+1);
	for (c2=0; c2<2; c2++) {
	  for (j=0; j<=maxRulesIndex; j++) {
	    currentRule2 = aRuleArray->get (j);
	    if ((int)currentRule2->getRuleClass () == c2) {
	      if (cooccurence[i][j] == 0) {
		printf ("|    ");
	      }
	      else if (cooccurence[i][j] < 10) {
		printf ("|  %d ", cooccurence[i][j]);
	      }
	      else {
		printf ("| %.2d ", cooccurence[i][j]);
	      }
	    }
	  }
	  printf ("|");
	}
	printf ("|\n");
      }
    }
    for (i=0; i<=maxRulesIndex; i++) {
      if (aRuleArray->get(i)->getRuleClass () != 2) {
	printf ("=====");
      }
    }
    printf ("========\n");
  }
  printf ("\n");
}

void Graph::displayCooccurenceMatrix ()
{
  int i, j;
  
  printf ("Co-occurence matrix (no classification)\n");
  printf ("    ");
  for (i=0; i<=maxRulesIndex; i++) {
    if (used[i])
      printf ("  %.2d ", i+1);
  }
  printf ("\n");
  
  for (i=0; i<=maxRulesIndex; i++) {
    if (used[i]) {
      printf (" %.2d ", i+1);
      for (j=0; j<=maxRulesIndex; j++) {
	if (used[j]) {
	  if (cooccurence[i][j] == 0) {
	    printf ("|    ");
	  }
	  else if (cooccurence[i][j] < 10) {
	    printf ("|  %d ", cooccurence[i][j]);
	  }
	  else {
	    printf ("| %.2d ", cooccurence[i][j]);
	  }
	}
      }
      printf ("|\n");
    }
  }
  printf ("\n");
}

void Graph::debug()
{
  int i;

  printf ("Debugging graph object\n");

  printf ("Nodes number: %d\n", nodesNumber);
  printf ("MaxCells: %d\n", maxCells);
  
  for (i=0; i<nodesNumber; i++) {
    (nodeAtIndex (i))->debug ();
  }
}

/* algo de kohonen (visualisation du graphe) */

int Graph::getCurrentIt ()
{
  return currentIteration;
}

int Graph::getMaxIt ()
{
  return maxIterations;
}

double Graph::getEpsilon1 ()
{
  return epsilon1;
}

double Graph::getEpsilon2 ()
{
  return epsilon2;
}
  
void Graph::setMaxIt (int aInt)
{
  maxIterations = aInt;
}

void Graph::setCurrentIt (int aInt)
{
  currentIteration = aInt;
}

void Graph::updateEpsilon1 ()
{
  epsilon1 = k1*(1-currentIteration/maxIterations);
}

void Graph::updateEpsilon2 ()
{
  epsilon2 = k2*(1-currentIteration/maxIterations);
}

void Graph::resetKohonen ()
{
  currentIteration = 0;
  maxIterations = defaultMaxIterations;
  k1 = defaultK1;
  k2 = defaultK2;
  D = defaultD;
  updateEpsilon1 ();
  updateEpsilon2 ();
}

void Graph::computeKohonenIteration ()
{
  RandomGenerator randomGenerator1 (PSEUDO, 0, 7999); /* newX */
  RandomGenerator randomGenerator2 (PSEUDO, 0, 7999); /* newY */
  double newX, newY;
  double x, y;
  int modifiedNodeId;
  double bestSqDistance;
  double SqDistance;
  int i;
  GraphNode* current;
  double d;
  
  /* on remet a jour epsilon1 et epsilon2 */
  
  currentIteration++;
  updateEpsilon1 ();
  updateEpsilon2 ();
  
  /* premiere etape: on tire un point au hasard */
  
  newX = (double)randomGenerator1.getNumber ();
  newY = (double)randomGenerator2.getNumber ();
  
  /* deuxieme etape: on cherche le noeud concerne */
  
  bestSqDistance = 1e20; /* c'est beaucoup, hein ??? */
  modifiedNodeId = -1; 
  
  for (i=0; i<nodesNumber; i++) {
    x = node[i]->getPx ();
    y = node[i]->getPy ();
    SqDistance = (newX-x)*(newX-x) + (newY-y)*(newY-y);
    if (SqDistance < bestSqDistance) { 
      bestSqDistance = SqDistance;
      modifiedNodeId = i;
    }
  }

  /* troisieme etape: on rapproche le noeud */

  x = node[modifiedNodeId]->getPx ();
  y = node[modifiedNodeId]->getPy ();
  node[modifiedNodeId]->setPx (x+epsilon1*(newX-x));
  node[modifiedNodeId]->setPy (y+epsilon1*(newY-y));
  
  newX = x+epsilon1*(newX-x);
  newY = y+epsilon1*(newY-y);
  
  /* quatrieme etape: on rapproche les voisins (input) */
  
  for (i=0; i<node[modifiedNodeId]->inputEdgesNumber (); i++) {
    /* current: noeud voisin courant */
    current = node[node[modifiedNodeId]->indexOfInputNeighbour (i)];
    x = current->getPx ();
    y = current->getPy ();
    d = sqrt ((x-newX)*(x-newX)+(y-newY)*(y-newY));
    if (3*currentIteration > maxIterations) {
      epsilon2 = k2*exp (-d/D);
    }
    current->setPx (x+epsilon2*(newX-x));
    current->setPy (y+epsilon2*(newY-y));
  }
  
  /* cinquieme etape: on rapproche les voisins (output) */
  
  for (i=0; i<node[modifiedNodeId]->outputEdgesNumber (); i++) {
    /* current: noeud voisin courant */
    current = node[node[modifiedNodeId]->indexOfOutputNeighbour (i)];
    x = current->getPx ();
    y = current->getPy ();
    d = sqrt ((x-newX)*(x-newX)+(y-newY)*(y-newY));
    if (3*currentIteration > maxIterations) {
     epsilon2 = k2*exp (-d/D);
    }
    current->setPx (x+epsilon2*(newX-x));
    current->setPy (y+epsilon2*(newY-y)); 
  }
  
}

void Graph::computeKohonen ()
{
  int i;
  
  resetKohonen ();
  for (i=0; i<defaultMaxIterations; i++) {
    computeKohonenIteration ();
  }
}


/* calcul des sequences de regles */

void Graph::initRulesSequences ()
{
  int i;

 // printf ("initRulesSequences ()\n");

  for (i=0; i<nodesNumber; i++) {
    node[i]->initSequences ();
  }

}

void Graph::computeRulesSequences ()
{
  int i,j,k;
/* USELESS  int currentRuleId; */
  Sequences sequences[maxRules][maxRules];

  initRulesSequences ();
  
  for (i=1; i<nodesNumber; i++) {
    // printf ("Computing node %d\n", i);
    for (j=0; j<maxRulesIndex; j++) {
      node[i]->initSequences ();
      node[i]->searchSequences (node[i], j, 0);
      for (k=0; k<node[i]->sequencesNumber; k++) {
	sequences[node[i]->getRuleId()][j].addData (node[i]->sequencesLength[k]);
      }
    }
  }
  
  /* for (currentRuleId=0; currentRuleId<maxRulesIndex; currentRuleId++) {
     if (sequences[currentRuleId][currentRuleId].numberOfDatas () != 0) {
     printf ("Rule %d - %d sequences mini: %d maxi: %d average: %f\n",
     currentRuleId,
     sequences[currentRuleId][currentRuleId].numberOfDatas (),
     sequences[currentRuleId][currentRuleId].mini (),
     sequences[currentRuleId][currentRuleId].maxi (),
     sequences[currentRuleId][currentRuleId].average ());
     }
  } */
  
  printf ("Matrice de contraintes (occurence):\n");
  printf ("    ");
  for (i=0; i<=maxRulesIndex; i++) {
    if (used[i])
      printf ("  %.2d  ", i+1);
  }
  printf ("\n");
  
  for (i=0; i<=maxRulesIndex; i++) {
    if (used[i]) {
      printf (" %.2d ", i+1);
      for (j=0; j<=maxRulesIndex; j++) {
	if (used[j]) {
	  if (sequences[i][j].numberOfDatas() == 0) {
	    printf ("|     ");
	  }
	  else if (sequences[i][j].numberOfDatas () < 10) {
	    printf ("|   %d ", sequences[i][j].numberOfDatas ());
	  }
	  else if (sequences[i][j].numberOfDatas () < 100) {
	    printf ("|  %d ", sequences[i][j].numberOfDatas ());
	  }
	  else {
	    printf ("| %d ", sequences[i][j].numberOfDatas ());
	  }
	}
      }
      printf ("|\n");
    }
  }
  printf ("\n");
  
  printf ("Matrice de contraintes (average):\n");
  printf ("    ");
  for (i=0; i<=maxRulesIndex; i++) {
    if (used[i])
      printf ("   %.2d   ", i+1);
  }
  printf ("\n");
  
  for (i=0; i<=maxRulesIndex; i++) {
    if (used[i]) {
      printf (" %.2d ", i+1);
      for (j=0; j<=maxRulesIndex; j++) {
	if (used[j]) {
	  if (sequences[i][j].average () == 0) {
	    printf ("|       ");
	  }
	  else if (sequences[i][j].average () < 10) {
	    printf ("|  %.2f ", sequences[i][j].average ());
	  }
	  else {
	    printf ("| %.2f ", sequences[i][j].average ());
	  }
	}
      }
      printf ("|\n");
    }
  }
  printf ("\n");
  
  /* for (i=0; i<nodesNumber; i++) {
     node[i]->debugSequences ();
     } */
}

void Graph::computePatterns (InfoBox* infoBox)
     /* si infoBox est NULL, pas d'affichage des stats */
{
  int i;
/* USELESS  int j,k; */
/* USELESS  int currentRuleId; */
  int total;
  int currentPatternId;
  int displayStats;
  
  if (arePatternsComputed) {
    return;
  }
  
  if (infoBox == NULL) {
    displayStats = 0;
  }
  else {
    displayStats = 1;
  }

  initRulesSequences ();
  
  /* on delete la liste des patterns existants (s'il y en a) */
  
  if (patternsNumber != 0) {
    delete[] patternList;
  }
  
  /* on cherche les patterns */
  
  total = 0;
  patternsNumber = 0;
  
  for (i=1; i<nodesNumber; i++) {
    if (displayStats) {
      infoBox->setStatWindow ((double)(i)/(double)(nodesNumber));
      if (infoBox->refreshStatWindow () == 1) {
	break;
      }
    }
    node[i]->initPattern (archSize);
    if (node[i]->searchPattern (node[i], 0)) {
      /* nouveau pattern trouve */
      total+=node[i]->pattern->cellsNumber ();
      patternsNumber++;
    }
    else {
      /* rien trouve */
    }
  }
  
  /* on a maintenant extrait patternsNumber patterns */

  /* on reserve de la memoire pour stocker les patterns */
  
  patternList = new Pattern*[patternsNumber];
  
  /* on cree la patternList */

  currentPatternId = 0;
  
  for (i=1; i<nodesNumber; i++) {
    if (node[i]->isPatternInitialized) {
      /* on positionne le pattern currentPatternId */
      patternList[currentPatternId] = node[i]->pattern;
      currentPatternId++;
    }
  }

  /* on code tous les patterns */

  for (i=0; i<patternsNumber; i++) {
    patternList[i]->setEncoding ();
  }
  
  /* mise a jour de la fitness */
  
  arePatternsComputed = 1;
  isComplexityComputed = 1;
  if (patternsNumber != 0) {
    complexityFitness = (double)total/(double)(patternsNumber)/40.0;
  }
  else {
    complexityFitness = 0.0;
  }
  
}

void Graph::freeAllPatterns ()
{
  int i;

  if (patternsNumber != 0) {
    delete[] patternList;
  }

  for (i=0; i<getNodesNumber (); i++) {
    (nodeAtIndex (i))->freePattern ();
  }

  arePatternsComputed = 0;
  patternsNumber = 0;
}

void Graph::computePatternMatching (InfoBox* infoBox)
{
  correlationType correlation;
  int angle, dx, dy, dz; 
  int i, j, k;
  int displayStats;
  RandomGenerator randomGenerator (PSEUDO, 0, 
				   (int)((float)patternsNumber*3/4));
  int maxTests;
  float total;
/* USELESS  float totalCells; */
/* USELESS  float differentCells; */
  
  if (infoBox == NULL)
      displayStats = 0;
  else
      displayStats = 1;

    
  if (patternsNumber != 0)
      {
	  maxTests = 200;
	  total = 0.0;
	  
	  for (k=0; k<maxTests; k++)
	      {
		  i = randomGenerator.getNumber ();
		  j = randomGenerator.getNumber ();
      
		  if (displayStats)
		      {
			  infoBox->setStatWindow 
			      ((double)(k)/(double)(maxTests));
			  if (infoBox->refreshStatWindow () == 1)
			      break;
		      }
      
		  total += patternList[i]->matchWithOtherPattern (patternList[j],
								  &correlation, 
								  &angle, &dx, &dy, &dz);
	      }
    
	  isCoherenceComputed = 1;
	  coherenceFitness = (double)total/(double)k;
	  isPatternMatchingDone = 1;
      }
  
  else
      {
	  // printf ("ERROR: unable to compute pattern matching: no patterns\n");
	  isCoherenceComputed = 1;
	  /* coherenceFitness = (double)total/(double)k; ??? */ 
      }
  
}

double Graph::getComplexityFitness ()
{
  if (!isComplexityComputed) {
    computePatterns (NULL);
  }
  
  return complexityFitness;
}

double Graph::getCoherenceFitness ()
{
  if (!isComplexityComputed) {
    computePatterns (NULL);
  }
  
  if (!isCoherenceComputed) {
    computePatternMatching (NULL);
  }
  
  return coherenceFitness;
}

/* Objet Sequences */

Sequences::Sequences ()
{
  total = 0;
  dataNumber = 0;
  miniValue = -1;
  maxiValue = -1;
  averageValue = 0;
}

Sequences::~Sequences ()
{
  /* rien a faire ici */
}

void Sequences::addData (int aData)
{
  total += aData;
  dataNumber++;
  if ((aData < miniValue) || (miniValue == -1)) {
    miniValue = aData;
  }
  if ((aData > maxiValue) || (maxiValue == -1)) {
    maxiValue = aData;
  }
  averageValue = (double)total/(double)dataNumber;
}

int Sequences::numberOfDatas ()
{
  return dataNumber;
}

int Sequences::mini ()
{
  return miniValue;
}

int Sequences::maxi ()
{
  return maxiValue;
}

double Sequences::average ()
{
  return averageValue;
}