QQ登录

只需要一步,快速开始

 注册地址  找回密码
查看: 9|回复: 0
打印 上一主题 下一主题

一群"鸟"帮你解数学题:聊聊粒子群算法

[复制链接]
字体大小: 正常 放大

1192

主题

4

听众

2946

积分

该用户从未签到

跳转到指定楼层
1#
发表于 2026-8-23 17:38 |只看该作者 |正序浏览
|招呼Ta 关注Ta

算法笔记 · 粒子群优化 · 附完整可运行代码

先讲一个真实的故事。

1995 年,两位研究者在模拟鸟群觅食:一群鸟在一片区域里随机飞,谁也不知道食物在哪。但每只鸟都记得自己飞过的最好位置,也能看到整个鸟群目前找到的最好位置。于是每只鸟的飞行策略变得很简单——一部分凭自己的经验继续飞,一部分朝群里公认的好方向靠拢。

结果让人惊讶:不需要任何"领导",不需要知道食物在哪,鸟群竟能自发地聚拢到食物周围。

这两位研究者一位是电子工程师 James Kennedy,一位是社会心理学家 Russell Eberhart。他们把这个模拟简化成几个公式,就成了今天要讲的主角——粒子群优化算法(Particle Swarm Optimization,PSO)。顺带一提,这个算法诞生 25 年后,他们拿到了 IEEE 的诺贝尔级别的荣誉——它至今仍是工程优化里最常用的算法之一。

如果你正在准备数学建模竞赛,这篇文章能帮你在"优化类"赛题里多一件趁手的兵器。

一、把优化问题想成"找食物"

先忘掉算法,想一个具体场景:

你要找一个函数的最小值,比如 $$f(x, y) = x² + y²$$,最小值显然在原点。但现实问题里,函数可能长得像一片丘陵——到处是小山包(局部最优),真正的谷底(全局最优)藏在某个角落。梯度下降这类方法从一点出发,很容易卡在最近的小山包里下不来。

PSO 的思路完全不同:与其派一个人去找,不如派一群鸟去找。

[td]
概念
在算法里是什么
一只鸟 / 粒子(Particle)一个候选解(一组变量值)
位置(Position)解本身,比如 (x, y) = (3.2, -1.5)
速度(Velocity)下一步往哪个方向飞、飞多远
食物的位置函数的最小值点
离食物多近适应度:把位置代入函数算出的值
自己的经验个体极值 pbest:我飞过的最好位置
群体的情报全局极值 gbest:全群找到的最好位置
二、灵魂公式:只有两行

PSO 全部精髓浓缩在两条更新规则里。每一代(每一轮迭代),每只鸟这样飞:


+ Y% t1 }" u2 J8 m: B* |+ h$$v_{k+1} = w \cdot v_k + c_1 r_1 \cdot (pbest - x_k) + c_2 r_2 \cdot (gbest - x_k)$$* {2 a, [$ P5 @, G9 A* d
) u2 L4 O2 D) i- T7 ?6 L3 I
$$x_{k+1} = x_k + v_{k+1}$$
- f  _) B0 J" P/ d7 I. k

第一条是速度更新,即"下一刻往哪飞"。拆开看,一共三股力量在拉扯这只鸟:

[td]
部分
含义
生活翻译
w·vk惯性:保留原来的飞行方向"我正往这边飞呢,别让我急转弯"
c₁r₁·(pbest − x)认知项:朝自己最好位置飞"我以前在那儿吃过食,回去看看"
c₂r₂·(gbest − x)社会项:朝全群最好位置飞"老王说他找到吃的了,跟过去!"

第二条是位置更新,简单到不像话:新位置 = 旧位置 + 速度。就是"按刚才决定的方向飞过去"。

这里的 r₁、r₂ 是 0 到 1 之间的随机数——这点很关键。如果每次都精确地朝 pbest 和 gbest 飞,所有鸟很快会挤成一条线,搜索就死了。加上随机扰动,鸟群才能保持"散着搜索、又大致向好看的方向去"的状态。

一句话记住 PSO: 每只鸟被三股力拉着——惯性(接着飞)、自己的记忆(回自己找到过的好地方)、群体的情报(去大家公认的好地方)。飞到新位置后,更新记忆,再飞。如此反复,鸟群整体就向最优解聚拢。

三、三个参数怎么调[td]
参数
名字
常用值
调大调小会怎样
w惯性权重0.6 ~ 0.8大→搜得广但收敛慢;小→收敛快但易陷局部最优
c₁认知因子1.5 ~ 2大→各飞各的,多样性好
c₂社会因子1.5 ~ 2大→一窝蜂跟风,收敛快但容易早熟

工程上还有个进阶技巧:惯性权重线性递减——从 0.9 慢慢降到 0.4。前期惯性大,鸟群撒开网全局搜索;后期惯性小,大家稳稳地精修局部。这被称为"先广搜、后精修",是 PSO 最常用的改良。

