/* check whether the spot is free	   and start search for near free position,            if necessary	*/            pos = env->searchdroppingsite(pos, dat); 	/* add item to grid, possibly remember it */	env->add(pos, dat);	if (memsize > 0) {	    mem->remember(dat, pos, alpha);	}	dat = FREE;      			/* now look for a data item to pick up */	while (dat == FREE) {	    /* try next candidate */	    position<int> temp_pos = env->next();	    int temp_data = (*env)(temp_pos[0],temp_pos[1]);	    /* examine local neighbourhood */	    double fvalue = env->f(temp_data, temp_pos, radius, int(nsize), scalefactor, par->clusterphase);	   	    /* probabilistic picking decision */	    if (ran0(&idum) < th_pick(fvalue)) {		/* pick up data item and remove it from the grid */		pos = temp_pos;		lastpos = pos;		dat = temp_data;		env->remove(pos);		return;	    }	}    }    else {	/* keep track of failures for alpha adaptation */	if (adapt == TRUE) {	    failure_ctr++;	}    }}/* computation of picking probability */inline double accl::ant::th_pick(double f) {      if (par->newprob == FALSE) {	return square(kp / (kp + f));    }    else {	if (f <= 1.0) return 1.0;	return 1.0 / (pow(f,2));    }}/* computation of dropping probability */inline double accl::ant::th_drop(double f) {     if (par->newprob == FALSE) {	return square(f / (kd + f));    }    else {	if (f >= 1.0) return 1.0;	return pow(f,4);    }}// well.... inline double accl::ant::square(double x) {    return x*x;}inline int accl::ant::ithmemory(int i) {    return mem->index[i];}inline int accl::ant::ithmemoryx(int i) {    return mem->pos[i][0];}inline int accl::ant::ithmemoryy(int i) {    return mem->pos[i][1];}/***************************************************     accl***************************************************//* constructor for an ACCL instance   memory allocation and computation of normalization factor*/accl::accl(conf * par, databin<USED_DATA_TYPE> * bin, evaluation * e):clalg(par, bin, e){    par->radius = par->initradius;    /* construction of the grid */    env = new grid(par, bin);      /* par->mu must have been computed before       the construction of the antcolony    */    colony = new antcolony(par, env);}/* destructor:   deletion of grid and ant colony*/accl:: ~accl() {    delete env;    delete colony;}/* initialisation of the algorithm:   requires the intialization of both the environment and the colony   environment muts be initialized first*/void accl::init() {       par->radius = par->initradius;    generation_ctr = 0;    env->init();}/* main routine of the algorithm   requires previous initialization*/void accl::run() {    par->clusterphase = TRUE;   /* let ants pick up a first data item */    colony->init();     generation_ctr++;        time_t start;    time_t end;     /* main loop        consisting of repeated phases of the algorithm       (subdivision used for evaluation purposes)    */    while (generation_ctr <= par->generations) {//	cout << generation_ctr << endl;	/* increase of radius */	if (generation_ctr % (int)(0.21*par->generations) == 0) {	    double nsize = 0.0;	    par->radius += 1;		    if (par->weighting == TRUE) {		for (int dx=0; dx <= par->radius; dx++) {		    for (int dy=1; dy <= par->radius; dy++) {			nsize += 4*exp(-double(max(dx, dy))/double(par->radius));		    }		}	    }	    if (par->increase == TRUE) {		for (int i=0; i<par->colonysize; i++) {   		    (*colony)[i].radius += 1;		    if (par->weighting == FALSE) {				    }		    else {			(*colony)[i].nsize = nsize;		    }		}	    }	}		/* clustering phases */	if (par->clustering == TRUE) {	    if (generation_ctr > 0.7*par->generations) {		if (generation_ctr > 0.55*par->generations) {		    par->clusterphase  = TRUE;		    		}		else if (generation_ctr > 0.45*par->generations) {		    par->clusterphase = TRUE;		}		else par->clusterphase = FALSE;	    }	}	// one generation consists of many individual actions	colony->cometolive();		generation_ctr++;    }     /* end of main loop */    colony->finish();}/* agglomerative construction of the clustering   starting from spatial mapping*/void accl::constructclustering() {    cout << "Cluster construction (using agglomerative clustering)" << endl;    /* distance matrixto reduce computational complexity */    tmatrix<USED_DATA_TYPE> dmatrix(par->binsize);    tmatrix<USED_DATA_TYPE> slinkmatrix(par->binsize);    int * currentsize = new int[par->binsize];    /* memory allocation */    int ** clusters = new int * [par->binsize];    position<double> ** docindex = new position<double>*[par->binsize];    for (int i=0; i<par->imax; i++) {	for (int j=0; j<par->jmax; j++) {	    if ((*env)(i,j) != FREE) {		docindex[(*env)(i,j)] = new position<double>(i, j);	    }	}    }    int num = par->binsize;    position<double> * centre = new position<double>[par->binsize];    /* initially each cluster contains one data element */    for (int i=0; i<par->binsize; i++) {	currentsize[i] = 1;	clusters[i] = new int[par->binsize+1];	clusters[i][0] = i;	for (int j=1; j<par->binsize+1; j++) {	    clusters[i][j] = FREE;	}    }    /* intialize distance matrix */    for (int i=0; i<par->binsize; i++) {	for (int j=0; j<i; j++) {	    dmatrix(i,j) = minlink(clusters[i], clusters[j], docindex);	    slinkmatrix(i,j) = dmatrix(i,j) / 2.0;	}    }    double minval = 0;    /* merge clusters until only one cluster is left       or stopping criteria is met */    while ( num > 1 ) {	minval = par->imax*par->jmax;	int minindex1 = 0;	int minindex2 = 0;    	/* determine the two closest clusters */	for (int i=0; i<par->binsize; i++) {	    if (clusters[i] == NULL) continue;	    for (int j=0; j<i; j++) {		if (clusters[j] == NULL) continue;			// This is the spot where we save computations, we just need to update the 		// matrix after each merging procedure (only one row though!!)		double d = dmatrix(i,j);		if (d < minval) {		    minval = d;		    minindex1 = i;		    minindex2 = j;		}	    }	}		//	if (minval >= 5) break;		if (minval >= par->radius) {	  break;	}	/* merge the two clusters */	int ptr=0;	while (clusters[minindex1][ptr] != FREE) {	    ptr++;	}	int i=0;	while (clusters[minindex2][i] != FREE) {	    clusters[minindex1][ptr++] = clusters[minindex2][i];	    i++;	}	delete [] clusters[minindex2];	clusters[minindex2] = NULL;	currentsize[minindex1] += currentsize[minindex2];	num--;	/* update row of distance matrix */	double quot; 	for (int i=0; i<par->binsize; i++) {	    if ( i!= minindex1 && clusters[i] != NULL) {		slinkmatrix(minindex1, i) = min(slinkmatrix(minindex1,i),slinkmatrix(minindex2,i));			if (currentsize[minindex1] < currentsize[i]) {		    quot = double(currentsize[minindex1]) / double(currentsize[i]);		}		else {		    quot = double(currentsize[i]) / double(currentsize[minindex1]);		}		double weight = 1.0+log10(1.0+9.0*quot);		dmatrix(minindex1,i) = slinkmatrix(minindex1,i)*weight;	    }	}    }     /* detect clusters that are too small (smaller or equal to 2) */    int * toosmall = new int[par->binsize];    int toosmallctr = 0;     for (int i=0; i<par->binsize; i++) {	if (clusters[i] == NULL) continue;	else {	    int j=0; 	    	 	    while (clusters[i][j] != FREE) {		j++;	    }	    if (j <= 2) { 		toosmall[toosmallctr] = i;		toosmallctr++;    		num--;		    }	    else cout << j << endl;	}    }    /* remove the small clusters by merging them with the closest cluster */    int minindex = 0;    double mind = par->imax*par->jmax;    double d;    for (int i=0; i<toosmallctr; i++) {	cout << "Small cluster " << toosmall[i] << " corrected" << endl;	/* detect closest cluster */	for (int j=0; j<par->binsize; j++) {	start:	    if (toosmall[i] == j) continue;	    for (int k=0; k<toosmallctr; k++) {		if (toosmall[k] == j) {		    if (j < par->binsize-1) {			j++;			goto start;		    }		    else goto ende;		}	    }	    if (clusters[j] != NULL) {		d = minlink(clusters[toosmall[i]], clusters[j], docindex);		if (d < mind) {		    minindex = j;		    mind = d;		}	    }	}    ende:	cout << "Merging cluster " << minindex << " and " << toosmall[i] << endl;	/* merge them */    	int ptr=0;	while (clusters[minindex][ptr] != FREE) {	    ptr++;	}	int j=0;	while (clusters[toosmall[i]][j] != FREE) {	    clusters[minindex][ptr++] = clusters[toosmall[i]][j];	    j++;	}	delete [] clusters[toosmall[i]];	clusters[toosmall[i]] = NULL;    }    /* generate clustering object */    clust = new clustering(par);    clust->init(num);    int ctr = 0;    for (int i=0; i<par->binsize; i++) {	if (clusters[i] == NULL) continue;	else {	    int j=0; 	    data<USED_DATA_TYPE> datacentre(par);	 	    while (clusters[i][j] != FREE) {		int index = clusters[i][j];		(*clust)[index] = ctr;		datacentre.add((*bin)[index]);		j++;	    }	    datacentre.div(j);	    clust->newcentre(ctr, &datacentre);	    ctr++;	}    }    cout << "Number of clusters: " << num << endl;    }/* modified single link linkage metric */double accl::minlink(int * cluster1, int * cluster2, position<double> ** docindex) {    double mind = par->imax*par->jmax;    int i=0;    int ctr1 = 0;    int ctr2 = 0;    int ctr = 0;     while (cluster1[i] != FREE) {	int j=0;	int d1 = cluster1[i];	ctr1++;	ctr2 = 0;	while (cluster2[j] != FREE) {	    int d2 = cluster2[j];	    double d = env->spatialdistance(docindex[d1], docindex[d2], TRUE);	    mind = min(mind, d);	    j++;	    ctr2++;	    ctr++;	}	i++;    }    /* weighting factor */    double quot = double(ctr1) / double(ctr2);    if (ctr1 > ctr2) {	quot = double(ctr2) / double(ctr1);    }    quot = 1.0 + log10(1.0+9.0*quot);	    return mind*quot;}