PSO也写了一个多月了吧，记得是今年3月底完成了
差不多是破釜沉舟式写出来的，当时写出来那个激动哇。。
因为就花一天，而且很简单，更主要的是，写过了JADE，这个写起来就思路清晰很多了。
我觉得这个问题是比较有意思的一个问题，它是模拟鸟群捕食行为的一种算法：
设想这样一个场景：一群鸟在随机搜索食物。在这个区域里只有一块食物。所有的鸟都不知道食物在那里。但是他们知道当前的位置离食物还有多远。那么找到食物的最优策略是什么呢。最简单有效的就是搜寻目前离食物最近的鸟的周围区域。
但是这个问题又继而演化成了一个物理问题了（这里得打个问号了，是我自己这么认为的）
PSO最主要的是有一个速度公式和一个位移公式，如下：
v[] = w * v[] + c1 * rand() * (pbest[] - present[]) + c2 * rand() * (gbest[] - present[])
present[] = persent[] + v[]
这里面的3个参数w,c1,c2是非常之重要的，w是惯性权重，如果它太大，相当于运动的惯性系数太大，会突然间超过目标（停不下来。。。），如果太小，那么飞行到目标解又得花很长的时间。更奇葩的应该是c1,c2的选取吧，等下看我程序里给的值吧（其实这里，标准算法的值是w=1.0,c1=c2=2.0的）
先看下PSO的基本流程吧
再贴一下百度上找的伪代码：
clear all;
clc;
format long;
%------给定初始化条件----------------------------------------------
c1=1.4962; %学习因子1
c2=1.4962; %学习因子2
w=0.7298; %惯性权重
MaxDT=1000; %最大迭代次数
D=10; %搜索空间维数（未知数个数）
N=40; %初始化群体个体数目
eps=10^(-6); %设置精度(在已知最小值时候用)
%------初始化种群的个体(可以在这里限定位置和速度的范围)------------
for i=1:N
for j=1:D
x(i,j)=randn; %随机初始化位置
v(i,j)=randn; %随机初始化速度
end
end
%------先计算各个粒子的适应度，并初始化Pi和Pg----------------------
for i=1:N
p(i)=f(x(i,:),D);
y(i,:)=x(i,:);
end
pg=x(1,:); %Pg为全局最优
for i=2:N
if f(x(i,:),D)<f(pg,D)
pg=x(i,:);
end
end
%------进入主要循环，按照公式依次迭代，直到满足精度要求------------
for t=1:MaxDT
for i=1:N
v(i,:)=w*v(i,:)+c1*rand*(y(i,:)-x(i,:))+c2*rand*(pg-x(i,:));
x(i,:)=x(i,:)+v(i,:);
if f(x(i,:),D)<p(i)
p(i)=f(x(i,:),D);
y(i,:)=x(i,:);
end
if p(i)<f(pg,D)
pg=y(i,:);
end
end
Pbest(t)=f(pg,D);
end
%------最后给出计算结果
disp('*************************************************************')
disp('函数的全局最优位置为：')
Solution=pg'
disp('最后得到的优化极值为：')
Result=f(pg,D)
disp('*************************************************************')
%------算法结束---DreamSun GL & HF-----------------------------------最后再贴下我自己写的代码吧（测试函数依旧用的JADE里的第一个，精度也非常精确了~~）
代码风格依旧很挫。。。。
#include<stdio.h>
#include<time.h>
#include<stdlib.h>
#include<iostream>
//#include<random>
#include<cmath>
#include<malloc.h>
using namespace std;
#define URAND rand()/(RAND_MAX+1.0)
const int P_num=100;
const double PI=acos(-1.0);
const int G=1500;
const double w=0.7162; //用0.72//// 崔文明用的0.7162，因此我也改用了，就相当准确了，不过标准的是0.729
//可以考虑线性变化的w的取值...
//w=0.9-t/T*0.5,不过测试之后，效果更差了
const double c1=1.49445;
const double c2=1.49445;
const int D=30;
int times;
struct individual
{
double p[D+10];
double v[D+10];
double fitness;
double pbest[D+10];
double val;
} gbest,P[P_num+10];
double rand_low_high(double low,double high)
{
return (high-low)*URAND+low;
}
double calc_val(individual P)
{
double val=0;
for(int i=1;i<=D;i++)
val+=P.p[i]*P.p[i];
return val;
}
void calc_P_val()
{
for(int i=1; i<=P_num; i++)
{
P[i].val=calc_val(P[i]);
}
}
void init_P()//初始化
{
srand((unsigned)time(NULL));
for(int i=1; i<=P_num; i++)
{
for(int j=1; j<=D; j++)
{
P[i].p[j]=rand_low_high(-100,100);
P[i].v[j]=rand_low_high(-100,100);
}
}
calc_P_val();
//找出第一代里的pbest和gbest
for(int i=1;i<=P_num;i++)
{
for(int j=1;j<=D;j++)
P[i].pbest[j]=P[i].p[j];
P[i].fitness=P[i].val;
}
gbest=P[1];
//for(int i=2;i<=P_num;i++)
//if(P[i].fitness>gbest.fitness)
//gbest=P[i];
}
int main()
{
srand((unsigned)time(NULL));
init_P();
// double w=0.7;
for(times=1; times<=G; times++)//其实更新应该是可以从times=2开始的
//for(times=2; times<=G; times++)
{
for(int j=1;j<=P_num;j++)
{
for(int k=1;k<=D;k++)
{
P[j].v[k]=w*P[j].v[k]+c1*URAND*(P[j].pbest[k]-P[j].p[k])+
c2*URAND*(gbest.p[k]-P[j].p[k]);
if(P[j].v[k]>100)P[j].v[k]=100.0;
// printf("%f\n",P[j].v[k]);
if(P[j].v[k]<-100)P[j].v[k]=-100.0;
P[j].p[k]+=P[j].v[k];
if(P[j].p[k]>100)P[j].p[k]=100;
// printf("%f\n",P[j].p[k]);
if(P[j].p[k]<-100)P[j].p[k]=-100;
}
P[j].val=calc_val(P[j]);
if(P[j].val<P[j].fitness)
{
for(int k=1;k<=D;k++)
P[j].pbest[k]=P[j].p[k];
P[j].fitness=P[j].val;
}
//收敛速度更快？？？？
if(P[j].fitness<gbest.fitness)gbest=P[j];
}
// w=0.7-times/G*0.5;
//标准算法是在每代开始更新gbest的
//for(int j=1;j<=P_num;j++)if(P[j].fitness<gbest.fitness)gbest=P[j];
}
gbest.val=calc_val(gbest);
cout<<"the best value: "<<gbest.val<<endl<<endl;
printf("time: %.6f\n\n\n\n",(double)clock()/CLOCKS_PER_SEC);
system("pause");
return 0;
}
本文章为转载内容，我们尊重原作者对文章享有的著作权。如有内容错误或侵权问题，欢迎原作者联系我们进行内容更正或删除文章。
