GOclasses/algorithms/MPSO.cpp

   1 /*
   2  *  MPSO.cpp
   3  *  PaGMO
   4  *
   5  *  Created by Dario Izzo on 10/23/08.
   6  *  Copyright 2008 __MyCompanyName__. All rights reserved.
   7  *
   8  */
   9
  10 #include <vector>
  11
  12 #include "MPSO.h"
  13
  14 using namespace std;
  15
  16 void MPSOalgorithm::initMPSO(int generationsInit, int SolDimInit, double omegaInit, double eta1Init, double eta2Init,double vcoeffInit, int nswarmsInit, uint32_t randomSeed){
  17         generations = generationsInit;
  18         SolDim = SolDimInit;
  19         omega = omegaInit;
  20         eta1 = eta1Init;
  21         eta2 = eta2Init;
  22         vcoeff = vcoeffInit;
  23         nswarms = nswarmsInit;
  24         drng.seed(randomSeed);
  25 }
  26
  27 Population MPSOalgorithm::evolve(Population deme, GOProblem& problem){
  28
  29     const std::vector<double>& LB = problem.getLB();
  30     const std::vector<double>& UB = problem.getUB();
  31
  32         int NP = deme.size()/nswarms;                                                   //potentially dangerous when the deme size is not divisible by the numberof swamrs
  33         int D = LB.size();
  34
  35         vector<double> dummy(D,0);                                                              //used for initialisation purposes
  36         vector<double> dummy2(NP,0);                                                    //used for initialisation purposes
  37         vector< vector<double> > dummy3(NP,dummy);                              //used for initialisation purposes
  38
  39         //Particle position: X[i][j][k] = i-th swarm, j-th individual, k-th phenotype
  40         vector< vector< vector<double> > > X(nswarms,dummy3);
  41         //Particle velocity: V[i][j][k] = i-th swarm, j-th individual, k-th phenotype
  42         vector< vector< vector<double> > > V(nswarms,dummy3);
  43         //Particle fitness: fit[i][j] = i-th swarm, j-th individual
  44         vector< vector<double> > fit(nswarms,dummy2);
  45         //Global best swarm fitness: gbfit[i] = i-th swarm
  46         vector<double> gbfit(nswarms,0);
  47         //Global best swarm individual: gbX[i][j] = i-th swarm, j-th phenotype
  48         vector< vector<double> > gbX(nswarms,dummy);
  49
  50         //Local best fitness: lbfit[i][j] = i-th swarm, j-th individual
  51         vector< vector <double> > lbfit(nswarms,dummy2);
  52         //Local best chromosome: lbX[i][j][k] =  i-th swarm, j-th individual, k-th phenotype
  53         vector< vector< vector<double> > > lbX(nswarms,dummy3);
  54
  55         double vwidth;                                                                          //Width of the search space
  56     vector<double> MINV(D),MAXV(D);                                             //Maximum and minumum velocity allowed
  57
  58         // Initialise the particles (class Individual) positions, their velocities and their fitness to that of the deme
  59         for ( int i = 0; i<nswarms; i++){
  60                 for ( int j = 0; j<NP; j++ ){
  61                                 X[i][j] =       deme[NP*i+j].getDecisionVector();
  62                                 V[i][j] =       deme[NP*i+j].getVelocity();
  63                                 fit[i][j]       =       deme[NP*i+j].getFitness();
  64                 }
  65         }
  66    // Initialise the minimum and maximum velocity
  67    for ( int i = 0; i<D; i++ ) {
  68                 vwidth = (UB[i]-LB[i])/vcoeff;
  69                 MINV[i] = -1.0*vwidth;
  70                 MAXV[i] = vwidth;
  71         }
  72
  73    for (int i=0;i<nswarms;i++){
  74         // Initialise the global and local bests
  75         gbX[i]=X[i][0];
  76         gbfit[i]=fit[i][0];
  77
  78         lbX[i]=X[i];                    //at the first generation the local best position is the particle position
  79         lbfit[i]=fit[i];                //same for the fitness
  80
  81         for (int j = 1; j<NP; j++){             //the int j = 1 jumps the first member as it is already set as the best
  82                 if (fit[i][j] < gbfit[i]){
  83                         gbfit[i] = fit[i][j];
  84                         gbX[i] = X[i][j];
  85                 }
  86         }
  87    }
  88
  89    // Main PSO loop
  90    for (int iter = 0; iter<generations; iter++){
  91         //loop through the swarms
  92         for (int i = 0;i<nswarms;i++){
  93                 //1 - move the particles and check that velocity and positions are in allowed range
  94                 for (int j = 0; j< NP; j++){
  95                         for (int k = 0; k< D; k++){
  96
  97                                 //new velocity
  98                                 V[i][j][k] = omega * V[i][j][k] + eta1 * drng() * (lbX[i][j][k] - X[i][j][k]) + eta2 * drng() * (gbX[i][k] - X[i][j][k]);
  99
 100                                 //check that it is within the allowed velocity range
 101                                 if ( V[i][j][k] > MAXV[k] )
 102                                         V[i][j][k] = MAXV[k];
 103
 104                                 else if ( V[i][j][k] < MINV[k] )
 105                                         V[i][j][k] = MINV[k];
 106
 107                                 //update position
 108                                 X[i][j][k] = X[i][j][k] + V[i][j][k];
 109
 110                                 if (X[i][j][k] < LB[k])
 111                                         X[i][j][k] = drng() * (UB[k] - LB[k]) + LB[k];
 112                                 else if (X[i][j][k] > UB[k])
 113                                         X[i][j][k] = drng() * (UB[k] - LB[k]) + LB[k];
 114                         }
 115
 116                         //We evaluate the new individual fitness now as to be able to update immediately the global best
 117                         //in case a better solution is found
 118                         fit[i][j] = problem.objfun(X[i][j]);
 119                         //update local and global best
 120                         if (fit[i][j] < lbfit[i][j]){
 121                                 lbfit[i][j] = fit[i][j];        //local best
 122                                 lbX[i][j] = X[i][j];
 123                                 if(fit[i][j] < gbfit[i]){
 124                                         gbfit[i] = fit[i][j];   //global best
 125                                         gbX[i]  = X[i][j];
 126                                 }
 127                         }
 128                 } //End of loop on the population members
 129         } //End of loop on the swarms
 130
 131         //exchanges two random elements from randomly selected swarms
 132         if (iter % (int)5 == 0){
 133                 int sw1 = (int)(drng()*nswarms);                //select 1st swarm
 134                 int sw2 = (int)(drng()*nswarms);
 135                 do{                                                                                     //endless loop if nswarms<2
 136                   sw2 = (int)(drng()*nswarms);          //selects 2nd swarm different from the first
 137                 } while (sw2 == sw1);
 138
 139                 int in1 = (int)(drng()*NP);             //select 1st individual
 140                 int in2;
 141                 do{                                                                                     //endless loop if nswarms<2
 142                   in2 = (int)(drng()*NP);                   //selects 2nd individual
 143                 } while (in2 == in1);
 144                 //swap position
 145                 dummy = X[sw1][in1];
 146                 X[sw1][in1]=X[sw2][in2];
 147                 X[sw2][in2] = dummy;
 148
 149                 //swap velocity
 150                 dummy = V[sw1][in1];
 151                 V[sw1][in1]=V[sw2][in2];
 152                 V[sw2][in2] = dummy;
 153
 154                 //swap local bests lbX
 155                 dummy = lbX[sw1][in1];
 156                 lbX[sw1][in1]=lbX[sw2][in2];
 157                 lbX[sw2][in2] = dummy;
 158
 159                 //swap local best lbfit;
 160                 dummy[0] = lbfit[sw1][in1];
 161                 lbfit[sw1][in1]=lbfit[sw2][in2];
 162                 lbfit[sw2][in2] = dummy[0];
 163         }//end of swap
 164    } // end of main PSO loop
 165
 166    //we end by constructing the object Population containing the final results
 167    Population popout;
 168    Individual temp;
 169    for (int i=0; i<nswarms;i++){
 170    for (int j=0; j<NP; j++){
 171         temp.setDecisionVector(lbX[i][j]);
 172         temp.setFitness(lbfit[i][j]);
 173         temp.setVelocity(V[i][j]);
 174         popout.addIndividual(temp);
 175    }
 176    }
 177    return popout;
 178 }
 179