chameleon.c

变色龙层次聚类算法

    printf("cannot allocate memory for node0->points. \n");
    return -1;
  }
  for(i=0; i<groupLength[group1]; i++){
    node1->points[i]=groups[group1][i];
  }/* end for */

  /* Now compute absolute inter-connectivity and 
   * absolute closeness.
   */
  computeAIC_AC(node0, node1, &aic, &ac);


  if(ac > 0.0 && aic > 0.0){
    /* now need to compute the internal connectivity of 
     * node0 and node1.
     */
    left0 = (struct node *)calloc(1, sizeof(struct node));
    right0 = (struct node *)calloc(1, sizeof(struct node));
    if(left0==NULL || right0==NULL){
      printf("cannot allocate memory for node0 and node1. \n");
      return -1;
    }/* end if */
    
    if(cutNode(node0, left0, right0, ubFactor)<0)
      {
	printf("cut node0 error! \n");
	return -1;
      }
    
    /* Now compute internal inter-connectivity and 
     * absolute closeness of node0.
     */
    computeAIC_AC(left0, right0, &aic0, &ac0);
    
    left1 = (struct node *)calloc(1, sizeof(struct node));
    right1 = (struct node *)calloc(1, sizeof(struct node));
    if(left1==NULL || right1==NULL){
      printf("cannot allocate memory for node0 and node1. \n");
      return -1;
    }/* end if */
    
    if(cutNode(node1, left1, right1, ubFactor)<0)
      {
	printf("cut node0 error! \n");
	return -1;
      }
    
    /* Now compute internal inter-connectivity and 
     * absolute closeness of node0.
     */
    computeAIC_AC(left1, right1, &aic1, &ac1);
    
    
    /* Now compute Relative Inter-Connectivity */
    ri = (aic*2)/(aic0+aic1);

    /* Now compute Relative Closeness */
    rc = ac*(node0->numberPoints+node1->numberPoints)/((node0->numberPoints)*ac0 + 
						       (node1->numberPoints)*ac1);
  
    /* 
       printf("aic is %f, aic0 is %f, aic1 is %f \n", aic, aic0, aic1);
       printf("ac is %f, ac0 is %f, ac1 is %f \n", ac, ac0, ac1);
       printf("node0->numberPoints is %d, node1->numberPoints is %d \n", 
       node0->numberPoints, node1->numberPoints);
    */

    if(ri>10000 || rc>10000){
      printf("ri or ic is wrong! ri= %f, rc= %f \n", ri, rc);
      return -1;
    }
    /* The final goodness according to p11 of the CHAMELEON paper */
    /* printf("ri is %f, rc is %f. \n", ri, rc);*/
    goodness = ri*pow((double)rc, ALPHA);


    free(node0->points);
    free(node1->points);
    free(left0->points);
    free(right0->points);
    free(left1->points);
    free(right1->points);
    free(left0);
    free(right0);
    free(left1);
    free(right1);
    free(node0);
    free(node1);
  }
  else{
    ri=0.0;
    rc=0.0;
    /* printf("ri is %f, rc is %f. \n", ri, rc);*/
    goodness = 0.0;

    free(node0->points);
    free(node1->points);
    free(node0);
    free(node1);
  }/* end if...else...*/

  return goodness;
}/* end computeGoodness */

/****************************************************************
 * Name    : merge                                              *
 * Function: merge 2 clusters.                                  *
 * Input   : int group0 -- group0's index.                      *
 *           int group1 -- group1's index.                      *
 * Output  :  int                                               *
 ****************************************************************/
static int merge(int group0, int group1){
  int i, j;

  groups[group0] = (int *)realloc(groups[group0], 
				  (groupLength[group0]+groupLength[group1])*sizeof(int)) ;
  if(groups[group0]==NULL){
    printf("cannot enlarge groups[group0]; \n");
    return -1;
  }
  /* move elements from group1 to group0 */
  for(i=groupLength[group0], j=0; i<(groupLength[group0]+groupLength[group1]); i++, j++){
    groups[group0][i] = groups[group1][j];
  }/* end for */

  free(groups[group1]);
  groups[group1]=NULL;

  groupLength[group0] = groupLength[group0]+groupLength[group1];
  groupLength[group1] = 0;

  return 1;
}/* end merge */


/****************************************************************
 * Name    : phase2                                             *
 * Function: phase 2, merging clusters.                         *
 * Input   : void                                               *
 * Output  :  int                                               *
 ****************************************************************/
