pso.cpp

pso 粒子智能优化算法

// PSO.cpp: implementation of the CPSO class.
//
//////////////////////////////////////////////////////////////////////
//群微粒算法：本算法涌群微粒算法求目标函数的最大值
//本算法使用步骤
// （1）派生自己的群微粒类，类中必须定义double GetFit(PARTICLE&)方法，用来计算每个微粒的适合度
// （2）生成派生类实例，并在构造函数中指明微粒坐标维数和群体个数
// （2）设置微粒坐标上界数组和下界数组，并用SetXup与SetXdown设置微粒坐标上下界
// （3）用SetVmax方法设置微粒最大速度
// （4）设置可选参数：C1，C2，W和通讯函数
// （5）采用Run方法进行优化运算，优化后用GetBest方法获得最优个体适合度和坐标

#include "stdafx.h"
#include "SBGA.h"
#include "PSO.h"


#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
#define new DEBUG_NEW
#endif

//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

//PSO构造函数
CPSO::CPSO()
{
	Particle = 0;
	PNum = 0;
	GBestIndex = 0;
	Xup = 0;
	Xdown = 0;
	W = 1;
	C1 = 2;
	C2 = 2;
	Com = 0;
}
CPSO::CPSO(int dim, int num)
{
	Particle = new PARTICLE[num];
	for(int i=0; i<num; i++) Particle[i].SetDim(dim);
	PNum = num;
	GBestIndex = 0;
	Xup = new double[dim];
	Xdown = new double[dim];
	Vmax = new double[dim];
	W = 1;
	C1 = 2;
	C2 = 2;
	Com = 0;
}

//析构函数
CPSO::~CPSO()
{
    if(Particle) delete []Particle;
	if(Xup) delete []Xup;
	if(Xdown) delete []Xdown;
	if(Vmax) delete []Vmax;
}
//设置坐标上界
void CPSO::SetXup(double *up)
{
	if(!Particle) return;
	for(int i=0; i<Particle[0].Dim; i++)
		Xup[i] = up[i];
}
//设置坐标下界
void CPSO::SetXdown(double *d)
{
	if(!Particle) return;
	for(int i=0; i<Particle[0].Dim;i++)
		Xdown[i] = d[i];
}
//设置最大速度
void CPSO::SetVmax(double *max)
{
	if(!Particle) return;
	for(int i=0; i<Particle[0].Dim;i++)
		Vmax[i] = max[i];
}
void CPSO::SetVmax(double p)
{
	if(!Particle) return;
	for(int i=0; i<Particle[0].Dim; i++)
		Vmax[i] = (Xup[i]-Xdown[i])*p;
}
//初始化群体
void CPSO::Initialize()
{
	if(!Particle) return;
	static int kk=(unsigned)time(NULL);
	srand((unsigned)time(NULL)+kk++);
	GBestIndex = 0;
	for(int i=0; i<PNum; i++)
	{
		for(int j=0; j<Particle[i].Dim; j++)
		{
			Particle[i].X[j] = rand()/(double)RAND_MAX*(Xup[j]-Xdown[j])+Xdown[j];//初始化坐标
			Particle[i].XBest[j] = Particle[i].X[j];
			Particle[i].V[j] = rand()/(double)RAND_MAX*Vmax[j]-Vmax[j]/2;//初始化速度
		}
		Particle[i].Fit = GetFit(Particle[i]); //计算该微粒适合度
		Particle[i].FitBest = Particle[i].Fit; //设最优适合度初值
		if(Particle[i].Fit<Particle[GBestIndex].Fit) GBestIndex = i;//查找群体最优微粒
	}
}

//计算群体各个微粒适合度
void CPSO::CalFit()
{
	if(!Particle) return;
	for(int i=0; i<PNum; i++)
		Particle[i].Fit = GetFit(Particle[i]);
}

