25 #include <boost/random/uniform_int.hpp>
26 #include <boost/random/uniform_real.hpp>
31 #include "pso_generational.h"
32 #include "../problem/base_stochastic.h"
36 namespace pagmo {
namespace algorithm {
58 pso_generational::pso_generational(
int gen,
double omega,
double eta1,
double eta2,
double vcoeff,
int variant,
int neighb_type,
int neighb_param):
base(),m_gen(gen),m_omega(omega),m_eta1(eta1),m_eta2(eta2),m_vcoeff(vcoeff),m_variant(variant),m_neighb_type(neighb_type),m_neighb_param(neighb_param) {
60 pagmo_throw(value_error,
"number of generations must be nonnegative");
63 if (m_omega < 0 || m_omega > 1) {
66 pagmo_throw(value_error,
"the particles' inertia must be in the [0,1] range");
69 pagmo_throw(value_error,
"the particles' constriction coefficient must be in the [0,1] range");
72 if (eta1 < 0 || eta2 < 0 || eta1 > 4 || eta2 > 4) {
73 pagmo_throw(value_error,
"the eta parameters must be in the [0,4] range");
76 if (vcoeff <= 0 || vcoeff > 1) {
77 pagmo_throw(value_error,
"fraction of variables' range in which velocity may vary should be in ]0,1]");
80 if (variant < 1 || variant > 6) {
81 pagmo_throw(value_error,
"algorithm variant must be one of 1 ... 6");
84 if (neighb_type < 1 || neighb_type > 4) {
85 pagmo_throw(value_error,
"swarm topology variant must be one of 1 ... 4");
116 pagmo_throw(value_error,
"There is no continuous part in the problem decision vector for PSO to optimise");
119 if( prob_c_dimension != 0 ){
120 pagmo_throw(value_error,
"The problem is not box constrained and PSO is not suitable to solve it");
123 if( prob_f_dimension != 1 ){
124 pagmo_throw(value_error,
"The problem is not single objective and PSO is not suitable to solve it");
128 if (swarm_size == 0 || m_gen == 0) {
133 std::vector<double> dummy(Dc,0);
135 std::vector<decision_vector> X(swarm_size,dummy);
136 std::vector<fitness_vector> fit(swarm_size);
138 std::vector<decision_vector > V(swarm_size,dummy);
140 std::vector<decision_vector> lbX(swarm_size,dummy);
141 std::vector<fitness_vector> lbfit(swarm_size);
144 std::vector< std::vector<int> > neighb(swarm_size);
148 bool best_fit_improved;
162 for( d = 0; d < Dc; d++ ){
163 vwidth = ( ub[d] - lb[d] ) * m_vcoeff;
164 minv[d] = -1.0 * vwidth;
169 for( p = 0; p < swarm_size; p++ ){
176 for( p = 0; p < swarm_size; p++ ){
182 switch( m_neighb_type ){
183 case 1: initialize_topology__gbest( pop, best_neighb, best_fit, neighb );
break;
184 case 3: initialize_topology__von( neighb );
break;
185 case 4: initialize_topology__adaptive_random( neighb );
186 best_fit = pop.champion().
f;
189 default: initialize_topology__lbest( neighb );
194 double acceleration_coefficient = m_eta1 + m_eta2;
203 for(
int g = 0; g < m_gen; ++g ){
206 for( p = 0; p < swarm_size; p++ ){
211 if( m_neighb_type != 1 && m_variant != 6)
212 best_neighb = particle__get_best_neighbor( p, neighb, lbX, lbfit, prob );
217 if( m_variant == 1 ){
218 for( d = 0; d < Dc; d++ ){
221 V[p][d] = m_omega * V[p][d] + m_eta1 * r1 * (lbX[p][d] - X[p][d]) + m_eta2 * r2 * (best_neighb[d] - X[p][d]);
228 else if( m_variant == 2 ){
229 for( d = 0; d < Dc; d++ ){
231 V[p][d] = m_omega * V[p][d] + m_eta1 * r1 * (lbX[p][d] - X[p][d]) + m_eta2 * r1 * (best_neighb[d] - X[p][d]);
237 else if( m_variant == 3 ){
240 for( d = 0; d < Dc; d++ ){
241 V[p][d] = m_omega * V[p][d] + m_eta1 * r1 * (lbX[p][d] - X[p][d]) + m_eta2 * r2 * (best_neighb[d] - X[p][d]);
248 else if( m_variant == 4 ){
250 for( d = 0; d < Dc; d++ ){
251 V[p][d] = m_omega * V[p][d] + m_eta1 * r1 * (lbX[p][d] - X[p][d]) + m_eta2 * r1 * (best_neighb[d] - X[p][d]);
268 else if( m_variant == 5 ){
269 for( d = 0; d < Dc; d++ ){
272 V[p][d] = m_omega * ( V[p][d] + m_eta1 * r1 * (lbX[p][d] - X[p][d]) + m_eta2 * r2 * (best_neighb[d] - X[p][d]) );
286 else if( m_variant == 6 ){
287 for( d = 0; d < Dc; d++ ){
289 for( n = 0; n < neighb[p].size(); n++ )
290 sum_forces +=
m_drng() * acceleration_coefficient * ( lbX[ neighb[p][n] ][d] - X[p][d] );
292 V[p][d] = m_omega * ( V[p][d] + sum_forces / neighb[p].size() );
298 for( p = 0; p < swarm_size; p++ ){
301 for( d = 0; d < Dc; d++ ){
303 if( V[p][d] > maxv[d] )
306 else if( V[p][d] < minv[d] )
310 new_x = X[p][d] + V[p][d];
322 else if( new_x > ub[d] ){
340 for( p = 0; p < swarm_size; p++ ){
342 prob.
objfun( fit[p], X[p] );
344 prob.
objfun( lbfit[p], lbX[p] );
346 pop.push_back(lbX[p]);
353 for( p = 1; p < swarm_size; p++ ){
361 catch (
const std::bad_cast& e)
364 for( p = 0; p < swarm_size; p++ ){
366 prob.
objfun( fit[p], X[p] );
375 best_fit_improved =
false;
376 for( p = 0; p < swarm_size; p++ ){
384 if( ( m_neighb_type == 1 || m_neighb_type == 4 ) && prob.
compare_fitness( fit[p], best_fit ) ){
387 best_fit_improved =
true;
393 if( m_neighb_type == 4 && !best_fit_improved )
395 initialize_topology__adaptive_random( neighb );
411 decision_vector pso_generational::particle__get_best_neighbor(
population::size_type pidx, std::vector< std::vector<int> > &neighb,
const std::vector<decision_vector> &lbX,
const std::vector<fitness_vector> &lbfit,
const problem::base &prob )
const
415 switch( m_neighb_type ){
418 pagmo_throw(value_error,
"particle__get_best_neighbor() invoked while using a gbest swarm topology");
425 bnidx = neighb[pidx][0];
426 for( nidx = 1; nidx < neighb[pidx].size(); nidx++ )
427 if( prob.
compare_fitness( lbfit[ neighb[pidx][nidx] ], lbfit[ bnidx ] ) )
428 bnidx = neighb[pidx][nidx];
453 void pso_generational::initialize_topology__gbest(
const population &pop,
decision_vector &gbX,
fitness_vector &gbfit, std::vector< std::vector<int> > &neighb )
const
457 gbX = pop.champion().x;
458 gbfit = pop.champion().f;
464 if( m_variant == 6 ){
466 for( i = 0; i < neighb.size(); i++ )
467 neighb[0].push_back( i );
468 for( i = 1; i < neighb.size(); i++ )
469 neighb[i] = neighb[0];
486 void pso_generational::initialize_topology__lbest( std::vector< std::vector<int> > &neighb )
const
488 int swarm_size = neighb.size();
492 int radius = m_neighb_param / 2;
494 for( pidx = 0; pidx < swarm_size; pidx++ ){
495 for( j = -radius; j <= radius; j++ ){
497 nidx = (pidx + j) % swarm_size;
498 if( nidx < 0 ) nidx = swarm_size + nidx;
499 neighb[pidx].push_back( nidx );
533 void pso_generational::initialize_topology__von( std::vector< std::vector<int> > &neighb )
const
535 int swarm_size = neighb.size();
541 rows = std::sqrt( swarm_size );
542 while( swarm_size % rows != 0 )
544 cols = swarm_size / rows;
546 for( pidx = 0; pidx < swarm_size; pidx++ ){
550 for( nidx = 0; nidx < 4; nidx++ ){
551 n_x = ( p_x + vonNeumann_neighb_diff[nidx][0] ) % cols;
if( n_x < 0 ) n_x = cols + n_x;
552 n_y = ( p_y + vonNeumann_neighb_diff[nidx][1] ) % rows;
if( n_y < 0 ) n_y = rows + n_y;
554 neighb[pidx].push_back( n_y * cols + n_x );
582 void pso_generational::initialize_topology__adaptive_random( std::vector< std::vector<int> > &neighb )
const
584 int swarm_size = neighb.size();
589 for( pidx = 0; pidx < swarm_size; pidx++ )
590 neighb[pidx].clear();
593 for( pidx = 0; pidx < swarm_size; pidx++ ){
596 neighb[pidx].push_back( pidx );
598 for( j = 1; j < m_neighb_param; j++ ){
599 nidx =
m_drng() * swarm_size;
600 neighb[nidx].push_back( pidx );
611 return "PSO - Generational";
621 std::ostringstream s;
622 s <<
"gen:" << m_gen <<
' ';
623 s <<
"omega:" << m_omega <<
' ';
624 s <<
"eta1:" << m_eta1 <<
' ';
625 s <<
"eta2:" << m_eta2 <<
' ';
626 s <<
"variant:" << m_variant <<
' ';
627 s <<
"topology:" << m_neighb_type <<
' ';
628 if( m_neighb_type == 2 || m_neighb_type == 4 )
629 s <<
"topology param.:" << m_neighb_param <<
' ';
boost::shared_ptr< base > base_ptr
Alias for shared pointer to base algorithm.
std::vector< double > decision_vector
Decision vector type.
fitness_vector cur_f
Current fitness vector.
Particle Swarm optimization generational.
const individual_type & get_individual(const size_type &) const
Get constant reference to individual at position n.
std::string human_readable_extra() const
Extra human readable algorithm info.
const int vonNeumann_neighb_diff[4][2]
Von Neumann neighborhood (increments on particles' lattice coordinates that produce the coordinates o...
pso_generational(int gen=1, double omega=0.7298, double eta1=2.05, double eta2=2.05, double vcoeff=0.5, int variant=5, int neighb_type=2, int neighb_param=4)
Constructor.
size_type get_dimension() const
Return global dimension.
bool compare_fitness(const fitness_vector &, const fitness_vector &) const
Compare fitness vectors.
fitness_vector f
Fitness vector.
std::string get_name() const
Algorithm name.
fitness_vector objfun(const decision_vector &) const
Return fitness of pagmo::decision_vector.
size_type get_i_dimension() const
Return integer dimension.
fitness_vector best_f
Best fitness vector so far.
std::vector< double > fitness_vector
Fitness vector type.
c_size_type get_c_dimension() const
Return global constraints dimension.
decision_vector cur_v
Current velocity vector.
const decision_vector & get_ub() const
Upper bounds getter.
Base Stochastic Optimization Problem.
container_type::size_type size_type
Population size type.
decision_vector cur_x
Current decision vector.
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.
base_ptr clone() const
Clone method.
rng_double m_drng
Random number generator for double-precision floating point values.
decision_vector best_x
Best decision vector so far.
decision_vector::size_type size_type
Problem's size type: the same as pagmo::decision_vector's size type.
void evolve(population &) const
Evolve implementation.