算法原理

粒子群算法思想來源于實際生活中鳥捕食的過程。假設在一個n維的空間中，有一群鳥（m只）在捕食，食物位于n維空間的某個點上，對于第i只鳥某一時刻來說，有兩個向量描述，一個是鳥的位置向量，第二個是鳥的速度。假設鳥能夠判斷一個位置的好壞，所謂“好壞”，就是離食物更近了還是更遠了。鳥在捕食的過程中會根據自己的經驗以及鳥群中的其他鳥的位置決定自己的速度，根據當前的位置和速度，可以得到下一刻的位置，這樣每只鳥通過向自己和鳥群學習不斷的更新自己的速度位置，最終找到食物，或者離食物足夠近的點。更新速度和位置的表達式如下。

算法流程

算例

語言python 3.7

求解下列函數的最小值

步驟1：繪制函數圖像

設置a=10

import numpy as np import matplotlib.pyplot as plt from matplotlib import cm from mpl_toolkits.mplot3d import Axes3D# 生成X和Y的數據 X = np.arange(-5, 5, 0.1) Y = np.arange(-5, 5, 0.1) X, Y = np.meshgrid(X, Y)# 目標函數 A = 10 Z = 2 * A + X ** 2 - A * np.cos(2 * np.pi * X) + Y ** 2 - A * np.cos(2 * np.pi * Y)# 繪圖 fig = plt.figure() ax = Axes3D(fig) surf = ax.plot_surface(X, Y, Z, cmap=cm.coolwarm) plt.show()

可以發現此函數有很多尖波 ，有很多極大值和極小值。尋找全局最小值比較困難。

算法實現

由于有2個變量（x,y）假設有N個粒子，則V和X都是N*2的矩陣

# 速度 # Vi+1 = w*Vi + c1 * r1 * (pbest_i - Xi) + c2 * r2 * (gbest_i - Xi) # 位置 # Xi+1 = Xi + Vi+1 # vi, xi 分別表示粒子第i維的速度和位置 # pbest_i, gbest_i 分別表示某個粒子最好位置第i維的值、整個種群最好位置第i維的值

首先定義3個函數，分別是適應度函數fitness_func，速度更新函數velocity_update、位置更新函數position_update

import numpy as np import matplotlib.pyplot as plt import matplotlib as mplmpl.rcParams['font.sans-serif'] = ['SimHei'] # 指定默認字體 mpl.rcParams['axes.unicode_minus'] = False # 解決保存圖像是負號'-'顯示為方塊的問題def fitness_func(X):"""計算粒子的的適應度值，也就是目標函數值，X 的維度是 size * 2 """A = 10pi = np.pix = X[:, 0]y = X[:, 1]return 2 * A + x ** 2 - A * np.cos(2 * pi * x) + y ** 2 - A * np.cos(2 * pi * y)def velocity_update(V, X, pbest, gbest, c1, c2, w, max_val):"""根據速度更新公式更新每個粒子的速度#種群大小20:param V: 粒子當前的速度矩陣，20*2 的矩陣:param X: 粒子當前的位置矩陣，20*2 的矩陣:param pbest: 每個粒子歷史最優位置，20*2 的矩陣:param gbest: 種群歷史最優位置，1*2 的矩陣"""size = X.shape[0]#種群大小r1 = np.random.random((size, 1))#生成size 個 0-1之間的隨機數r2 = np.random.random((size, 1))#生成size 個 0-1之間的隨機數V = w * V + c1 * r1 * (pbest - X) + c2 * r2 * (gbest - X) # 直接對照公式寫就好了# 防止越界處理V[V < -max_val] = -max_valV[V > max_val] = max_valreturn Vdef position_update(X, V):"""根據公式更新粒子的位置:param X: 粒子當前的位置矩陣，維度是 20*2:param V: 粒子當前的速度舉著，維度是 20*2"""return X + V

這是PSO 主函數，首先定義粒子群算法的參數，然后在每一輪的迭代中，更新局部最優和全局最優，直到滿足最大迭代次數。

def pso():# PSO的參數w = 1 # 慣性因子，一般取1c1 = 2 # 學習因子，一般取2c2 = 2 #r1 = None # 為兩個（0,1）之間的隨機數r2 = Nonedim = 2 # 維度的維度size = 20 # 種群大小，即種群中小鳥的個數iter_num = 1000 # 算法最大迭代次數max_val = 0.5 # 限制粒子的最大速度為0.5best_fitness = float(9e10) # 初始的適應度值，在迭代過程中不斷減小這個值fitneess_value_list = [] # 記錄每次迭代過程中的種群適應度值變化# 初始化種群的各個粒子的位置# 用一個 20*2 的矩陣表示種群，每行表示一個粒子X = np.random.uniform(-5, 5, size=(size, dim))#2維 表示x、y# 初始化種群的各個粒子的速度V = np.random.uniform(-0.5, 0.5, size=(size, dim))##2維 表示x、y的速度# 計算種群各個粒子的初始適應度值p_fitness = fitness_func(X)# 計算種群的初始最優適應度值g_fitness = p_fitness.min()# 講添加到記錄中fitneess_value_list.append(g_fitness)# 初始的個體最優位置和種群最優位置pbest = Xgbest = X[p_fitness.argmin()]# 接下來就是不斷迭代了for i in range(1, iter_num):V = velocity_update(V, X, pbest, gbest, c1, c2, w, max_val) # 更新速度X = position_update(X, V) # 更新位置p_fitness2 = fitness_func(X) # 計算各個粒子的適應度g_fitness2 = p_fitness2.min() # 計算群體的最優適應度# 更新每個粒子的歷史最優位置for j in range(size):if p_fitness[j] > p_fitness2[j]:pbest[j] = X[j]p_fitness[j] = p_fitness2[j]# 更新群體的最優位置if g_fitness > g_fitness2:gbest = X[p_fitness2.argmin()]g_fitness = g_fitness2# 記錄最優迭代結果fitneess_value_list.append(g_fitness)# 迭代次數+1i += 1# 打印迭代的結果print("最優值是：%.5f" % fitneess_value_list[-1])print("最優解是：x=%.5f, y=%.5f" % (gbest[0],gbest[1]))# 最優值是：0.00000# 最優解是：x=0.00000, y=-0.00000# 繪圖plt.plot(fitneess_value_list, color='r')plt.title('迭代過程')plt.show() pso()

作者：電氣 余登武

總結

以上是生活随笔為你收集整理的粒子群算法求解无约束优化问题 源码实现的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。