static int phase2(){
  int clusterNO;
  float goodness;
  int group0, group1;
  int i, j;

  /* 
   * groupIndex at the end of phase 1 gives the
   * number of initial clusters after phase 1.
   */
  clusterNO = groupIndex;

#ifdef Debug_phase2
  printf("At the beginning of phase 2, clusterNO is %d, stopClusterNO is %d \n", 
	 clusterNO, stopClusterNO);
#endif

  if(clusterNO<stopClusterNO){
    printf("at the beginning of phase 2, clusterNO is already smaller than stopClusterNO! error\n");
    return -1;
  }/* end if */

  bestGoodness = 0.0;
  while(clusterNO>stopClusterNO){
    /* the for loop decide which two groups need to be merged */
    for(i=0; i<(groupIndex-1); i++){
      if(groupLength[i]!=0){
	group0 = i;
	for(j=(group0+1); j<groupIndex; j++){
	  if(groupLength[j]!=0){
	    group1 = j;
	    
	    /* Now compute the goodness for merging the two clusters */
	    goodness = computeGoodness(group0, group1);
	    if(goodness > 10000 || goodness < 0){
	      printf("group0 is %d, group1 is %d,!!!\n", 
		     group0, group1);
	      printf("goodness is %f, some thing wrong!\n", goodness);

	      return -1;
	    }
	    /* printf("group0 is %d, group1 is %d, goodness is %f !!! \n", 
	       group0, group1, goodness);
	    */
	    if(goodness > bestGoodness)
	      {
		bestGoodness = goodness;
		mergeGroups[0] = group0;
		mergeGroups[1] = group1;
	      }
	  }/* end if */
	}/* end for j*/
      }/* end if */
    }/* end for i*/
    
    /* now merge the two selected groups */
    if(mergeGroups[0]==mergeGroups[1]){
      printf("mergeGroups[0]==mergeGroups[1]==%d, something wrong! \n", 
	     mergeGroups[0]);
      return -1;
    }

    if(bestGoodness == 0.0){
      printf("bestGoodness reached 0 \n");
      return 1;
    }

    printf("best goodness is %f \n", bestGoodness);

    if(merge(mergeGroups[0], mergeGroups[1])<0)
      {
	printf("merge error! \n");
	return -1;
      }

    clusterNO--;
    
#ifdef Debug_phase2
    printf("In the while loop,  clusterNO is %d, stopClusterNO is %d, mergeGroups[0] is %d, mergeGroups[1] is %d. \n", 
	   clusterNO, stopClusterNO, mergeGroups[0], mergeGroups[1]);
#endif

    bestGoodness = 0.0;
    mergeGroups[0]=0;
    mergeGroups[1]=0;

  }/* end while */

  return 1;

}/* end phase2 */


/****************************************************************
 * Name    : belongs                                            *
 * Function: decide whether a point belongs to a node.          *
 * Input   : struct node * sourceNode --  the node              *
 *           int point -- the point                             *
 * Output  :  int                                               *
 ****************************************************************/
static int belongs(struct node * sourceNode, int point){
  int i, tmp;

  tmp = 0;

  for(i=0; i<sourceNode->numberPoints; i++){
    if(sourceNode->points[i] == point)
      {
	tmp = 1;
	break;
      }/* end if */
  }

  return tmp;

}/* end belongs */

/****************************************************************
 * Name    : cutNode                                            *
 * Function: cut a node of the tree into two parts, and thus    *
 *           return the left pointer and the right pointer.     *
 * Input   : struct node * source --  the source node to be cut.*
 *           struct node * left -- the left resulting node.     *
 *           struct node * right -- the right resulting node.   *
 * Output  :  int                                               *
 ****************************************************************/
/* May 30th, Weinan: one problem I have just found is that I need to 
 * establish a kind of checking table to connect my global vertex 
 * number with local vertex number to be used in the HMETIS_PartRecursive
 * function calling.
 */
