0 賽題思路

（賽題出來以后第一時間在群內(nèi)分享）

2023高教社杯全國大學生數(shù)學建模競賽

資料思路分享Q群：575833716

1 什么是粒子群算法？

粒子群算法（Particle Swarm Optimization,PSO）是一種模仿鳥群、魚群覓食行為發(fā)展起來的一種進化算法。其概念簡單易于編程實現(xiàn)且運行效率高、參數(shù)相對較少，應用非常廣泛。粒子群算法于1995年提出，距今（2019）已有24年歷史。 　　 　　粒子群算法中每一個粒子的位置代表了待求問題的一個候選解。每一個粒子的位置在空間內(nèi)的好壞由該粒子的位置在待求問題中的適應度值決定。每一個粒子在下一代的位置有其在這一代的位置與其自身的速度矢量決定，其速度決定了粒子每次飛行的方向和距離。在飛行過程中，粒子會記錄下自己所到過的最優(yōu)位置 Ｐ，群體也會更新群體所到過的最優(yōu)位置Ｇ 。粒子的飛行速度則由其當前位置、粒子自身所到過的最優(yōu)位置、群體所到過的最優(yōu)位置以及粒子此時的速度共同決定。

2 舉個例子

在一個湖中有兩個人他們之間可以通信，并且可以探測到自己所在位置的最低點。初始位置如上圖所示，由于右邊比較深，因此左邊的人會往右邊移動一下小船。

現(xiàn)在左邊比較深，因此右邊的人會往左邊移動一下小船

一直重復該過程，最后兩個小船會相遇

得到一個局部的最優(yōu)解

將每個個體表示為粒子。每個個體在某一時刻的位置表示為，x(t),方向表示為v(t)

p（t）為在t時刻x個體的自己的最優(yōu)解，g(t)為在t時刻所有個體的最優(yōu)解，v(t)為個體在t時刻的方向，x(t)為個體在t時刻的位置

下一個位置為上圖所示由x,p,g共同決定了

種群中的粒子通過不斷地向自身和種群的歷史信息進行學習，從而可以找到問題的最優(yōu)解。

3 還是一個例子

粒子群算法是根據(jù)鳥群覓食行為衍生出的算法?，F(xiàn)在，我們的主角換成是一群鳥。

小鳥們的目標很簡單，要在這一帶找到食物最充足的位置安家、休養(yǎng)生息。它們在這個地方的搜索策略如下： 　　1． 每只鳥隨機找一個地方，評估這個地方的食物量。 　　2． 所有的鳥一起開會，選出食物量最多的地方作為安家的候選點G。 　　3． 每只鳥回顧自己的旅程，記住自己曾經(jīng)去過的食物量最多的地方P。 　　4． 每只鳥為了找到食物量更多的地方，于是向著G飛行，但是呢，不知是出于選擇困難癥還是對P的留戀，或者是對G的不信任，小鳥向G飛行時，時不時也向P飛行，其實它自己也不知道到底是向G飛行的多還是向P飛行的多。 　　5． 又到了開會的時間，如果小鳥們決定停止尋找，那么它們會選擇當前的G來安家；否則繼續(xù)2->3->4->5來尋找它們的棲息地。

上圖描述的策略4的情況，一只鳥在點A處，點G是鳥群們找到過的食物最多的位置，點P是它自己去過的食物最多的地點。V是它現(xiàn)在的飛行速度（速度是矢量，有方向和大小)，現(xiàn)在它決定向著P和G飛行，但是這是一只佛系鳥，具體飛多少隨緣。如果沒有速度V，它應該飛到B點，有了速度V的影響，它的合速度最終使它飛到了點C，這里是它的下一個目的地。如果C比P好那么C就成了下一次的P，如果C比G好，那么就成了下一次的G。

算法流程

算法實現(xiàn)

這里學長用python來給大家演示使用粒子群解函數(shù)最優(yōu)解

import numpy as np
import matplotlib.pyplot as plt
import random


# 定義“粒子”類
class parti(object):
? ?def __init__(self, v, x):
? ? ? ?self.v = v ? ? ? ? ? ? ? ? ? ?# 粒子當前速度
? ? ? ?self.x = x ? ? ? ? ? ? ? ? ? ?# 粒子當前位置
? ? ? ?self.pbest = x ? ? ? ? ? ? ? ?# 粒子歷史最優(yōu)位置

class PSO(object):
? ?def __init__(self, interval, tab='min', partisNum=10, iterMax=1000, w=1, c1=2, c2=2):
? ? ? ?self.interval = interval ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?# 給定狀態(tài)空間 - 即待求解空間
? ? ? ?self.tab = tab.strip() ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?# 求解最大值還是最小值的標簽: 'min' - 最小值；'max' - 最大值
? ? ? ?self.iterMax = iterMax ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?# 迭代求解次數(shù)
? ? ? ?self.w = w ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?# 慣性因子
? ? ? ?self.c1, self.c2 = c1, c2 ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? # 學習因子
? ? ? ?self.v_max = (interval[1] - interval[0]) * 0.1 ? ? ? ? ? ? ? ? ? ? ?# 設置最大遷移速度
? ? ? ?#####################################################################
? ? ? ?self.partis_list, self.gbest = self.initPartis(partisNum) ? ? ? ? ? ? ? ? # 完成粒子群的初始化，并提取群體歷史最優(yōu)位置
? ? ? ?self.x_seeds = np.array(list(parti_.x for parti_ in self.partis_list)) ? ?# 提取粒子群的種子狀態(tài) ###
? ? ? ?self.solve() ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?# 完成主體的求解過程
? ? ? ?self.display() ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?# 數(shù)據(jù)可視化展示

