pso.cpp

标准PSO算法

#include"stdio.h"
#include"stdlib.h"
#include"time.h"
#include"math.h"

//定义随机数
#define rdint(i) (rand()%(int)(i))  // 随机产生小于i的整数
#define rdft() (float)((double)rdint(16384)/(16383.0))  //随机产生小于1的浮点型数据 
#define rnd(a,b) (rdint((int)(b)-(int)(a)+1)+(int)(a))  //[a]+小于b-a+1的随机数

//宏定义
#define POPSIZE 20
#define DIMENSION 2  //2维空间 x y 轴

//全局变量宏定义
float W = 1.0 ;   //速度惯性权重
float C1 = 1.8 ;   
float C2 = 1.8 ;
float VMAX = 2 ;
float XMIN = -2.0;
float XMAX = 2.0;
float P[DIMENSION];
float PBEST;

struct indi //个体微粒属性结构
{
	float number[DIMENSION]; //微粒i的第j维位置（自变量的个数）
	float best[DIMENSION];   //微粒i的第j维的最好位置
	float bestfitness;       //微粒i历史最好适应值
	float fitness;           //微粒i当前适应值
	float speed[DIMENSION];  //微粒i的j维速度

}individual[POPSIZE];        //20个微粒结构体

void initiate(void);
void calculation(int number);
void globalbest(int number);
void localbest(int number);

//程序初始化定义
void initiate()
{
//	printf("开始初始化");
	int i,j;
	//所有微粒各维位置初始化
	for(i=0;i<POPSIZE;i++)
		for(j=0;j<DIMENSION; j++)
			individual[i].number[j] = rdft()*(XMAX - XMIN) + XMIN; 

		//所有微粒各维速度初始化
	for(i=0;i<POPSIZE;i++)
		for(j=0;j<DIMENSION;j++)
			individual[i].speed[j] = VMAX*rdft();

		//所有微粒各维历史最好位置初始化
	for(i=0; i<POPSIZE; i++)
		for(j=0;j<DIMENSION;j++)
			individual[i].best[j] = individual[i].number[j];

		//计算各微粒当前适应值（即函数值）fitness
	for(i=0;i<POPSIZE; i++)
		calculation(i);
	
	//初始化最好适应值为当前适应值
	for(i=0;i<POPSIZE;i++)
		individual[i].bestfitness = individual[i].fitness;
//	printf("\n这里也执行了");

	//
	globalbest(0);
		
//	printf("初始化结束");
}

//微粒群历史最优位置修改程序
void localbest(int number)
{
	int i;
	if(individual[number].bestfitness > individual[number].fitness)
	{
		for(i=0;i<DIMENSION;i++)
			individual[number].best[i] = individual[number].number[i];
	}
	individual[number].bestfitness=individual[number].fitness;

}

//种群历史最优位置修改程序
void globalbest(int number)
{
	int i,j;
	float s=0;
	int flag = 0;
	if(number ==0 )
	{
		s=individual[0].fitness;
		flag=0;

		for(i=1;i<POPSIZE; i++)
		{
			if(individual[i].fitness < s)
			{
				s=individual[i].fitness;
				flag=i;
			}
		}
		for (i=0;i<DIMENSION;i++)
		{
			P[i] = individual[flag].number[i];
		}
		
		PBEST=individual[flag].fitness;
		

	}
	else 
	{
		for (i=0;i<POPSIZE;i++)
		{
			if (individual[i].bestfitness<PBEST)
			{
				for (j=0;j<DIMENSION;j++)
				{
					P[j]=individual[i].best[j];
				}
				PBEST=individual[i].bestfitness;
			}
		}
	}

	
}


//适应值函数计算程序
void calculation(int num)
{
	//printf("\n计算适应值开始");
	int i;
	float s=0.0,h=0.0;
	s=pow(individual[num].number[0]+individual[num].number[1]+1.0,2.0)
		*(19.0-14.0*individual[num].number[0]
		+3.0*pow(individual[num].number[0],2.0)
		-14.0*individual[num].number[1]
		+6.0*individual[num].number[0]*individual[num].number[1]
		+3.0*pow(individual[num].number[1],2.0))+1.0;

	h=30.0+pow(2.0*individual[num].number[0]-3.0*individual[num].number[1],2.0)
		*(18.0-32.0*individual[num].number[0]
		+12.0*pow(individual[num].number[0],2.0)
		+48.0*individual[num].number[1]
		-36.0*individual[num].number[0]*individual[num].number[1]
		+27.0*pow(individual[num].number[1],2.0));

	individual[num].fitness=s*h;
	//printf("计算适应值结束");
}


//主程序
void main()
{
	
		
	int i,j,k,t,total=0;
	float sum=0;
	for (j=0;j<50;j++)//50个子域
	{
        
		initiate();
		for (i=0;i<500;i++) //500次迭代进化
		{
			W=1.0-i*0.6/499; //权值逐渐减小以保证收敛于全局
			if ((PBEST-3.0)<0.0001)  //误差足够小的时候退出
			{
				total++;
				sum=sum+i;
				break;
			} 
			else
			{
				
				for (k=0;k<POPSIZE;k++)  //群体规模中每个粒子都进行进化
				{
					for (t=0;t<DIMENSION;t++)
					{
						individual[k].speed[t]=W*individual[k].speed[t]
							+C1*rdft()*(individual[k].best[t]-individual[k].number[t] )
							+C2*rdft()*(P[t]-individual[k].number[t] );
						if(individual[k].speed[t]>VMAX )
						{
							individual[k].speed[t]=VMAX;
						}

						individual[k].number[t]=individual[k].number[t]+individual[k].speed[t];

						if (individual[k].number[t]<XMIN )//为什么这样处理？
						{
							individual[k].number[t]=2*XMIN - individual[k].number[t];

						}
						if (individual[k].number[t]>XMAX  )
						{
							individual[k].number[t]=2*XMAX-individual[k].number[t];
						}
					}
					
					calculation(k);  //计算微粒k本次适应值
					
					localbest(k);  //计算微粒k历史最佳适应值
				//	printf("程序执行到这里");
					
				}
				
				globalbest(1);  //计算微粒群最佳适应值
			}
		}
		printf("%d,%f \t",i,PBEST);
	}

	printf("total number is %d \n",total);
	printf("sum is %f\n",sum);
	printf("Average is %f\n",sum/total);

	return;
}