用MATLAB编写PSO算法及实例.docx
《用MATLAB编写PSO算法及实例.docx》由会员分享，可在线阅读，更多相关《用MATLAB编写PSO算法及实例.docx（7页珍藏版）》请在冰点文库上搜索。
用MATLAB编写PSO算法及实例
1.1粒子群算法
PSO从这种模型中得到启示并用于解决优化问题。
PSO中，每个优化问题的潜在解都是搜索空间中的一只鸟，称之为粒子。
所有的粒子都有一个由被优化的函数决定的适值（fitnessvalue），每个粒子还有一个速度决定它们飞翔的方向和距离。
然后粒子们就追随当前的最优粒子在解空间中搜索。
PSO初始化为一群随机粒子（随机解），然后通过迭代找到最优解。
在每一次迭代中，粒子通过跟踪两个极值来更新自己；第一个就是粒子本身所找到的最优解，这个解称为个体极值；另一个极值是整个种群目前找到的最优解，这个极值是全局极值。
另外也可以不用整个种群而只是用其中一部分作为粒子的邻居，那么在所有邻居中的极值就是局部极值。
假设在一个维的目标搜索空间中，有个粒子组成一个群落，其中第个粒子表示为一个维的向量
，。
第个粒子的“飞行”速度也是一个维的向量，记为
，。
第个粒子迄今为止搜索到的最优位置称为个体极值，记为
，。
整个粒子群迄今为止搜索到的最优位置为全局极值，记为
在找到这两个最优值时，粒子根据如下的公式（1.1）和（1.2）来更新自己的速度和位置：
（1.1）
（1.2）
其中：
和为学习因子，也称加速常数（accelerationconstant），和为[0，1]范围内的均匀随机数。
式（1.1）右边由三部分组成，第一部分为“惯性（inertia）”或“动量（momentum）”部分，反映了粒子的运动“习惯（habit）”，代表粒子有维持自己先前速度的趋势；第二部分为“认知（cognition）”部分，反映了粒子对自身历史经验的记忆（memory）或回忆（remembrance），代表粒子有向自身历史最佳位置逼近的趋势；第三部分为“社会（social）”部分，反映了粒子间协同合作与知识共享的群体历史经验。
二、算法设计
2.1算法流程图
2.2算法实现
算法的流程如下：
①初始化粒子群，包括群体规模，每个粒子的位置和速度
②计算每个粒子的适应度值；%它的适应度就是指目标函数的值。
一般来说，目标函数的选择由具体问题来决定，假如是背包问题，适应度即放入包中物体的总价格。
初始粒子位置和速度的位置一般随机产生。
但是在某些领域，如果已有其他的算法可以产生可行解的话，可以用这个可行解来初始化，这样更容易得到最优的解
③对每个粒子，用它的适应度值和个体极值比较，如果，则用替换掉；
④对每个粒子，用它的适应度值和全局极值比较，如果则用替；
⑤根据公式（1.1），（1.2）更新粒子的速度和位置；
⑥如果满足结束条件（误差足够好或到达最大循环次数）退出，否则返回②。
2.3参数选择
本算法中主要的参数变量为（惯性权值），，（加速因子），N（种群数），M（迭代次数），D（粒子维数）。
（1）种群规模
通常，种群太小则不能提供足够的采样点，以致算法性能很差；种群太大尽管可以增加优化信息，阻止早熟收敛的发生，但无疑会增加计算量，造成收敛时间太长，表现为收敛速度缓慢。
种群规模一般设为100~1000。
本文选择种群规模为100。
（2）最大迭代次数
迭代次数越多能保证解的收敛性，但是影响运算速度，本文选1000次。
（3）惯性权值
惯性权重表示在多大程度上保留原来的速度。
较大，全局收敛能力强，局部收敛能力弱；较小，局部收敛能力强，全局收敛能力弱。
本文选0.6。
（4）加速因子
加速常数和分别用于控制粒子指向自身或邻域最佳位置的运动。
文献[20]建议，并通常取。
本文也取。
（5）粒子维数
本文中粒子维数取决于待优化函数的维数。
需要说明的是，本文的程序允许改变这些参数，因为本文编写的程序参照matlab工具箱，留给用户解决这类问题一个接口函数，上述的各个参数正是接口函数的参数，因此允许改变。
另外对于和c也可采用变参数法，即随迭代次数增加，利用经验公式使它们动态调整，本文采用固定值。
3.1求三维函数f=x
（1）.^2+x
（2）.^2+x（3）.^2的最小值
步骤：
1.初始化x,v;2.求出每个粒子的适应值；3.初始化pb,pg个体最优和全局最优；4.根据式子更新x,v;5.是否满足条件，满足跳出循环，否则重复2-4步
尝试编码：
（1）pso.m文件
%此算法是PSO算法，汪汪的20161024号版本
function[xm,fv]=PSO（fitness,N,c1,c2,w,M,D）
%{
xm,fv算法最后得到的最优解时的x及最优解，fitness为适应度，即要优化的目标函数，
N为种群数量，c1,c2为学习因子，w为惯性权重，M为迭代次数，D为粒子的维数
%}
formatlong;
%初始化种群
fori=1:
N
forj=1:
D
x（i,j）=randn;%随机初始化位置
v（i,j）=randn;%随机初始化速度
end
end
%先计算各个粒子的适应度pi，并初始化y-粒子个体极值，pg-全局极值
fori=1:
N
p（i）=fitness（x（i,:
））;%适应度问题：
将x（i,:
）改成x（i,j）是否可以，答不能
y（i,:
）=x（i,:
）;%个体极值
end
pg=x（N,:
）;%初始化全局极值/最优
fori=1:
N-1
iffitness（x（i,:
））pg=x（i,:
）;%替换并选出全局极值
end
end
%进入粒子群算法主要循环，更新v及x
fort=1:
M
fori=1:
N
v（i,:
）=w*v（i,:
）+c1*rand*（y（1,:
）-x（1,:
））+c2*rand*（pg-x（i,:
））;
x（i,:
）=x（i,:
）+v（i,:
）;
iffitness（x（i,:
））p（i）=fitness（x（i,:
））;
y（i,:
）=x（i,:
）;
end
ifp（i）pg=y（i,:
）;
end
end
pbest（t）=fitness（pg）;%M次迭代后最优解
end
xm=pg';%为何要共轭转置？
fv=fitness（pg）;
（2）目标函数fitness.m文件
functionf=fitness（x）
f=x
（1）.^2+x
（2）.^2+x（3）.^2;
end
需要说明的是，针对不同的函数优化，只需要改变目标函数就可以。
（3）在命令行输入或建立调用m文件
在命令行先后输入[xm,fv]=PSO（@fitness,100,2,2,0.6,1000,3），或建立包涵该语句的m文件，运行即可得到结果。
四、结果与分析
xm=1.0e-04*
-0.285730140229565
0.676783696397148
-0.250529540096653
fv=
6.024429352056337e-09
fv是最优值，xm为最优值对应的自变量值。
3.2高斯函数
[xy]=meshgrid（-100:
100,-100:
100）;%生成xy坐标系且选取范围
sigma=50;%高斯函数中的sigma为50
img=（1/（2*pi*sigma^2））*exp（-（x.^2+y.^2）/（2*sigma^2））;
%目标函数，高斯函数
mesh（img）;%三维画图
holdon;
升级会员 