? ?def initPartis(self, partisNum):
? ? ? ?partis_list = list()
? ? ? ?for i in range(partisNum):
? ? ? ? ? ?v_seed = random.uniform(-self.v_max, self.v_max)
? ? ? ? ? ?x_seed = random.uniform(*self.interval)
? ? ? ? ? ?partis_list.append(parti(v_seed, x_seed))
? ? ? ?temp = 'find_' + self.tab
? ? ? ?if hasattr(self, temp): ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? # 采用反射方法提取對應的函數(shù)
? ? ? ? ? ?gbest = getattr(self, temp)(partis_list)
? ? ? ?else:
? ? ? ? ? ?exit('>>>tab標簽傳參有誤："min"|"max"<<<')
? ? ? ?return partis_list, gbest

? ?def solve(self):
? ? ? ?for i in range(self.iterMax):
? ? ? ? ? ?for parti_c in self.partis_list:
? ? ? ? ? ? ? ?f1 = self.func(parti_c.x)
? ? ? ? ? ? ? ?# 更新粒子速度，并限制在最大遷移速度之內(nèi)
? ? ? ? ? ? ? ?parti_c.v = self.w * parti_c.v + self.c1 * random.random() * (parti_c.pbest - parti_c.x) + self.c2 * random.random() * (self.gbest - parti_c.x)
? ? ? ? ? ? ? ?if parti_c.v > self.v_max: parti_c.v = self.v_max
? ? ? ? ? ? ? ?elif parti_c.v < -self.v_max: parti_c.v = -self.v_max
? ? ? ? ? ? ? ?# 更新粒子位置，并限制在待解空間之內(nèi)
? ? ? ? ? ? ? ?if self.interval[0] <= parti_c.x + parti_c.v <=self.interval[1]:
? ? ? ? ? ? ? ? ? ?parti_c.x = parti_c.x + parti_c.v
? ? ? ? ? ? ? ?else:
? ? ? ? ? ? ? ? ? ?parti_c.x = parti_c.x - parti_c.v
? ? ? ? ? ? ? ?f2 = self.func(parti_c.x)
? ? ? ? ? ? ? ?getattr(self, 'deal_'+self.tab)(f1, f2, parti_c) ? ? ? ? ? ? # 更新粒子歷史最優(yōu)位置與群體歷史最優(yōu)位置

? ?def func(self, x): ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? # 狀態(tài)產(chǎn)生函數(shù) - 即待求解函數(shù)
? ? ? ?value = np.sin(x**2) * (x**2 - 5*x)
? ? ? ?return value

? ?def find_min(self, partis_list): ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? # 按狀態(tài)函數(shù)最小值找到粒子群初始化的歷史最優(yōu)位置
? ? ? ?parti = min(partis_list, key=lambda parti: self.func(parti.pbest))
? ? ? ?return parti.pbest

? ?def find_max(self, partis_list):
? ? ? ?parti = max(partis_list, key=lambda parti: self.func(parti.pbest)) ? # 按狀態(tài)函數(shù)最大值找到粒子群初始化的歷史最優(yōu)位置
? ? ? ?return parti.pbest

? ?def deal_min(self, f1, f2, parti_):
? ? ? ?if f2 < f1: ? ? ? ? ? ? ? ? ? ? ? ? ?# 更新粒子歷史最優(yōu)位置
? ? ? ? ? ?parti_.pbest = parti_.x
? ? ? ?if f2 < self.func(self.gbest):
? ? ? ? ? ?self.gbest = parti_.x ? ? ? ? ? ?# 更新群體歷史最優(yōu)位置

? ?def deal_max(self, f1, f2, parti_):
? ? ? ?if f2 > f1: ? ? ? ? ? ? ? ? ? ? ? ? ?# 更新粒子歷史最優(yōu)位置
? ? ? ? ? ?parti_.pbest = parti_.x
? ? ? ?if f2 > self.func(self.gbest):
? ? ? ? ? ?self.gbest = parti_.x ? ? ? ? ? ?# 更新群體歷史最優(yōu)位置

? ?def display(self):
? ? ? ?print('solution: {}'.format(self.gbest))
? ? ? ?plt.figure(figsize=(8, 4))
? ? ? ?x = np.linspace(self.interval[0], self.interval[1], 300)
? ? ? ?y = self.func(x)
? ? ? ?plt.plot(x, y, 'g-', label='function')
? ? ? ?plt.plot(self.x_seeds, self.func(self.x_seeds), 'b.', label='seeds')
? ? ? ?plt.plot(self.gbest, self.func(self.gbest), 'r*', label='solution')
? ? ? ?plt.xlabel('x')
? ? ? ?plt.ylabel('f(x)')
? ? ? ?plt.title('solution = {}'.format(self.gbest))
? ? ? ?plt.legend()
? ? ? ?plt.savefig('PSO.png', dpi=500)
? ? ? ?plt.show()
? ? ? ?plt.close()


if __name__ == '__main__':
? ?PSO([-9, 5], 'max')

效果

建模資料

資料分享: 最強建模資料

2023高教社杯全國大學生數(shù)學建模競賽

資料思路分享Q群：575833716