粒子群算法Python实现

东明山庄 2015-04-13 12:49:52

5835

分类专栏：软计算算法列传文章标签：粒子群

最后发布:2015-04-13 12:49:52首次发布:2015-04-13 12:49:52

本文链接：https://blog.csdn.net/golden1314521/article/details/45024307

版权

根据适应度值来计算局部最优和全局最优。前者表示迭代以来，该粒子的最优解；后者是局部最优的最优。
局部最优对应了局部学习率，全局最优对应了全局学习率，用来计算更新速度。更新速度负责对粒子位置进行更新。
粒子的速度有个范围，例子的位置也有个范围。

寻找最优解的函数为
$\frac{\sin\sqrt{x^2+y^2}}{\sqrt{x^2+y^2}} +e^{\frac{\cos2\pi x +cos2\pi y}{2}}-2.71289$
其函数图像如下
这里写图片描述
该三维图像代码如下：

from mpl_toolkits.mplot3d import Axes3D

fig = figure()
ax = Axes3D(fig)
X = np.arange(-2, 2, 0.05)
Y = np.arange(-2, 2, 0.05)
X, Y = np.meshgrid(X, Y)
Z = np.sin(np.sqrt(X**2+Y**2))/np.sqrt(X**2+Y**2)+np.exp((np.cos(2*np.pi*X)+np.cos(2*np.pi*Y))/2)-2.71289
ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap='hot')

粒子群算法Python实现源码：

import numpy as np

def getweight():
    # 惯性权重
    weight = 1
    return weight

def getlearningrate():
    # 分别是粒子的个体和社会的学习因子，也称为加速常数
    lr = (0.49445,1.49445)
    return lr

def getmaxgen():
    # 最大迭代次数
    maxgen = 300
    return maxgen

def getsizepop():
    # 种群规模
    sizepop = 20
    return sizepop

def getrangepop():
    # 粒子的位置的范围限制,x、y方向的限制相同
    rangepop = (-2,2)
    return rangepop

def getrangespeed():
    # 粒子的速度范围限制
    rangespeed = (-0.5,0.5)
    return rangespeed

def func(x):
    # x输入粒子位置
    # y 粒子适应度值
    if (x[0]==0)&(x[1]==0):
        y = np.exp((np.cos(2*np.pi*x[0])+np.cos(2*np.pi*x[1]))/2)-2.71289
    else:
        y = np.sin(np.sqrt(x[0]**2+x[1]**2))/np.sqrt(x[0]**2+x[1]**2)+np.exp((np.cos(2*np.pi*x[0])+np.cos(2*np.pi*x[1]))/2)-2.71289
    return y

def initpopvfit(sizepop):
    pop = np.zeros((sizepop,2))
    v = np.zeros((sizepop,2))
    fitness = np.zeros(sizepop)

    for i in range(sizepop):
        pop[i] = [(np.random.rand()-0.5)*rangepop[0]*2,(np.random.rand()-0.5)*rangepop[1]*2]
        v[i] = [(np.random.rand()-0.5)*rangepop[0]*2,(np.random.rand()-0.5)*rangepop[1]*2]
        fitness[i] = func(pop[i])
    return pop,v,fitness

def getinitbest(fitness,pop):
    # 群体最优的粒子位置及其适应度值
    gbestpop,gbestfitness = pop[fitness.argmax()].copy(),fitness.max()
    #个体最优的粒子位置及其适应度值,使用copy()使得对pop的改变不影响pbestpop，pbestfitness类似
    pbestpop,pbestfitness = pop.copy(),fitness.copy()
    
    return gbestpop,gbestfitness,pbestpop,pbestfitness  

w = getweight()
lr = getlearningrate()
maxgen = getmaxgen()
sizepop = getsizepop()
rangepop = getrangepop()
rangespeed = getrangespeed()

pop,v,fitness = initpopvfit(sizepop)
gbestpop,gbestfitness,pbestpop,pbestfitness = getinitbest(fitness,pop)

result = np.zeros(maxgen)
for i in range(maxgen):
        
        #速度更新
        for j in range(sizepop):
            v[j] += lr[0]*np.random.rand()*(pbestpop[j]-pop[j])+lr[1]*np.random.rand()*(gbestpop-pop[j])
        v[v<rangespeed[0]] = rangespeed[0]
        v[v>rangespeed[1]] = rangespeed[1]
        
        #粒子位置更新
        for j in range(sizepop):
            pop[j] += 0.5*v[j]
        pop[pop<rangepop[0]] = rangepop[0]
        pop[pop>rangepop[1]] = rangepop[1]
        
        #适应度更新
        for j in range(sizepop):
            fitness[j] = func(pop[j])
            
        for j in range(sizepop):
            if fitness[j] > pbestfitness[j]:
                pbestfitness[j] = fitness[j]
                pbestpop[j] = pop[j].copy()
        
        if pbestfitness.max() > gbestfitness :
            gbestfitness = pbestfitness.max()
            gbestpop = pop[pbestfitness.argmax()].copy()
        
        result[i] = gbestfitness
                
            
plot(result)

全局最优解随迭代次数增加的变化情况如下图所示：

这里写图片描述