现在的位置: 首页 > 综合 > 正文

C语言实现标准PSO算法

2013年08月08日 ⁄ 综合 ⁄ 共 4789字 ⁄ 字号 小 中 大 ⁄ 评论关闭

简介以及基本PSO的实现：http://blog.csdn.net/fovwin/article/category/1256709

相对于基本PSO，标准PSO加入了惯性权重系数W，W代表者粒子群对全局空间的搜索能力和局部收敛速度的权衡。

也就是说，若惯性权重W越大，速度更新式子的第一项所占的比重，即速度惯性项，比较大。粒子就会比较不受约束，可以冲来冲去。不受世俗的约束，也就不容易变“俗”（陷入局部极小点出不来）；若惯性权重W越小，速度更新式子的后两项所占的比重，即个体历史极值和全局最优解，比较大。粒子被自己的经验和社会的世俗所约束了，听“妈妈”的话，当然人生走的顺畅，但是你一般在你20岁就知道你后辈子的路，即容易陷入局部极小点出不来——过早收敛早熟。

所以一般在一开始一般设置较大的W值，让它充分搜索全局空间（就像小时候不听话，有时候老爸老妈说了这样不行，我们也会自己试了再说，呵呵），到粒子迭代的后半阶段，设置W较小，可以提高收敛速度（后半辈子大部分按经验活着）。

本文的W按照如下设置：

/*惯性权重系数*/
#define W_START 1.4
#define W_END	0.4
#define W_FUN	(W_START-(W_START-W_END)*pow((double)cur_n/ITE_N,2))

当然，我们可以有的别的设置公式，但是大体呈现下降的趋势，也可以用模糊控制规则（下下步实现目标）来设置～～～

至于作用效果么，比基本的能更快早到最优解～～～不过木有pp，1.X版再说～～～

相对与Verson 0.0版减少好几个函数，第一版初级阶段，主要为了学习原理，Version 1.0版代码更加简洁，将pso.c和pso.h独立出来，好以后扩展，通过define来设置你的待优化函数（下一步尝试多种诡异的函数，有的还蛮漂亮的），但是还要完善，但是加入了扩展适应度函数的框架～～～将初始化第一次计算P，X和之后的合并为一体～～～

main.c

#include "pso.h"
#include "stdio.h"


int main(int argc, const char *argv[])
{
	cur_n=0;
	RandInitofSwarm();             //初始化粒子群
	while(cur_n++ != ITE_N)
	{
		UpdatePandGbest();    //更新个体历史最优解P和全局最优解GBest
		UpdateofVandX();        //速度和位置更新，即飞翔
	}

	getchar();
	return 0;
}

pso.h

#ifndef _PSO_H_
#define _PSO_H_

/*各种适应度函数选择，要用哪个，就设置为1,但只能有一个为1*/
#define TEST 0  
#define GOLDSTEIN_PRICE	0
#define SCHAFFER 1
#define HANSEN	0
#define NEEDLE	0

#define Dim 2		//粒子维度
#define PNum 20		//种群规模
#define ITE_N  50  	//最大迭代次数
int cur_n;			//当前迭代次数

/*惯性权重函数*/
#define W_START 1.4
#define W_END	0.4
#define W_FUN	(W_START-(W_START-W_END)*pow((double)cur_n/ITE_N,2))
/*个体和种群结构体*/
struct PARTICLE
{
	double X[Dim];
	double P[Dim];
	double V[Dim];
	double Fitness;
} particle;

struct SWARM
{
	struct PARTICLE Particle[PNum];
	int GBestIndex;
	double GBest[Dim];
	double C1;
	double C2;
	double Xup[Dim];
	double Xdown[Dim];
	double Vmax[Dim];
} swarm;
/*是的，只要三个就好了，更好理解点感觉*/
void RandInitofSwarm(void);
void UpdateofVandX(void);
void UpdatePandGbest(void);

#endif

pso.c

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


