PaGMO  1.1.5
bee_colony.cpp
1 /*****************************************************************************
2  * Copyright (C) 2004-2015 The PaGMO development team, *
3  * Advanced Concepts Team (ACT), European Space Agency (ESA) *
4  * *
5  * https://github.com/esa/pagmo *
6  * *
7  * act@esa.int *
8  * *
9  * This program is free software; you can redistribute it and/or modify *
10  * it under the terms of the GNU General Public License as published by *
11  * the Free Software Foundation; either version 2 of the License, or *
12  * (at your option) any later version. *
13  * *
14  * This program is distributed in the hope that it will be useful, *
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
17  * GNU General Public License for more details. *
18  * *
19  * You should have received a copy of the GNU General Public License *
20  * along with this program; if not, write to the *
21  * Free Software Foundation, Inc., *
22  * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
23  *****************************************************************************/
24 
25 #include <string>
26 #include <vector>
27 #include <boost/random/uniform_int.hpp>
28 #include <boost/random/uniform_real.hpp>
29 
30 #include "bee_colony.h"
31 #include "../exceptions.h"
32 #include "../population.h"
33 #include "../problem/base.h"
34 #include "../types.h"
35 
36 
37 
38 
39 namespace pagmo { namespace algorithm {
40 
42 
49 bee_colony::bee_colony(int gen, int limit):base(),m_iter(gen), m_limit(limit) {
50  if (gen < 0) {
51  pagmo_throw(value_error,"number of generations must be nonnegative");
52  }
53 
54  if (limit < 0) {
55  pagmo_throw(value_error,"limit value must be nonnegative");
56  }
57 
58 }
59 
62 {
63  return base_ptr(new bee_colony(*this));
64 }
65 
67 
74 {
75  // Let's store some useful variables.
76  const problem::base &prob = pop.problem();
77  const problem::base::size_type prob_i_dimension = prob.get_i_dimension(), D = prob.get_dimension(), Dc = D - prob_i_dimension, prob_c_dimension = prob.get_c_dimension();
78  const decision_vector &lb = prob.get_lb(), &ub = prob.get_ub();
79  const population::size_type NP = (int) pop.size();
80 
81  //We perform some checks to determine whether the problem/population are suitable for ABC
82  if ( Dc == 0 ) {
83  pagmo_throw(value_error,"There is no continuous part in the problem decision vector for ABC to optimise");
84  }
85 
86  if ( prob.get_f_dimension() != 1 ) {
87  pagmo_throw(value_error,"The problem is not single objective and ABC is not suitable to solve it");
88  }
89 
90  if ( prob_c_dimension != 0 ) {
91  pagmo_throw(value_error,"The problem is not box constrained and ABC is not suitable to solve it");
92  }
93 
94  if (NP < 2) {
95  pagmo_throw(value_error,"for ABC at least 2 individuals in the population are needed");
96  }
97 
98  // Get out if there is nothing to do.
99  if (m_iter == 0) {
100  return;
101  }
102 
103  // Some vectors used during evolution are allocated here.
104  fitness_vector fnew(prob.get_f_dimension());
105  decision_vector dummy(D,0); //used for initialisation purposes
106  std::vector<decision_vector > X(NP,dummy); //set of food sources
107  std::vector<fitness_vector> fit(NP); //food sources fitness
108 
109  decision_vector temp_solution(D,0);
110 
111  std::vector<int> trial(NP,0);
112 
113  std::vector<double> probability(NP);
114 
115  population::size_type neighbour = 0;
116 
117  decision_vector::size_type param2change = 0;
118 
119  std::vector<double> selectionfitness(NP), cumsum(NP), cumsumTemp(NP);
120  std::vector <population::size_type> selection(NP);
121 
122 
123  double r = 0;
124 
125  // Copy the food sources position and their fitness
126  for ( population::size_type i = 0; i<NP; i++ ) {
127  X[i] = pop.get_individual(i).cur_x;
128  fit[i] = pop.get_individual(i).cur_f;
129  }
130 
131  // Main ABC loop
132  for (int j = 0; j < m_iter; ++j) {
133  //1- Send employed bees
134  for (population::size_type ii = 0; ii< NP; ++ii) {
135  //selects a random component (only of the continuous part) of the decision vector
136  param2change = boost::uniform_int<decision_vector::size_type>(0,Dc-1)(m_urng);
137  //randomly chose a solution to be used to produce a mutant solution of solution ii
138  //randomly selected solution must be different from ii
139  do{
140  neighbour = boost::uniform_int<population::size_type>(0,NP-1)(m_urng);
141  }
142  while(neighbour == ii);
143 
144  //copy local solution into temp_solution (the whole decision_vector, also the integer part)
145  for(population::size_type i=0; i<D; ++i) {
146  temp_solution[i] = X[ii][i];
147  }
148 
149  //mutate temp_solution
150  temp_solution[param2change] = X[ii][param2change] + boost::uniform_real<double>(-1,1)(m_drng) * (X[ii][param2change] - X[neighbour][param2change]);
151 
152  //if generated parameter value is out of boundaries, it is shifted onto the boundaries*/
153  if (temp_solution[param2change]<lb[param2change]) {
154  temp_solution[param2change] = lb[param2change];
155  }
156  if (temp_solution[param2change]>ub[param2change]) {
157  temp_solution[param2change] = ub[param2change];
158  }
159 
160  //Calling void prob.objfun(fitness_vector,decision_vector) is more efficient as no memory allocation occur
161  //A call to fitness_vector prob.objfun(decision_vector) allocates memory for the return value.
162  prob.objfun(fnew,temp_solution);
163  //If the new solution is better than the old one replace it with the mutant one and reset its trial counter
164  if(prob.compare_fitness(fnew, fit[ii])) {
165  X[ii][param2change] = temp_solution[param2change];
166  pop.set_x(ii,X[ii]);
167  prob.objfun(fit[ii], X[ii]); //update the fitness vector
168  trial[ii] = 0;
169  }
170  else {
171  trial[ii]++; //if the solution can't be improved increase its trial counter
172  }
173  } //End of loop on the population members
174 
175  //2 - Send onlooker bees
176  //We scale all fitness values from 0 (worst) to absolute value of the best fitness
177  fitness_vector worstfit=fit[0];
178  for (pagmo::population::size_type i = 1; i < NP;i++) {
179  if (prob.compare_fitness(worstfit,fit[i])) worstfit=fit[i];
180  }
181 
182  for (pagmo::population::size_type i = 0; i < NP; i++) {
183  selectionfitness[i] = fabs(worstfit[0] - fit[i][0]) + 1.;
184  }
185 
186  // We build and normalise the cumulative sum
187  cumsumTemp[0] = selectionfitness[0];
188  for (pagmo::population::size_type i = 1; i< NP; i++) {
189  cumsumTemp[i] = cumsumTemp[i - 1] + selectionfitness[i];
190  }
191  for (pagmo::population::size_type i = 0; i < NP; i++) {
192  cumsum[i] = cumsumTemp[i]/cumsumTemp[NP-1];
193  }
194 
195  for (pagmo::population::size_type i = 0; i < NP; i++) {
196  r = m_drng();
197  for (pagmo::population::size_type j = 0; j < NP; j++) {
198  if (cumsum[j] > r) {
199  selection[i]=j;
200  break;
201  }
202  }
203  }
204 
205  for(pagmo::population::size_type t = 0; t < NP; ++t) {
206  r = m_drng();
207  pagmo::population::size_type ii = selection[t];
208  //selects a random component (only of the continuous part) of the decision vector
209  param2change = boost::uniform_int<decision_vector::size_type>(0,Dc-1)(m_urng);
210  //randomly chose a solution to be used to produce a mutant solution of solution ii
211  //randomly selected solution must be different from ii
212  do{
213  neighbour = boost::uniform_int<population::size_type>(0,NP-1)(m_urng);
214  }
215  while(neighbour == ii);
216 
217  //copy local solution into temp_solution (also integer part)
218  for(population::size_type i=0; i<D; ++i) {
219  temp_solution[i] = X[ii][i];
220  }
221 
222  //mutate temp_solution
223  temp_solution[param2change] = X[ii][param2change] + boost::uniform_real<double>(-1,1)(m_drng) * (X[ii][param2change] - X[neighbour][param2change]);
224 
225  /*if generated parameter value is out of boundaries, it is shifted onto the boundaries*/
226  if (temp_solution[param2change]<lb[param2change]) {
227  temp_solution[param2change] = lb[param2change];
228  }
229  if (temp_solution[param2change]>ub[param2change]) {
230  temp_solution[param2change] = ub[param2change];
231  }
232 
233  //Calling void prob.objfun(fitness_vector,decision_vector) is more efficient as no memory allocation occur
234  //A call to fitness_vector prob.objfun(decision_vector) allocates memory for the return value.
235  prob.objfun(fnew,temp_solution);
236  //If the new solution is better than the old one replace it with the mutant one and reset its trial counter
237  if(prob.compare_fitness(fnew, fit[ii])) {
238  X[ii][param2change] = temp_solution[param2change];
239  pop.set_x(ii,X[ii]);
240  prob.objfun(fit[ii], X[ii]); //update the fitness vector
241  trial[ii] = 0;
242  }
243  else {
244  trial[ii]++; //if the solution can't be improved increase its trial counter
245  }
246  }
247 
248  //3 - Send scout bees
249  int maxtrialindex = 0;
250  for (population::size_type ii=1; ii<NP; ++ii)
251  {
252  if (trial[ii] > trial[maxtrialindex]) {
253  maxtrialindex = ii;
254  }
255  }
256  if(trial[maxtrialindex] >= m_limit)
257  {
258  //select a new random solution
259  for(problem::base::size_type jj = 0; jj < Dc; ++jj) {
260  X[maxtrialindex][jj] = boost::uniform_real<double>(lb[jj],ub[jj])(m_drng);
261  }
262  trial[maxtrialindex] = 0;
263  pop.set_x(maxtrialindex,X[maxtrialindex]);
264  }
265 
266  } // end of main ABC loop
267 
268 }
269 
271 std::string bee_colony::get_name() const
272 {
273  return "Artificial Bee Colony optimization";
274 }
275 
276 
278 
282 {
283  std::ostringstream s;
284  s << "gen:" << m_iter << ' ';
285  s << "limit:" << m_limit << ' ';
286  return s.str();
287 }
288 
289 }} //namespaces
290 
291 BOOST_CLASS_EXPORT_IMPLEMENT(pagmo::algorithm::bee_colony)
boost::shared_ptr< base > base_ptr
Alias for shared pointer to base algorithm.
bee_colony(int gen=1, int limit=20)
Constructor.
Definition: bee_colony.cpp:49
std::string get_name() const
Algorithm name.
Definition: bee_colony.cpp:271
Root PaGMO namespace.
std::vector< double > decision_vector
Decision vector type.
Definition: types.h:40
fitness_vector cur_f
Current fitness vector.
Definition: population.h:89
const individual_type & get_individual(const size_type &) const
Get constant reference to individual at position n.
Definition: population.cpp:277
Base algorithm class.
Base problem class.
Definition: problem/base.h:148
base_ptr clone() const
Clone method.
Definition: bee_colony.cpp:61
Population class.
Definition: population.h:70
size_type get_dimension() const
Return global dimension.
bool compare_fitness(const fitness_vector &, const fitness_vector &) const
Compare fitness vectors.
void evolve(population &) const
Evolve implementation.
Definition: bee_colony.cpp:73
fitness_vector objfun(const decision_vector &) const
Return fitness of pagmo::decision_vector.
size_type get_i_dimension() const
Return integer dimension.
std::vector< double > fitness_vector
Fitness vector type.
Definition: types.h:42
c_size_type get_c_dimension() const
Return global constraints dimension.
const decision_vector & get_ub() const
Upper bounds getter.
container_type::size_type size_type
Population size type.
Definition: population.h:192
The Artificial Bee Colony Solver (ABC)
Definition: bee_colony.h:57
decision_vector cur_x
Current decision vector.
Definition: population.h:83
rng_uint32 m_urng
Random number generator for unsigned integer values.
f_size_type get_f_dimension() const
Return fitness dimension.
const decision_vector & get_lb() const
Lower bounds getter.
rng_double m_drng
Random number generator for double-precision floating point values.
std::string human_readable_extra() const
Extra human readable algorithm info.
Definition: bee_colony.cpp:281
decision_vector::size_type size_type
Problem's size type: the same as pagmo::decision_vector's size type.
Definition: problem/base.h:160