//微粒飞翔,产生新一代微粒
void CPSO::ParticleFly()
{
	static double FitBak[100];
	if(!Particle) return;
	static int tt=(unsigned)time(NULL);
	srand((unsigned)time(NULL)+tt++);
	//整个群体飞向新的位置
	for(int i=0; i<PNum; i++)
	{
		for(int j=0; j<Particle[i].Dim; j++)
			Particle[i].V[j] = W*Particle[i].V[j]+// (FitBak[i]-Particle[i].Fit)+//修改速度
			rand()/(double)RAND_MAX*C1*(Particle[i].XBest[j]-Particle[i].X[j])+
			rand()/(double)RAND_MAX*C2*(Particle[GBestIndex].XBest[j]-Particle[i].X[j]);
		for(j=0; j<Particle[i].Dim; j++) //检查速度最大值
		{
			if(Particle[i].V[j]>Vmax[j]) Particle[i].V[j] = Vmax[j];
			if(Particle[i].V[j]<-Vmax[j]) Particle[i].V[j] = -Vmax[j];
		}
		for(j=0; j<Particle[i].Dim; j++)
		{
			Particle[i].X[j] += Particle[i].V[j]; //修改坐标
			if(Particle[i].X[j]>Xup[j]) Particle[i].X[j]=Xup[j];//保护
			if(Particle[i].X[j]<Xdown[j]) Particle[i].X[j]=Xdown[j];
		}
	}
	//计算各微粒适合度
	CalFit();
	for(i=0; i<PNum; i++) FitBak[i] = Particle[i].Fit;
	//设置新的个体最好位置
	for(i=0; i<PNum; i++)
		if(Particle[i].Fit<=Particle[i].FitBest)
		{
			Particle[i].FitBest = Particle[i].Fit;
			for(int j=0; j<Particle[i].Dim; j++)
				Particle[i].XBest[j] = Particle[i].X[j];
		}
		//设置新的最优个体
		GBestIndex = 0;
		for(i=0; i<PNum; i++)
			if(Particle[i].FitBest<=Particle[GBestIndex].FitBest && i!=GBestIndex) GBestIndex = i;
}

//运行群粒算法
PARTICLE& CPSO::Run(int n)
{
	Initialize();
	double *opt_p = new double[Particle[0].Dim]; //通讯用数组，最优点坐标
	double **opt_a = new double*[PNum]; //通讯用数组，所有点坐标
	for(int i=0; i<n; i++)
	{
		ParticleFly();
		if(Com) //通讯函数存在，完成通讯
		{
			for(int k=0; k<Particle[0].Dim; k++) opt_p[k]=Particle[GBestIndex].XBest[k];//拷贝最优点坐标
			for(k=0; k<PNum; k++) opt_a[k] = Particle[k].X; //指向所有点坐标
			if(!Com(Particle[GBestIndex].FitBest,opt_p,opt_a,GBestIndex)) break;
		}
	}
	delete []opt_p;
	delete []opt_a;
	return Particle[GBestIndex];
}
PARTICLE& CPSO::Run(double fit)
{
	double *opt_p = new double[Particle[0].Dim]; //通讯用数组，最优点坐标
	double **opt_a = new double*[PNum]; //通讯用数组，所有点坐标
	Initialize();
	do
	{
		ParticleFly();
		if(Com) //通讯函数存在，完成通讯
		{
			for(int k=0; k<Particle[0].Dim; k++) opt_p[k]=Particle[GBestIndex].XBest[k];//拷贝最优点坐标
			for(k=0; k<PNum; k++) opt_a[k] = Particle[k].X; //指向所有点坐标
			if(!Com(Particle[GBestIndex].FitBest,opt_p,opt_a,GBestIndex)) break;
		}
	}while(Particle[GBestIndex].FitBest<fit);
	delete []opt_p;
	delete []opt_a;
	return Particle[GBestIndex];
}

//返回最佳个体
double CPSO::GetBest(double *r)
{
	for(int i=0; i<Particle[GBestIndex].Dim; i++)
		r[i] = Particle[GBestIndex].XBest[i];
	return Particle[GBestIndex].FitBest;
}