四、上代码:70 行实现完整 PSO

下面这段代码零依赖,保存成 pso_demo.py 直接运行。我们拿两个经典函数开刀:Sphere(单峰,考收敛速度)和 Rastrigin(满地局部最优,考跳出陷阱的能力)。

  1. # -*- coding: utf-8 -*-- y, e$ _7 X' R6 ?% C6 c
  2. """粒子群优化算法(PSO)完整演示,零依赖可直接运行"""
      w1 z3 t! g3 c: v* n  K- Z
  3. import random* B7 F( k\" j\" T
  4. import math
    6 K5 n  D3 {+ P: w; a: ?
  5. ; B0 c- m* J9 l( z, ]& Z
  6. def sphere(x):; G# O: B+ ^5 H7 o& C! f3 w$ C
  7.     """Sphere 函数:f(0,0,...,0)=0,单峰"""0 d$ s, _7 k% R2 K. O/ ?
  8.     return sum(xi ** 2 for xi in x)
    & g* F. i# ^. u* i
  9. ! W! i8 T3 [. ]
  10. def rastrigin(x):
    \" N% i9 `; |4 F4 v! o# P' x
  11.     """Rastrigin 函数:大量局部最优点,全局最优仍在原点"""  H: R6 M+ N& _' S9 P
  12.     return sum(xi**2 - 10*math.cos(2*math.pi*xi) + 10 for xi in x)
    & ~' W& A' e; q0 s8 ?\" N  ], b' t$ h

  13. , I* z  s6 C0 Z- \7 r6 s) }
  14. class Particle:! p* c3 n\" B% I  {\" i# T
  15.     """一只『鸟』:有位置、速度,记得自己飞过的最好地方"""# R' [. l% U7 A
  16.     def __init__(self, dim, bounds):- `& d4 G/ K2 ]
  17.         self.pos = [random.uniform(bounds[i][0], bounds[i][1]) for i in range(dim)]; t6 F0 g7 c8 {5 t/ P
  18.         self.vel = [random.uniform(-1, 1) for _ in range(dim)]% S+ g& Y; }7 m4 b
  19.         self.best_pos = list(self.pos)      # 个体极值 pbest6 ?6 c! k/ Y\" k0 E2 L
  20.         self.best_val = float("inf")0 m. Y, n! x5 f3 p+ a+ y
  21. $ l) R, _* F6 R9 B1 D) i& k/ N
  22.     def update_velocity(self, gbest, w, c1, c2):
    & {+ }& m$ y  k8 V: O) O+ J5 s
  23.         # 速度更新:v = w*v + c1*r1*(pbest-x) + c2*r2*(gbest-x)
    ( g% K1 _5 h: U. s$ k( |8 G+ g
  24.         for i in range(len(self.vel)):, o) c% y( X5 b. {\" b7 o6 r
  25.             r1, r2 = random.random(), random.random(), Y+ c/ ^* t; X, w0 `
  26.             cognitive = c1 * r1 * (self.best_pos[i] - self.pos[i])  # 自己的经验5 ^- t( q1 F/ p7 h
  27.             social    = c2 * r2 * (gbest[i] - self.pos[i])          # 群体的情报
    + e( U0 {+ \7 n
  28.             self.vel[i] = w * self.vel[i] + cognitive + social) l* r7 u4 r+ O& f# T: h

  29. 4 Y4 D: Y4 L6 @9 K5 Z0 {
  30.     def update_position(self, bounds):% p: ?1 z. m/ h\" G* I% P
  31.         # 位置更新:x = x + v,并夹回边界* ~. N. N6 v% p  b* L! v  c
  32.         for i in range(len(self.pos)):
    ' [# y7 {! p& y2 V, v! R7 w, U
  33.             self.pos[i] += self.vel[i]1 w+ F* ~' z' I+ b
  34.             lo, hi = bounds[i]
    6 @9 p3 z5 y2 [$ L% C9 ^( \
  35.             if self.pos[i] < lo:   self.pos[i] = lo7 P/ @/ E+ L5 k- ~& |& Y
  36.             elif self.pos[i] > hi: self.pos[i] = hi
      D, o% X* y  d5 b8 f
  37. 7 r- c  ?) N1 l
  38. def pso(func, dim=2, bounds=None, n_particles=30, max_iter=100,$ N\" z0 I' z; Y( X8 F% i  b
  39.         w=0.7, c1=1.5, c2=1.5, verbose_every=10):
    * [$ L- O# l; `+ y  h- i
  40.     """标准 PSO 主流程""") s( D7 w- h. p2 u. L( Z
  41.     if bounds is None:  A9 h' F- d% n5 \3 u) I8 b) r% F
  42.         bounds = [(-5.12, 5.12)] * dim/ q; t0 z* h% g) o: F! Y
  43.     particles = [Particle(dim, bounds) for _ in range(n_particles)]
    ( b$ x  m* o! L% N9 Z0 Y& z
  44.     gbest, gbest_val = list(particles[0].pos), float("inf")
    ; r* v- V: o8 }
  45.     history = []
    ! L% |  |- f* W  Y$ K2 P

  46. 9 H3 {/ R. H) Q) Y% g\" |
  47.     for it in range(max_iter):
    ' f7 P* F' W1 e0 [4 M+ q. P, w
  48.         # 第 1 步:评估每只鸟,更新 pbest 和 gbest! X! j* s$ h* O& k) R& ~0 G
  49.         for p in particles:$ `2 M) t0 b' `) o( u
  50.             val = func(p.pos)
    2 k. C# b4 Y6 d
  51.             if val < p.best_val:& Q+ ~, H\" h& ?3 @/ b
  52.                 p.best_val, p.best_pos = val, list(p.pos)
    8 o6 Z6 t0 y/ ^
  53.             if val < gbest_val:
    , |\" w9 m4 T. u* V; e% S* ?
  54.                 gbest_val, gbest = val, list(p.pos)
    8 X& [6 o- z  B/ N! Q
  55.         history.append(gbest_val)
    ! }* D6 f+ O5 ?4 B% q

  56. \" f, V- ^* n$ A9 }# [2 r. x
  57.         # 第 2 步:更新速度和位置
    ! T7 T/ z( d% V, A0 F% w, t
  58.         for p in particles:
    / n! a$ Z( K3 Z
  59.             p.update_velocity(gbest, w, c1, c2)+ G0 D$ b5 b; E! [$ I9 O% a
  60.             p.update_position(bounds)
    ; V! n3 `, i% q8 D

  61. 2 s4 i# @/ X( D7 N5 c7 [
  62.         if verbose_every and (it + 1) % verbose_every == 0:
    \" D9 p( f9 ~# T0 a8 ^\" z- C
  63.             print(f"    第 {it+1:3d} 代 | 全局最优值 = {gbest_val:.6f}")
    $ b9 A4 \! b, b1 d) y+ G; @6 a9 w& n
  64.     return gbest, gbest_val, history
    # H\" C5 g8 r2 _3 _( Y4 o

  65. 9 X) C1 U0 |2 k1 s# [1 F. I# a# l
  66. if __name__ == "__main__":
    - {$ @\" b/ W: Y# k0 o
  67.     # 实验 1:单峰函数,看收敛速度
    1 j7 \2 ~\" ?5 U
  68.     print("[实验 1] Sphere 函数(理论最优值 = 0)")$ i; [2 {/ I5 j; ~
  69.     best_pos, best_val, _ = pso(sphere, dim=2, bounds=[(-10, 10)]*2,
    ) r& F1 ^3 E- g' v  A& d0 a) h
  70.                                 n_particles=30, max_iter=100)' B4 w* i) \7 F7 Z8 J: y
  71.     print(f"  >>> 最优位置: ({best_pos[0]:.4f}, {best_pos[1]:.4f})")( L/ _4 x! w0 E& M
  72.     print(f"  >>> 最优值: {best_val:.8f}")2 z9 B! r, ]7 U# f3 e
  73. 0 U% U$ w. W, b6 O6 h\" t
  74.     # 实验 2:多峰函数,看跳出局部最优, H# R' d% I7 {3 _2 \3 ^8 u
  75.     print("[实验 2] Rastrigin 函数(局部最优遍地都是)")
    ! g6 e4 m- J% V8 ]0 o( J
  76.     best_pos, best_val, _ = pso(rastrigin, dim=2,
    + d% O; z% [$ Q9 T
  77.                                 n_particles=50, max_iter=200)
    \" @$ I9 S# Z8 q  y, o8 y5 d
  78.     print(f"  >>> 最优值: {best_val:.6f}")
    4 W# f9 y6 ~; r& ?0 ~. c

  79. $ `$ D6 t; @* Q3 _$ k! H
  80.     # 实验 3:惯性权重对比(10 维)+ [; {7 t# B1 Q  C9 o
  81.     print("[实验 3] 惯性权重 w 的影响(10 维 Rastrigin)"). X$ S- ?' L8 o/ C7 F
  82.     for w in (0.2, 0.7, 1.2):
    5 Z2 t9 _\" ~# n. a+ W5 r
  83.         _, val, hist = pso(rastrigin, dim=10, bounds=[(-5.12, 5.12)]*10,
    3 m- f6 h; X+ a9 E
  84.                            n_particles=50, max_iter=200, w=w, verbose_every=0)
    4 G) ]8 W% K8 n$ V
  85.         print(f"  w = {w:<4} | 最终最优值 = {val:10.5f}")
复制代码

三个实验分别说明:

实验 1:单峰函数太简单,50 只鸟迭代 50 代就精确命中原点,误差小于百万分之一。

实验 2:Rastrigin 函数像个布满小坑的月球表面,梯度下降几乎必卡坑里,PSO 靠"一群鸟分散搜索"硬是跳了出来。

实验 3最有意思——升到 10 维后参数的差别立刻显现:w=1.2 时惯性太大,鸟群一直"飘着"落不下来(101.6,基本没搜到);w=0.2 时收敛太快,早熟地卡在局部最优(24.9);w=0.7 恰到好处(7.0)。这就是参数敏感性的现场教学。

五、什么时候用 PSO,什么时候别用适合的场景[td]
场景
举例
函数不可导 / 不光滑目标函数是查表、仿真得出的,没有解析式
多峰、怕陷入局部最优神经网络调参、TSP 路径优化、portfolio 投资组合
混合整数 / 离散问题稍加改造即可处理 0-1 变量(粒子位置取整)
需要"能跑就行"的快速方案数模竞赛限时三天,PSO 半小时就能调通
要小心的地方

1. 不是全局最优保证。 PSO 是启发式算法,理论上可能错过全局最优。论文里别写"证明了最优解",要写"搜索到的近似最优解"。

2. 高维会退化。 维度上百时收敛变慢,可考虑增加粒子数或换差分进化(DE)。

3. 结果有随机性。 每次运行结果略有不同。写论文时要多次运行取平均/最优,并报告标准差——评委很看重这一点。

4. 约束处理要自己写。 常见做法是罚函数:越界就往目标函数上加一个很大的惩罚。

和其他算法的关系[td]
算法
思想
和 PSO 的区别
梯度下降沿坡往下走单点出发、需可导;PSO 群体出发、不需要导数
遗传算法 GA选择、交叉、变异有"淘汰"机制;PSO 不淘汰,所有粒子活到最后、共享信息
模拟退火 SA以一定概率接受坏解单点搜索;PSO 是群体协作
差分进化 DE向量差扰动与 PSO 各有胜负,高维时 DE 常更稳
六、写给数模竞赛的你

回到竞赛场景。PSO 在数模里的出场率非常高,典型的用法:

题型一:优化决策。"如何分配 X 使总成本最小/总收益最大",目标函数复杂、变量多、不可导——PSO 直接上。

题型二:模型调参。 先用回归/机器学习建了模型,再有一堆超参数要调。可以用 PSO 包在外面做参数寻优,论文里就是个漂亮的"双层优化"结构。

题型三:拟合非线性模型。 最小二乘拟合非线性曲线时,初值选不好就拟合歪。用 PSO 先粗搜一轮找初值,再交给 LM 算法精修——这是很扎实的组合拳写法。

论文里怎么写才加分: 不要只丢一个结果。写清楚三件事——① 粒子数、迭代数、参数 w/c₁/c₂ 怎么定的(最好说明试过几组);② 收敛曲线(画"迭代代数 vs 最优值"的图,评委一眼看到收敛行为);③ 多次运行的结果稳定性(跑 10 次报告均值和标准差)。做到这三点,PSO 就从"套模型"升级成了"用模型"。

写在最后

PSO 之所以流传三十年,恰恰因为它简单得"不像话"——两条公式,没有 leaders,没有淘汰,只有一群朴素的个体,各自带着一点点记忆和一点点从众,互相分享情报,就把一片茫茫的搜索空间梳了个遍。

有人说这是对"群体智慧"最好的计算演绎:没有一只鸟知道答案,但鸟群知道。

下次遇到解不出的优化题,不妨也撒一把粒子出去——记得先广搜,后精修。


如果这篇文章帮到了你,点个「在看」让更多备赛的朋友看到


+ _! a2 j" M" U9 [- L
  a3 O: [$ U; [

9 X2 S1 P+ j2 B. G
! Q. J' W2 N5 F8 r: i7 i; T( k, O$ r, x8 a) s

: }7 `  Q& H% q! H9 X: L2 x% N
zan
转播转播0 分享淘帖0 分享分享0 收藏收藏0 支持支持0 反对反对0 微信微信
您需要登录后才可以回帖 登录 | 注册地址

qq
收缩
  • 电话咨询

  • 04714969085
fastpost

关于我们| 联系我们| 诚征英才| 对外合作| 产品服务| QQ

手机版|Archiver| |繁體中文 手机客户端  

蒙公网安备 15010502000194号

Powered by Discuz! X2.5   © 2001-2013 数学建模网-数学中国 ( 蒙ICP备14002410号-3 蒙BBS备-0002号 )     论坛法律顾问:王兆丰

GMT+8, 2026-8-23 18:27 , Processed in 0.504640 second(s), 52 queries .

回顶部