粒子群算法

粒子群算法

PSO算法，是一种基于迭代的优化算法。系统初始化为一组随机解，通过迭代搜寻最优值。在粒子群算法中，每个个体称为一个粒子，粒子是待求解问题的潜在可能的解，粒子种群由若干个粒子组成。粒子群算法的基本原理是粒子种群在搜索空间以一定的速度飞行, 每个粒子在搜索时，考虑自己搜索到的历史最优位置和种群内其他粒子的历史最优位置, 在此基础上进行位置的变化。

更新速度和位置采用以下两个公式：
version[i]=w_newversion[i-1]+c1r1*(p_x_bestposition[i]-x_position[i])+c2r2(px[i]-x_position[i]);//新的速度
X[i]=x[i-1]+v[i]//新的位置

式中：w为惯性权重，非负数，调节对解空间的搜索范围，c1、c2为学习因子，也称为加速常数，r1、r2增加随机搜索性。
第一部分：惯性部分，反映了粒子的运动习惯，代表粒子有维持自己先前速度的趋势；
第二部分：自我认知，反映了粒子对自身历史经验的记忆，代表粒子有向自身最佳位置逼近的趋势；
第三部分：社会认知，反映了粒子间协同与知识共享的群体历史经验，代表粒子有向群体或领域历史最佳位置逼近的趋势。`

Optimizer.h文件
#ifndef _OPTIMIZER_H
#define _OPTIMIZER_H
#define c1max   1.5
#define c1min   0.9//c1学习因子1，自我认知
#define c2max   1.5
#define c2min   0.9 //c2学习因子2，社会认知
#define w_max  0.9 //最大惯性权重因子，影响着全局搜索
#define w_min  0.6 //最小惯性权重因子，局部深度搜索能力
#define M_particle   2000//最大迭代次数，
#define N_particle  500  //初始化群体的个体，N很小容易陷入局部优化,N很大优化能力很好，优化的寻优精度和运行时间
#define x_min -1000//变量的下限
#define x_max 1000//变量的上限
#define MINIMUM_ACCURACY 1e-20
#include <iostream>
#include<ctime>
#include<cstdlib>
#include<cmath>
using namespace std;

class optimizer
{
public:
    virtual double setOptimizer(void (*foo)(double *,double*),double *,double *,int,int)=0;//计算误差值




protected:
    double *px;//指针指向自变量
    double *py;//指针指向因变量
    int xd;//自变量个数
    int yd;//因变量个数
    void (*fun)(double *,double *);

};

class slove_x:public optimizer
{
public:
    virtual double setOptimizer(void (*foo)(double *,double*),double *,double *,int,int);//计算误差值
    slove_x();
    void Init(void (*foo)(double *,double*),double* a,double*b,double*m,int c,int d);
    void Initialize_vel_position();//给定初始化粒子群速度和位置
    void Move_fit(double * pt);//计算粒子群移动
    ~slove_x();//析构函数
    double error_value;//误差值
    double x_position[4];   //储存当期自变量的横坐标



private:

    double p_x_bestposition[4];//储存局部每代最优变量
    double version[4];   //储存速度
    double *p_current_calculate;
    int num=0;//该粒子已经迭代次数
    double r1, r2;//随机参数
    double w_new;//权重更新,初始化为最大值，利于全局搜索


};
slove_x::slove_x()
{

}
slove_x::~slove_x()
{

}



void slove_x::Init(void (*foo)(double *,double*),double* a,double*b,double*m,int c,int d)
{
    fun=foo;
    px=a;
    py=b;
    xd=c;
    yd=d;
    p_current_calculate=m;


}
void slove_x::Initialize_vel_position()
{
    for(int i=0; i<xd; i++)
    {


        double k=rand()/double(RAND_MAX);

        x_position[i]=k*20-10;//局部最优解初始化


        double k1=rand()/double(RAND_MAX);
        version[i]=(k1*20-10)/8;//局部漂移速度初始化


    }

    error_value=setOptimizer(fun,x_position,p_current_calculate,xd,yd);

}

void slove_x::Move_fit(double * pt)
{
    w_new=w_max-(w_max-w_min)*(double(num)/double(M_particle));
    double c1=(c1max-c1min)*(double(num)/double(M_particle))+c1min;
    double c2=(c2max-c2min)*(double(num)/double(M_particle))+c2min;

    r1=rand()/double(RAND_MAX);
    r2=rand()/double(RAND_MAX);

    for(int i=0; i<xd; i++)
    {

        version[i]=w_new*version[i]+c1*r1*(p_x_bestposition[i]-x_position[i])+c2*r2*(px[i]-x_position[i]);//新的速度
        x_position[i]+=version[i];
        if(x_position[i]>=x_max)
            x_position[i]=x_max;
        else if(x_position[i]<=x_min)
            x_position[i]=x_min;
    }

    if(setOptimizer(fun,x_position,p_current_calculate,xd,yd)<error_value)
    {
        error_value=setOptimizer(fun,x_position,p_current_calculate,xd,yd);

        for(int i=0; i<xd; i++)
        {
            p_x_bestposition[i]=x_position[i];
        }

    }

    num++;

}
double  slove_x::setOptimizer(void (*foo)(double *,double*),double *a,double *b,int c,int d)
{
    fun(a,b);
    double sum=0;
    for(int i=0; i<d; i++)
    {
        sum+=pow((b[i]-py[i]),2) ;

    }
    return sum;
}


#endif // _OPTIMIZER_H
main.cpp文件：
#include <iostream>
#include<ctime>
#include"optimizer.h"

void Initialization(slove_x * a,int x);//初始化函数
void search_for_answer(slove_x * a,int x);//粒子漂移
void fit_function(double *x,double *y);
using namespace std;
double error;     //误差
double x[6]= {0};  //全局最优解的自变量坐标
double y[24]= {0};  //因变量或方程的值
double calculated_value[24]= {0};
int dimension=2;//解空间维数，即未知数个数
int equation_num=1;//方程个数
int main()
{
    srand(time(NULL));
    cout<<"Please enter the function in the fit_function  "<<endl<<"Please enter the result and dimension!"<<endl;
    cout<<"Dimension:"<<endl;
    cin>>dimension;//输入要求解的自变量个数
    cout<<"The num of the equation"<<endl;
    cout<<"equation_num:"<<endl;//方程个数
    cin>>equation_num;
    for(int i=0; i<equation_num; i++)
    {
        cin>>y[i];
    }
    slove_x *elementswarm=new slove_x[N_particle];//生成多个粒子对象
    Initialization(elementswarm,N_particle);//初始化位置和速度
    search_for_answer(elementswarm,N_particle);//求解函数
    cout<<"The result is :"<<endl<<x[0]<<endl<<x[1]<<endl<<x[2]<<endl<<x[3]<<endl;
    cout<<"The error is : "<<error<<endl;

    return 0;
}
void Initialization(slove_x * a,int b)
{
    for(int i=0; i<b; i++)
    {
        a[i].Init(fit_function,x,y,calculated_value,dimension,equation_num) ;
        a[i].Initialize_vel_position();
    }







    double mini=a[0].error_value;
    int flag=0;
    for(int i=1; i<b; i++)
    {



        if(a[i].error_value<mini)
        {
            flag=i;
            mini=a[i].error_value;

        }

    }
    error=a[flag].error_value;

          for(int j=0; j<dimension; j++ )
    {
        x[j]=a[flag].x_position[j];//全局最优值
    }

}
void search_for_answer(slove_x * a,int b)
{
    int flag=0;
    for(int j=0; j<M_particle; j++)
    {


        for(int i=0; i<b; i++)
        {
            a[i].Move_fit(x);
            if(a[i].error_value<error)
            {
                flag=i;
            }
        }
        error=a[flag].error_value;
        for(int i=0; i<dimension; i++)
        {
            x[i]=a[flag].x_position[i];
        }
        if(error<MINIMUM_ACCURACY)
        {

            j=M_particle;

        }

    }

}

void fit_function(double *x,double *y)
{

                                            y[0]= -(x[0]+47)*sin(pow(fabs(x[1]+0.5*x[0]+47),0.5))-x[0]*sin(pow(fabs(x[0]-x[1]-47),0.5));
    //自定义函数

}