粒子群算法

原理：
http://www.cnblogs.com/hxsyl/p/4521778.html
http://www.cnblogs.com/maybe2030/p/5043356.html

代码参考：
http://blog.csdn.net/breezedust/article/details/12378519

#include "iostream"
#include "time.h"
using namespace std; 


/*  计算y=-x(x-2)的最大值 
 *  取值范围x--[-3,3] 
 * 
 */  
class PSOTest 
{  
public:
    int n;   //粒子个数  
    double *y;  
    double *x;  
    double *v;  
    double c1;  
    double c2;  
    double *pbest;  
    double gbest;  
    double vmax;  

    //适应函数
    void fitnessFunction()  
    {  
        for(int i=0; i<n; i++)
        {  
            y[i] = -1.0 * x[i] * (x[i] - 2.0);  
        }  
    }  

    //初始化
    void init()             
    {   
        n = 2;   //粒子个数，这里为了方便演示，我们只取两个，观察其运动方向  
        c1 = 2;  
        c2 = 2; 
        vmax = 0.1;
        x = new double[n];  
        v = new double[n];  
        y = new double[n];  
        pbest = new double[n];  
        /*** 
         * 本来是应该随机产生的，为了方便演示，我这里手动随机落两个点，分别落在最大值两边 
         */  
        x[0] = 0;  
        x[1] = 2;  
        v[0] = 0.01;  
        v[1] = 0.02;  

        fitnessFunction();  

        //初始化当前个体极值，并找到群体极值  
        for(int i=0;i<n;i++)
        {  
            pbest[i] = y[i];  
            gbest = max(gbest, pbest[i]);  
        }  
        cout << "start gbest: " << gbest << endl;  
    }  

    //粒子群算法  
    void PSO(int maxGens)
    {  
        for(int i=0; i<maxGens; i++)
        {  
            double w = 0.4;  
            int j;
            for(j=0; j<n; j++)
            {  
                //更新位置和速度  
                v[j] = w * v[j] + c1 * randval(0,1) * (pbest[j]-x[j]) + c2 * randval(0,1) * (gbest-x[j]);  
                if(v[j] > vmax) 
                    v[j] = vmax;  
                x[j] += v[j];  

                //越界判断  
                if(x[j] > 3) 
                    x[j] = 3;  
                if(x[j] < -3) 
                    x[j] = -3;  

            }  

            fitnessFunction();  

            //更新个体极值和群体极值  
            for(j=0; j<n; j++)
            {  
                pbest[j] = max(y[j], pbest[j]);  
                gbest = max(gbest, pbest[j]);
                //cout << x[j] << "  " << v[j] << endl;  
            }  
            cout << "======" << (i+1) << "======gbest:" << gbest << endl;  
        } 
    }

    double max(double a,double b)
    {  
        return a > b ? a : b;  
    }  

    //随机数产生函数 [low, high) 
    double randval(double low, double high)  
    {    
        return ((double)(rand()%RAND_MAX) / (double)RAND_MAX)*(high - low) + low;  
    }

    ~PSOTest()
    {
        delete y;  
        delete x;
        delete v;
        delete pbest;
    }
};

int main()
{  
    srand(time(NULL));
    PSOTest ts;  
    ts.init();  
    ts.PSO(100);  
    return 0;
} 

但是，如果将初始值设为
x[0] = -0.5
x[1] = 2.6
就找不到最优解了。