//初始化种群
void RandInitofSwarm(void)
{
	int i, j;
        //学习因子C1,C2
	swarm.C1 = 2.0;
	swarm.C2 = 2.0;
	for(j = 0; j < Dim; j++)
	{
		swarm.Xdown[j] = -10;    //搜索空间范围
		swarm.Xup[j] = 10;
		swarm.Vmax[j] = 0.1;       //粒子飞翔速度最大值
	}

	srand((unsigned)time(NULL));
	for(i = 0; i < PNum; i++)
	{
		for(j = 0; j < Dim; j++)
		{
			swarm.Particle[i].X[j] = rand() / (double)RAND_MAX * (swarm.Xup[j] - swarm.Xdown[j]) + swarm.Xdown[j];	//-Xdown~Xup
			swarm.Particle[i].V[j] = rand() / (double)RAND_MAX * swarm.Vmax[j] * 2 - swarm.Vmax[j];					//-Vmax~Vmax
		}
	}
}
/*计算适应度函数的适应度，待进一步完善*/
static double ComputAFitness(double X[])
{
	/*
		OK
		min-f(0,0)=3
	*/
#if TEST
	return X[0]*X[0]+X[1]*X[1]+3;
#endif

	/*
		
		min-f(0,-1)=3.x,y-(-2,2)
	*/
#if GOLDSTEIN_PRICE
	return (1+pow(X[0]+X[1]+1,2)*(19-14*X[0]+3*pow(X[0],2)-14*X[1]+6*X[0]*X[1]+ 3*pow(X[1],2))*(30+pow((2*X[0]-3*X[1]),2)*\
		(18-32*X[0]+12*pow(X[0],2)+48*X[1]-36*X[0]*X[1]	+ 27*pow(X[1],2))));
#endif

	/*
		
		min-f(0,0)=-1.x,y-(-4,4)
	*/
#if SCHAFFER
	//函数:Schaffer's F6

	return -1*(0.5-(sin(sqrt(X[0]*X[0]+X[1]*X[1]))*\
		sin(sqrt(X[0]*X[0]+X[1]*X[1]))-0.5)/pow(( 1+0.001*(X[0]*X[0]+X[1]*X[1])),2));

#endif

	/*
		此函数局部极值点有760个。
		min-f(x,y)=-176.541793.x,y-(-10,10)
		(-7.589893,-7.708314),(-7.589893,-1.425128)
		(-7.5898893,4.858057),(-1.306708,-7.708314)
		(-1.306707,-1.425128),(-1.306708,4.858057)
		(4.976478,-7.708314),(4.976478,-1.425128)
		(4.976478,4.858057)
	*/
#if HANSEN
	int i;
	double temp1=0,temp2=0;
	double hansenf=0;
	for (i=1;i <= 5;i++)
	{
		temp1+=i*cos((i-1)*X[0]+i);
		temp2+=i*cos((i-1)*X[1]+i);
	}
	hansenf=-1*temp1*temp2;

	return hansenf;
#endif

	/*
		
		min-f(0,0)=-3600,x,y-(-5.12,5.12)
		该问题的全局最优解被最差解包围,
		局部极值点：(-5.12,5.12),(-5.12,-5.12)
					(5.12,-5.12),(5.12,5.12)
					f=-2748.78
	*/
#if NEEDLE
	return -1*pow((3/(0.05+pow(X[0]-1,2)+pow(X[1]-1,2))),2);
#endif

}

//update  V and X
void UpdateofVandX(void)
{
	int i, j;
	srand((unsigned)time(NULL));
	for(i = 0; i < PNum; i++)
	{
		for(j = 0; j < Dim; j++)
			swarm.Particle[i].V[j] = W_FUN * swarm.Particle[i].V[j] +
			rand() / (double)RAND_MAX * swarm.C1 * (swarm.Particle[i].P[j] - swarm.Particle[i].X[j]) +
			rand() / (double)RAND_MAX * swarm.C2 * (swarm.GBest[j] - swarm.Particle[i].X[j]);
		for(j = 0; j < Dim; j++)
		{
			if(swarm.Particle[i].V[j] > swarm.Vmax[j])
				swarm.Particle[i].V[j] = swarm.Vmax[j];
			if(swarm.Particle[i].V[j] < -swarm.Vmax[j])
				swarm.Particle[i].V[j] = -swarm.Vmax[j];
		}

		for(j = 0; j < Dim; j++)
		{
			swarm.Particle[i].X[j] += swarm.Particle[i].V[j];
			if(swarm.Particle[i].X[j] > swarm.Xup[j])
				swarm.Particle[i].X[j] = swarm.Xup[j];
			if(swarm.Particle[i].X[j] < swarm.Xdown[j])
				swarm.Particle[i].X[j] = swarm.Xdown[j];
		}
	}
}

/*更新个体极值P和全局极值GBest*/
void UpdatePandGbest(void)
{
	int i, j;
	//update of P if the X is bigger than current P
	for (i = 0; i < PNum; i++)
	{
		if (swarm.Particle[i].Fitness < ComputAFitness(swarm.Particle[i].P))
		{
			for(j = 0; j < Dim; j++)
			{
				swarm.Particle[i].P[j] = swarm.Particle[i].X[j];
			}
		}
	}
	//update of GBest
	for(i = 0; i < PNum; i++)
		if(ComputAFitness(swarm.Particle[i].P) < ComputAFitness(swarm.Particle[swarm.GBestIndex].P))
			swarm.GBestIndex = i;
	for(j = 0; j < Dim; j++)
	{
		swarm.GBest[j] = swarm.Particle[swarm.GBestIndex].P[j];
	}
	printf("The %dth iteraction.\n",cur_n);
	printf("GBestIndex:%d \n",swarm.GBestIndex );
	printf("GBest:" );
	for(j=0;j<Dim;j++)
	{
		printf("%.4f ,",swarm.GBest[j]);
	}
	printf("\n" );
	printf("Fitness of GBest: %f \n\n",ComputAFitness(swarm.Particle[swarm.GBestIndex].P));
}

Version 1.0 标准PSO

返回
【上篇】
【下篇】

作者:

抱歉!评论已关闭.

返回首页

Copyright © 2013-2018 学步园  保留所有权利.
软文销售 QQ客服：2265327166 点击这里给我发消息（其他合作也可洽谈）
必威体育
必威电竞