介绍看此文:http://blog.csdn.net/fovwin/article/details/8069606
粒子群算法的流程图如上,看了好多版本,这个最靠谱,我的main函数完全按照这个来,好理解过程:
int main(int argc, const char *argv[])
{
int n=0;
//printf("Random Initialization of the swarm:\n\n");
RandInitofSwarm();
//printf("Computation of the fitness of each particle:\n");
ComputFitofSwarm();
//printf("FirstComputPandGbest:\n");
FirstComputPandGbest();
while(n++!=N)
{
printf("The %dth time to calculate .\n",n );
//printf("Updated of the swarm:\n\n");
UpdateofVandX();
//printf("Updated of the swarm's Fitness:\n");
ComputFitofSwarm();
//printf("Replaced of P and Gbest:\n\n");
UpdatePandGbest();
}
getchar();
return 0;
}
各种printf用来调试自己的信息,可以去掉看看中间值~~~
先来pso.h
现在执行最简单的功能,用粒子群搜索二维空间:y=x^2+y^+3的最大值~~~x,y的范围属于(-100,100),那最大值就是20003,伪随机粒子数20个,最大迭代次数为40次,最大搜索速度为2,惯性权重W设为1.4,学习因子c1=c2=2,看能不能找到最大值~~~
#ifndef _PSO_H_
#define _PSO_H_
#define Dim 2
#define PNum 20
#define N 40
typedef struct PARTICLE{
double X[Dim];
double P[Dim];
double V[Dim];
double Fitness;
}particle;
typedef struct SWARM{
particle Particle[PNum];
int GBestIndex;
double GBest[Dim];
double W;
double C1;
double C2;
double Xup[Dim];
double Xdown[Dim];
double Vmax[Dim];
}swarm;
void RandInitofSwarm(void);
void ComputFitofSwarm(void);
void FirstComputPandGbest(void);
void UpdateofVandX(void);
void UpdatePandGbest(void);
double InertWeight(void);
int CriteriaofStop(void);
#endif
pso.c
#include "stdio.h"
#include "stdlib.h"
#include "time.h"
#include "math.h"
#include "pso.h"
#include "stdio.h"
swarm *s=(swarm *)malloc(sizeof(swarm));
particle *p=(particle *)malloc(sizeof(particle));
int main(int argc, const char *argv[])
{
int n=0;
//printf("Random Initialization of the swarm:\n\n");
RandInitofSwarm();
//printf("Computation of the fitness of each particle:\n");
ComputFitofSwarm();
//printf("FirstComputPandGbest:\n");
FirstComputPandGbest();
while(n++!=N)
{
printf("The %dth time to calculate .\n",n );
//printf("Updated of the swarm:\n\n");
UpdateofVandX();
//printf("Updated of the swarm's Fitness:\n");
ComputFitofSwarm();
//printf("Replaced of P and Gbest:\n\n");
UpdatePandGbest();
}
getchar();
return 0;
}
//Random Initialization of the swarm
void RandInitofSwarm(void)
{
int i,j;
s->W=1.4;
s->C1=2.0;
s->C2=2.0;
for(j=0;j<Dim;j++)
{
s->Xdown[j] = -100;
s->Xup[j] = 100;
s->Vmax[j] = 2;
}
srand((unsigned)time(NULL));
for(i=0; i<PNum; i++)
{
//printf(" The %dth of X is: ",i);
for(j=0; j<Dim; j++)
{
s->Particle[i].X[j] = rand()/(double)RAND_MAX*(s->Xup[j]-s->Xdown[j])+s->Xdown[j]; //-100~100
s->Particle[i].V[j] = rand()/(double)RAND_MAX*s->Vmax[j]*2-s->Vmax[j]; //-2~2
//printf(" %.2f \n ",s->Particle[i].X[j]);
}
}
}
//Computation of the fitness of each particle
void ComputFitofSwarm(void)
{
int i;
srand((unsigned)time(NULL));
for(i=0; i<PNum; i++)
{
//printf(" The Fitness of %dth Particle: ",i);
s->Particle[i].Fitness = s->Particle[i].X[0]*s->Particle[i].X[0]+s->Particle[i].X[1]*s->Particle[i].X[1]+3;
//printf(" %.2f\n",s->Particle[i].Fitness);
}
}
void FirstComputPandGbest(void)
{
int i,j;
//P=X;
for(i=0; i<PNum; i++)
{
for(j=0;j<Dim;j++)
{
s->Particle[i].P[j]=s->Particle[i].X[j];
}
}
//Computation of GBest
s->GBestIndex = 0;
for(i=0; i<PNum; i++)
if(s->Particle[i].Fitness>=s->Particle[s->GBestIndex].Fitness)
s->GBestIndex = i;
for(j=0;j<Dim;j++)
{
s->GBest[j]=s->Particle[s->GBestIndex].P[j];
}
printf("GBestIndex , GBest , Fitness of GBest:%d ,%.2f ,%.2f ,%.2f \n",
s->GBestIndex ,s->GBest[0],s->GBest[1],
s->Particle[s->GBestIndex].Fitness);
}
//update V and X
void UpdateofVandX(void)
{
int i,j;
srand((unsigned)time(NULL));
for(i=0; i<PNum; i++)
{
//printf(" The %dth of X is: ",i);
for(j=0; j<Dim; j++)
s->Particle[i].V[j] = s->W*s->Particle[i].V[j]+
rand()/(double)RAND_MAX*s->C1*(s->Particle[i].P[j] - s->Particle[i].X[j])+
rand()/(double)RAND_MAX*s->C2*(s->GBest[j] - s->Particle[i].X[j]);
for(j=0; j<Dim; j++)
{
if(s->Particle[i].V[j]>s->Vmax[j])
s->Particle[i].V[j] = s->Vmax[j];
if(s->Particle[i].V[j]<-s->Vmax[j])
s->Particle[i].V[j] = -s->Vmax[j];
}
for(j=0; j<Dim; j++)
{
s->Particle[i].X[j] += s->Particle[i].V[j];
if(s->Particle[i].X[j]>s->Xup[j])
s->Particle[i].X[j]=s->Xup[j];
if(s->Particle[i].X[j]<s->Xdown[j])
s->Particle[i].X[j]=s->Xdown[j];
}
//printf(" %.2f %.2f \n",s->Particle[i].X[0],s->Particle[i].X[1]);
}
}
static double ComputAFitness(double X[])
{
return X[0]*X[0]+X[1]*X[1]+3;
}
void UpdatePandGbest(void)
{
int i,j;
//update of P if the X is bigger than current P
for (i = 0; i < PNum; i++)
{
//printf(" The %dth of P is: ",i);
if (s->Particle[i].Fitness > ComputAFitness(s->Particle[i].P))
{
for(j=0;j<Dim;j++)
{
s->Particle[i].P[j] = s->Particle[i].X[j];
}
}
//printf(" %.2f %.2f \n",s->Particle[i].P[0],s->Particle[i].P[1]);
}
for (i = 0; i < PNum; i++)
{
//printf("The %dth of P's Fitness is : %.2f \n",i,ComputAFitness(s->Particle[i].P));
}
//update of GBest
for(i=0; i<PNum; i++)
if(ComputAFitness(s->Particle[i].P) >= s->Particle[s->GBestIndex].Fitness)
s->GBestIndex = i;
for(j=0;j<Dim;j++)
{
s->GBest[j]=s->Particle[s->GBestIndex].P[j];
}
printf("GBestIndex , GBest , Fitness of GBest:%d ,%.2f ,%.2f ,%.2f \n",
s->GBestIndex ,s->GBest[0],s->GBest[1],
ComputAFitness(s->GBest));
}
double InertWeight(void)
{
return 1.0;
}
int CriteriaofStop()
{
int n=N;
return (n--==0);
}
C语言实现部分没有复杂的语法,要优化下,不过现在我用printf一步一步打印出来的,基本上实现了,关键是理解流程图,幸亏流程图写的简单明了~~~小bug估计还不少~~~
至于内存申请木有释放,就先放着吧,也不多,呵呵~~~但是要释放得需要先释放particle,再释放swarm~~~不然就内存泄漏了~~~
实验结果:
前天刚看到几个非常变态的函数,什么有276个局部最小点的,最大值在最小值周围的~~~可以各种虐啊~~~
Version0.0——基本PSO
论文:在google学术里面搜particle swarm optimization 被引用20000+的论文就是鼻祖了~~~是的,本文就是实现最基本的~~~