D:/Radagast_D/Projekt/ECF_trunk/ECF/AlgParticleSwarmOptimization.cpp

#include "AlgParticleSwarmOptimization.h"

ParticleSwarmOptimization::ParticleSwarmOptimization()
{
    // define algorithm name
    name_ = "ParticleSwarmOptimization";

    // create selection operators needed
    // in this case, SelBestOp
    selBestOp = static_cast<SelectionOperatorP> (new SelBestOp);
}


// register algorithm parameters
void ParticleSwarmOptimization::registerParameters(StateP state)
{
    registerParameter(state, "weightType", (voidP) new InertiaWeightType(CONSTANT), ECF::INT);
    registerParameter(state, "weight", (voidP) new double(0.8), ECF::DOUBLE);
    registerParameter(state, "maxVelocity", (voidP) new double(10.0), ECF::DOUBLE);
    registerParameter(state, "particleCount", (voidP) new uint(10), ECF::INT);
}


bool ParticleSwarmOptimization::initialize(StateP state)
{
    // initialize all operators
    selBestOp->initialize(state);

    voidP weightType = getParameterValue(state, "weightType");
    m_weightType = *((InertiaWeightType*) weightType.get());

    voidP weight = getParameterValue(state, "weight");
    m_weight = *((double*) weight.get());

    voidP maxV = getParameterValue(state, "maxVelocity");
    m_maxV = *((double*) maxV.get());

    // test if inertia weight type is time variant and if so, check if max iterations specified
    if(m_weightType == TIME_VARIANT) {
        if(state->getRegistry()->isModified("term.maxgen")) {
            // read maxgen parameter
            m_maxIter = *(boost::static_pointer_cast<int>( state->getRegistry()->getEntry("term.maxgen") ));
        }
        else {
            state->getLogger()->log(1, "Error: term.maxgen has to be specified in order to use time variant inertia eight in PSO algorithm");
            throw("");
        }
    }

    // algorithm accepts a single Binary Genotype (for now)
    BinaryP bin = (BinaryP) new Binary::Binary;
    if(state->getGenotypes()[0]->getName() != bin->getName()
        || state->getGenotypes().size() != 1) {
        state->getLogger()->log(1, "Error: PSO algorithm accepts only a single Binary genotype!");
        throw ("");
    }

    voidP sptr = state->getGenotypes()[0]->getParameterValue(state, "dimension");
    uint numDimension = *((uint*) sptr.get());

    BinaryP binary[4];
    for(uint iGen = 1; iGen < 4; iGen++) {

        binary[iGen] = (BinaryP) new Binary::Binary;
        state->setGenotype(binary[iGen]);

        if(iGen == 3)
            binary[iGen]->setParameterValue(state, "dimension", (voidP) new uint(1));
        else
            binary[iGen]->setParameterValue(state, "dimension", (voidP) new uint(numDimension));

        // other parameters are proprietary (ignored by the algorithm)
        binary[iGen]->setParameterValue(state, "lbound", (voidP) new double(0));
        binary[iGen]->setParameterValue(state, "ubound", (voidP) new double(1));
        binary[iGen]->setParameterValue(state, "precision", (voidP) new uint(1));
    }
    state->getLogger()->log(1, "PSO algorithm: added 3 Binary genotypes (particle velocity, best-so-far postition, best-so-far fitness value)");

    return true;
}


bool ParticleSwarmOptimization::advanceGeneration(StateP state, DemeP deme)
{
//       a) For each particle:
//          1) If the fitness value is better than the best fitness value (pBest) in  history
//          2) Set current value as the new pBest
//          End
//       b) For each particle:
//          1) Find in the particle neighborhood, the particle with the best fitness
//          2) Calculate particle velocity according to the velocity equation (1)
//          3) Apply the velocity constriction
//          4) Update particle position according to the position equation (2)
//          5) Apply the position constriction
//          6) Calculate fitness value
//          End

    // a)
    for( uint i = 0; i < deme->getSize(); i++ ) { // for each particle 
        IndividualP particle = deme->at(i);

        // the whole point of this section is to compare fitness and pbest
        BinaryP bin = boost::dynamic_pointer_cast<Binary::Binary> (particle->getGenotype(3));
        double &particlePbestFitness = bin->realValue[0];
        double fitness = particle->fitness->getValue();

        bin = boost::dynamic_pointer_cast<Binary::Binary> (particle->getGenotype(0));
        std::vector< double > &positions = bin->realValue;

        bin = boost::dynamic_pointer_cast<Binary::Binary> (particle->getGenotype(2));
        std::vector< double > &pbestx = bin->realValue;

        // set particle pbestx-es
        if( /*iter == 0 ||*/ fitness < particlePbestFitness ) { // minimize error
            particlePbestFitness = fitness;

            // set pbestx-es
            for( uint j = 0;j<pbestx.size();j++ ) {
                pbestx[j] = positions[j];
            }
        }
        // NOTE store best particle index?
    }

    // b)
    for( uint i = 0; i < deme->getSize(); i++ ) { // for each particle
        IndividualP particle = deme->at(i);

        IndividualP bestParticle = selBestOp->select( *deme );

        BinaryP bin = boost::dynamic_pointer_cast<Binary::Binary> (particle->getGenotype(0));
        std::vector< double > &positions = bin->realValue;

        bin = boost::dynamic_pointer_cast<Binary::Binary> (particle->getGenotype(1));
        std::vector< double > &velocities = bin->realValue;

        bin = boost::dynamic_pointer_cast<Binary::Binary> (particle->getGenotype(2));
        std::vector< double > &pbestx = bin->realValue;


        double weight_up;

        switch( m_weightType )
        {
            //time variant weight, linear from weight to 0.4
            case TIME_VARIANT:
            weight_up = ( m_weight - 0.4 ) * ( m_maxIter - state->getGenerationNo() ) / m_maxIter + 0.4;
            break;

            // constant inertia weight
            case CONSTANT:
            default:
            weight_up = m_weight;
            break;
        }
        // calculate particle velocity according to the velocity equation (1)
        bin = boost::dynamic_pointer_cast<Binary::Binary> (bestParticle->getGenotype(2));
        std::vector< double > &bestParticlesPbestx = bin->realValue;
        for( uint j = 0; j < velocities.size(); j++ ) {
            double velocity;

            velocity = weight_up * velocities[j] +
               2 * (rand()/(float)RAND_MAX) * (pbestx[j] - positions[j]) +
               2 * (rand()/(float)RAND_MAX) * (bestParticlesPbestx[j] - positions[j]);
            if( velocity > m_maxV ) velocity = m_maxV;
            velocities[j] = velocity;

            positions[j] += velocities[j];
            // TODO apply position constriction
        }

        // determine new particle fitness
        evaluate( particle );
    }

    return true;
}