From 2723624a67e454540cff46760f2080c6fbb45058 Mon Sep 17 00:00:00 2001
From: Dario Izzo <dario.izzo@googlemail.com>
Date: Wed, 17 Dec 2008 00:24:45 +0100
Subject: [PATCH] DiGMO simulateor case added

---
 SolversThreads/SolversThreads.cpp |   6 +-
 main.cpp                          | 320 +++++++++++++++++++++-----------------
 2 files changed, 179 insertions(+), 147 deletions(-)

diff --git a/SolversThreads/SolversThreads.cpp b/SolversThreads/SolversThreads.cpp
index 88ed117..7e558f2 100644
--- a/SolversThreads/SolversThreads.cpp
+++ b/SolversThreads/SolversThreads.cpp
@@ -254,7 +254,7 @@ void *ASAthread(void *data)
    vector <int> picks;
    vector<double> LB,UB;
    ASAalgorithm ASA;
-   Pk::Random32 rng;
+   Pk::Random32 rng;
    GOProblem* problem;
 
 
@@ -264,8 +264,8 @@ void *ASAthread(void *data)
    pthread_mutex_lock (PtrTP->TPmutex);
         rng = Pk::Random32(PtrTP->randomSeed);
 		deme=PtrTP->Ptr_pop->extractRandomDeme(PtrTP->NP,picks, rng);
-		problem = PtrTP->problem;
-		problem->getBounds(LB, UB);
+		problem = PtrTP->problem;
+		problem->getBounds(LB, UB);
 		unsigned int temp;
 		temp = rng.next();
 		ASA.initASA(PtrTP->generations,LB.size(),PtrTP->Ts,PtrTP->Tf, rng.next());