static int cutNode(struct node * source, struct node * left, struct node * right, int ub){
  int nvtxs, nhedges;
  int * vwgts = NULL;
  int * eptr = NULL;
  int * eind = NULL;
  int * hewgts = NULL;
  int nparts;
  int ubfactor;
  int * options = NULL;
  int * part = NULL;
  int * edgecut = NULL;

  int tmp, i, j, k, l, tmpNO, flag, found, group0, group1, index0, index1;
  /* this is the table used to correspond global vertice number with 
   * local vertice number.
   */
  int * checkTable;

  /* number of vertices */
  nvtxs = source->numberPoints;
  checkTable = (int *)calloc(nvtxs, sizeof(int));
  if(checkTable == NULL)
    {
      printf("cannot allocate memory for checkTable! \n");
      return -1;
    }
  /* load the vertice's number into the checkTable */
  for(i=0; i<nvtxs; i++){
    checkTable[i] = source->points[i];
  }/* end for */

  /* number of edges */
  tmp=0;
  for(i=0; i<nvtxs; i++){
    for(j=0; j<points[source->points[i]].length; j++){
      /* decide a point whether belongs to a node */
      if(belongs(source, points[source->points[i]].edges[j].pointNO)==1) /* !!!!!! */
	tmp++;
    }/* end for j*/
  }/* end for i*/
  nhedges = tmp;

  /* weight of the vertices, because 
   * the vertices are not weighted, so
   * according to p13 of the hMETIS manual,
   * vwgts can be NULL.
   */
  vwgts = NULL;

  /* eptr and eind */
  /* because all my edges are of
   * size 2, so it makes things 
   * easier.
   */
  tmp = nhedges+1;
  eptr = (int *)calloc(tmp, sizeof(int));
  if(eptr == NULL)
    {
      printf("cannot allocate memory for eptr! \n");
      return -1;
    }
  for(i=0; i<tmp; i++)
    eptr[i]=i*2;

  /* when loading eind, need to check the checkTable */
  eind = (int *)calloc(2*tmp, sizeof(int));
  if(eind == 0)
    {
      printf("cannot allocate memory for eind! \n");
      return -1;
    }
  k = 0;
  for(i=0; i<nvtxs; i++){
    for(j=0; j<points[source->points[i]].length; j++){
      /* decide a point whether belongs to a node */
      if(belongs(source, points[source->points[i]].edges[j].pointNO)==1){
	/* eind[k] = source->points[i];*/
	eind[k] = i;
	k++;
	/* eind[k] = points[source->points[i]].edges[j].pointNO; */
	tmpNO =  points[source->points[i]].edges[j].pointNO;
	flag = 0;
	for(l=0; l<nvtxs; l++){
	  if (tmpNO == checkTable[l]){
	    flag = 1;
	    found = l;
	    break;
	  }
	}
	if(flag == 1)
	  eind[k] = found;
	else
	  {
	    printf("some thing wrong with checkTable! \n");
	    return -1;
	  }/* end if...else... */
	k++;      
      }/* end if */
    }/* end for j*/
  }/* end for i*/
  if(k != 2*nhedges){
    printf("some thing wrong with eind! \n");
    return -1;
  }/* end if */

  /* hewgts: an array of size nhedges that stores the weight
   * of the hyperedges.
   */
  hewgts = (int *)calloc(nhedges, sizeof(int));
  if(hewgts == 0)
    {
      printf("cannot allocate memory for hewgts! \n");
      return -1;
    }
  k = 0;
  for(i=0; i<nvtxs; i++){
    for(j=0; j<points[source->points[i]].length; j++){
      /* decide a point whether belongs to a node */
      if(belongs(source, points[source->points[i]].edges[j].pointNO)==1){
	/*!!!!!! here I have to do a cast becasue now similarity is a double */
	hewgts[k] = (int)(points[source->points[i]].edges[j].similarity*DIAG);
	k++;      
      }/* end if */
    }/* end for j*/
  }/* end for i*/
  if(k != nhedges){
    printf("some thing wrong with hewgts! \n");
    return -1;
  }/* end if */

  /* nparts: number of desired partitions */
  nparts = 2;

  /* ubfactor: relative imbalance factor */
  ubfactor = ub;

  /* options */
  /* note that is options[0]=0,
   * then default values are used.
   */
  options =  (int *)calloc(9, sizeof(int));
  if(options == 0)
    {
      printf("cannot allocate memory for options! \n");
      return -1;
    }
  options[0] = 0;

  /* part */
  part =  (int *)calloc(nvtxs, sizeof(int));
  if(part == 0)
    {
      printf("cannot allocate memory for part! \n");
      return -1;
    }

  /* edgecut */
  edgecut = (int *)calloc(1, sizeof(int));
  if(edgecut == 0)
    {
      printf("cannot allocate memory for edgecut! \n");
      return -1;
    }


  HMETIS_PartRecursive(nvtxs, nhedges, vwgts, eptr, eind, hewgts, 
		       nparts, ubfactor, options, part, edgecut);

  group0 = 0;
  group1 = 0;
  for(i=0; i<nvtxs; i++){
    if(part[i] == 0)
      group0++;
    else if(part[i] == 1)
      group1++;
    else
      {
	printf("something wrong with part. \n");
	return -1;
      }/* end if..else...*/
  }/* end for */

  left->numberPoints = group0;
  right->numberPoints = group1;
  left->points = (int *)calloc(group0, sizeof(int));
  if(left->points == NULL)
    {
      printf("cannot allocate memory for left->points \n");
      return -1;
    }
  right->points = (int *)calloc(group1, sizeof(int));
  if(right->points == NULL)
    {