diff --git a/main.cpp b/main.cpp
index 5ebc81c..f24e312 100644
--- a/main.cpp
+++ b/main.cpp
@@ -53,7 +53,7 @@ while (choice != -1) {
 
 		//we choose the algorithm
 
-		cout << "Choose: 1-ASA, 2-PSO, 3-MPSO, 4-DE, 5-SGA, 6-IslandModel(ring), Default-DiGMO: ";
+		cout << "Choose: 1-ASA, 2-PSO, 3-MPSO, 4-DE, 5-SGA, 6-IslandModel(ring), 7-DiGMO(simulator): ";
 		cin >> choice;
 
 		switch (choice){
@@ -696,190 +696,222 @@ while (choice != -1) {
 			}
 			break;
 
-            case 9:
-            {
-                Individual x;
-                x.createRandomIndividual(LB,UB,rng);
-                cout << x << endl;
-                x[0]=666;
-                cout << x << endl;
-
-                Population pop;
-                pop.createRandomPopulation(LB,UB,5,rng);
-                pop.evaluatePopulation(problem);
-                cout << pop[3] << endl;
-                pop[3][0]=666;
-                cout << pop[3] << endl;
-                cout << pop << endl;
-
-            }
-            break;
-
-			default:
-			{
-			int NUM_THREADS=2;
-			//We create a random population and evaluate its fitness (this is done by the main thread - i.e. no parallelization here)
-			pop.createRandomPopulation(LB,UB,20*5, rng);
-			pop.evaluatePopulation(problem);
+            case 7:
+           {
 
-			//We instanciate the objects needed for pthreads allocating memory for the threads array
-			pthread_t *threads;
-			threads = new pthread_t [NUM_THREADS];
-			//pthread_attr_t attr;
-			//pthread_attr_init(&attr);
-			//pthread_attr_setdetachstate(&attr, PTHREAD_CREATE_JOINABLE);
-			pthread_mutex_t mutex;
-			pthread_cond_t exitcond;
-			pthread_mutex_init(&mutex, NULL);
-			pthread_cond_init (&exitcond, NULL);
+			//Experiment Settings
+			int NP = 200;					//main population size
+			int trials = 20;				//number of trials
+			int gen = 500;					//generations per algorithm call
+			double F = 0.8;					//F in DE
+			double CR = 0.8;				//CR in DE
+			int strategy = 2;				//DE startegy
+			int islandsN  = 1;				//Number of Islands
 
+			//stopping criteria
+			int itermax = 120;				//Maximum number of iterations allowed (i.e. output printed on the screen)
 
+			//Experiment Outputs
+			double mean = 0, dev= 0, max = 0, min = 200000000;
+			vector <double> results;
 
+			//Main cycle creating threads
+			for (int i=0;i<trials;i++){
+				cout << "\nTrial number #" << i+1 << endl;
+				//We create a random population in each island (this is done by the main thread - i.e. no parallelization here)
+				Population IslandPop;
+				vector<GOProblem*> parallelProblems(islandsN);
+				for (int i=0;i<islandsN;i++) { parallelProblems[i] = new problem_type(); } //It is necessary to create a different object per CPU as to make sure to access different memory
+																							//locations when the problem is mga or mga-1dsm for th r,v,and DV variables in the mgaproblem structures
+				IslandPop.createRandomPopulation(LB,UB,NP, rng);
+				IslandPop.evaluatePopulation(*parallelProblems[0]);							//all the problems are identical.... evaluation is done for the [0] one.
 
-			//We allocate memory for the isActive array containing the status ofeach thread	and we initialise it to false (no threads opened)
-			bool *isActive;
-			isActive = new bool[NUM_THREADS];		//all threads are inactive
-			for (int i=0;i<NUM_THREADS;i++) {isActive[i]=false;}
 
-			//We allocate memory and initialise the data threadParam array containing the information to be passed to the different threads
-			threadParam *data;
-			data = new threadParam [NUM_THREADS];
-			for (int i=0;i<NUM_THREADS;i++){
-				data[i].threadID = i;
+				//We instanciate the objects needed for pthreads allocating memory for the threads array
+				pthread_t *threads;
+				threads = new pthread_t [islandsN];
+				pthread_mutex_t mutex;
+				pthread_cond_t exitcond;
+				pthread_mutex_init(&mutex, NULL);
+				pthread_cond_init (&exitcond, NULL);
 
-				data[i].problem = &problem;
-				data[i].Ptr_pop = &pop;
+				//We allocate memory for the isActive array containing the status ofeach thread	and we initialise it to false (no threads opened)
+				bool *isActive;
+				isActive = new bool[islandsN];
+				for (int i=0;i<islandsN;i++) isActive[i]=false;  //all threads are inactive
 
-				data[i].generations = 500;
+				//We allocate memory and initialise the data threadParam array containing the information to be passed to the different threads
+				threadParam* data;
+				data = new threadParam [islandsN];
+				for (int i=0;i<islandsN;i++){
+					data[i].threadID = i;
 
-				//Initialise DE specific data
-				data[i].strategy = 2;
-				data[i].F = 0.8;
-				data[i].CR = 0.8;
+					data[i].problem = parallelProblems[i];
+					data[i].Ptr_pop = &IslandPop;			//Only one population in DiGMO
 
-				//Initialise PSO specific data
-				data[i].omega = 0.5;
-				data[i].eta1 = 2.0;
-				data[i].eta2 = 2.0;
-				data[i].vcoeff = 1000;
 
-				data[i].isActive = &isActive[i];
-				data[i].TPmutex = &mutex;
-				data[i].exit = &exitcond;
+					data[i].NP = 20;
+					data[i].generations = gen;
 
-				data[i].Ptr_log = &logfile;
-			}
+					//Initialise DE specific data
+					data[i].strategy = 1;
+					data[i].F = F;
+					data[i].CR = CR;
 
-			//Main cycle creating threads
-			int globalCount=0;
-			int count1=0,count2=0;
-			int rc;
-			int algorithmchoice = 0;
-
-			pthread_mutex_lock (&mutex);
-			time(&start);
-			while (globalCount<120){
-			//while(pop.extractBestIndividual().getFitness()!=0){
-					for (int i=0;i<NUM_THREADS;i++){
-					if ( !*(data[i].isActive) ){
-							cout << "main():" << endl <<  "\t\t\tcreating thread, ID " << i << endl;
-							algorithmchoice = 4;//r250()%10;
-							if (algorithmchoice == 0){
-								count1++;
+					//Initialise PSO/MPSO specific data
+					data[i].omega = 0.65;
+					data[i].eta1 = 2.0;
+					data[i].eta2 = 2.0;
+					data[i].vcoeff = 1;
+					data[i].nswarms = 4;
+
+					//Initialise SGA specific data
+					data[i].CRsga = 0.95;
+					data[i].M = 0.2;
+					data[i].insert_best =1;
+
+					//Initialise ASA specific data
+					data[i].Ts = 1;
+					data[i].Tf = data[i].Ts/1000;
+					data[i].isActive = &isActive[i];
+					data[i].TPmutex = &mutex;
+					data[i].exit = &exitcond;
+
+					data[i].Ptr_log = &logfile;
+				}
+
+				int iter=0;
+				int rc;
+				int IslandType = 0;
+				double loweststd = 1000000;
+
+				pthread_mutex_lock (&mutex);
+				time(&start);
+				while(iter < itermax){
+					for (int i=0;i<islandsN;i++){
+						if ( !*(data[i].isActive) ){
+							//Create again the i-th thread to simulate an Island
+							IslandType = 4;
+							if (IslandType == 0){
+								data[i].randomSeed = rng.next();
 								cout << "\t\t\tPSO:\t\t omega: "<< data[i].omega  <<  "\t\teta1: " << data[i].eta1 <<  "\t\teta2: " << data[i].eta2 << "\t\tVcoeff: " << data[i].vcoeff << "\t\tGenerations: " << data[i].generations << endl;
 								rc = pthread_create(&threads[i], NULL, PSOthread, (void *)&data[i]);
 							}
+							else if (IslandType == 1){
+								data[i].randomSeed = rng.next();
+								cout << "\t\t\tMPSO:\t\t omega: "<< data[i].omega  <<  "\t\teta1: " << data[i].eta1 <<  "\t\teta2: " << data[i].eta2 << "\t\tVcoeff: " << data[i].vcoeff  << "\t\tNswamrs: " << data[i].nswarms<< "\t\tGenerations: " << data[i].generations << endl;
+								rc = pthread_create(&threads[i], NULL, MPSOthread, (void *)&data[i]);
+							}
+							else if (IslandType == 2){
+								data[i].randomSeed = rng.next();
+								cout << "\t\t\tSGA:\t\t CR: "<< data[i].CRsga  <<  "\t\tM: " << data[i].M <<  "\t\tInsertBest: " << data[i].insert_best << "\t\tGenerations: " << data[i].generations << endl;
+								rc = pthread_create(&threads[i], NULL, SGAthread, (void *)&data[i]);
+							}
+							else if (IslandType == 3){
+								data[i].randomSeed = rng.next();
+								data[i].generations=10000;
+								data[i].NP = 1;
+								cout << "\t\t\tASA:\t\t Ts: "<< data[i].Ts  <<  "\t\tTf: " << data[i].Tf << "\t\tGenerations: " << data[i].generations << endl;
+								rc = pthread_create(&threads[i], NULL, ASAthread, (void *)&data[i]);
+							}
 							else {
-								data[i].strategy = i+1;
+								data[i].randomSeed = rng.next();
+								data[i].strategy = strategy;
 								cout << "\t\t\tDE: \t\t F: "<< data[i].F  <<  "\t\tCR: " << data[i].CR << "\t\tStrategy: " << data[i].strategy << "\t\tGenerations: " << data[i].generations << endl;
 								rc = pthread_create(&threads[i], NULL, DEthread, (void *)&data[i]);
-								count2++;
 							}
-
 							if (rc){
+							//Problem creating the thread
 								printf("ERROR; return code from pthread_create() is %d\n", rc);
 								exit(-1);
 							}
 							else{
-								globalCount += 1;
+							//Thread Successfully Created
+								iter += 1;
 								*(data[i].isActive) = true;
+								//thread is detached as it will never be joined
+								//(on OSx Leopard this avoids that at a certain point
+								//threads cannot be created anymore resulting in rc=35)
+								pthread_detach(threads[i]);
+							}
+
+							//evaluate standard deviation in main population
+							//loweststd = IslandPop[0].evaluateStd();
+							loweststd = IslandPop.evaluateStd();
+
+							//log islands best and std
+							for (int i2=0;i2<islandsN;i2++){
+							cout << "CPU#:" << i2 << ": " << IslandPop.extractBestIndividual().getFitness() << " " << IslandPop.evaluateStd() << endl;
 							}
+							cout << "HighestStd: " << loweststd << endl;
+							cout << "Iter: " << iter+1 << endl;
+							cout << "\nmain():" << endl <<  "\t\t\tcreating thread, ID " << i << endl;
+
+							//ring topology migration
+							
+
+						}
 					}
+					pthread_cond_wait(&exitcond, &mutex);
 				}
-				pthread_cond_wait(&exitcond, &mutex);
+				pthread_mutex_unlock (&mutex);
+				//The main cycle has finished: we wait for all threads to finish
+				for (int i=0; i<islandsN;i++){
+					while (*data[i].isActive); //infinite loop if a thread never ends
+				}
+
+				//we clean the thread objects and deallocate memory
+				pthread_mutex_destroy(&mutex);
+				pthread_cond_destroy(&exitcond);
+				delete[] data;
+				delete[] threads;
+				delete[] isActive;
+				time(&end);
+				dif = difftime(end,start);
+				cout << "\nSeconds elapsed: " << dif << endl<<endl;
+
+				double res=IslandPop.extractBestIndividual().getFitness();
+				//if (iter<itermax){
+				results.push_back(res);
+				//}
+
+				for (int i = 0; i < islandsN; i++) {
+					delete parallelProblems[i];
+				}
+			}
+
+			//evaluate experiment results
+			for (unsigned int i=0;i<results.size();i++){
+				mean += results[i];
+				if (results[i] > max) max = results[i];
+				if (results[i] < min) min = results[i];
 			}
-			pthread_mutex_unlock (&mutex);
+			mean /= results.size();
 
-			//The main cycle has finished: we wait for all threads to finish
-			for (int i=0; i<NUM_THREADS;i++){
-				pthread_join(threads[i],NULL);
+			for (unsigned int i=0;i<results.size();i++){
+				dev += (mean-results[i])*(mean-results[i]);
 			}
+			dev = sqrt(dev/results.size());
 
-			//we clean the thread objects and deallocate memory
-			pthread_mutex_destroy(&mutex);
-			//pthread_attr_destroy(&attr);
-			pthread_cond_destroy(&exitcond);
-			delete[] data;
-			delete[] threads;
-			delete[] isActive;
-			time(&end);
-			dif = difftime(end,start);
-			cout << "\nSeconds elapsed: " << dif << endl<<endl;
-			cout << "\nCount1: " <<  count1 << endl;
-		    cout << "\nCount2: " <<  count2 << endl;
+			for (unsigned int i=0;i<results.size();i++){
+				cout << "\nTrial #" << i << ": " << results[i] << endl;
+			}
 
+			//print results
+			cout << "\nMean: " << mean << endl;
+			cout << "Std: " << dev << endl;
+			cout << "Max: " << max << endl;
+			cout << "Min: " << min << endl;
+			cout << "Success: " << results.size() << endl;
 			}
 			break;
+           
 		}//end case
 
 } //end while
 
 
-			//The following lines are commented out as they implement exactly the DiGMO random strategy in a random
-			//way. This creates a very unefficient cooperation between the algorithms as it does not have the implicit
-			//learning of the best strategy. This is the mechanism that allow, in the distributed version of the algorithm,
-			//the same individual to be evolved by different algorithms and only the most succesful one to be selected
-			//
-			//It seems that the reinsertion strategy (only improved individuals substitute their old versions) together with
-			//the contemporary evolution of the same individuals is key to the efficiency of algorithm cooperation (mimicking co-evolution)
-			//In a sequential version this can be  enforced by evolving with DE, PSO and SA and only at the end perform reinsertion
-			//as in the uncommented lines above
-
-			/*algorithmchoice = r250()%3;
-			if ( algorithmchoice == 0 ){
-				deme=pop.extractRandomDeme(20,picks);
-				//we print the best
-				cout << "Using DE" << endl;
-				cout << "Initial fitness: " << deme.extractBestIndividual().getFitness() << endl;
-				//DE.initDE(50,LB.size(),dr250(),dr250(),(int)dr250()*11+1);
-				deme = DE.evolve(deme,&rosenbrock,LB,UB);
-			}
-
-			else if ( algorithmchoice == 1 ){
-				deme=pop.extractRandomDeme(20,picks);
-				//deme.resetVelocities(LB,UB);
-			    //we print the best
-				cout << "Using PSO" << endl;
-				cout << "Initial fitness: " <<  deme.extractBestIndividual().getFitness() << endl;
-				//PSO.initPSO(500,LB.size(),(dr250()+1)/2,4.0*dr250(),4.0*dr250(),(int)dr250()*10000);
-				deme = PSO.evolve(deme,&rosenbrock,LB,UB);
-			}
-
-			else if ( algorithmchoice == 2 ){
-				deme=pop.extractRandomDeme(1,picks);
-				//deme.resetVelocities(LB,UB);
-				//we print the best
-				cout << "Using ASA" << endl;
-				cout << "Initial fitness: " << deme.extractBestIndividual().getFitness() << endl;
-				//ASA.initASA(10000,1,8,LB.size(),dr250(),dr250(),1.0);
-				deme = ASA.evolve(deme[0],&rosenbrock,LB,UB);
-			}
-
-			//we print the result
-			cout << "Final fitness: " <<  deme.extractBestIndividual().getFitness() << endl;
-			//we insert it in the main population
-			pop.insertDeme(deme,picks);
-			picks.clear();*/
+
 
 
 		cout << "\nBest fitness found: " <<  pop.extractBestIndividual().getFitness() << endl;
-- 
2.11.4.GIT