Python实现粒子群优化算法：原理、流程与代码详解

南乔

发布时间：2026-10-14 | 6人浏览过

原创

本文详细解析粒子群优化算法（PSO）的核心原理，包括个体极值与全局极值的更新机制。通过Python代码演示如何初始化种群、计算适应度、更新速度与位置，并展示求解二元函数最小值的完整流程与结果分析。

1.引言

粒子群优化算法灵感来源于对鸟类觅食行为的观察。在自然状态下，鸟群在寻找食物时会表现出聚集、分散及轨迹调整等行为，同时个体间能保持默契的距离。该算法模拟这一过程，将待求解问题的搜索空间映射为鸟类飞行的空间，并将每个可行解抽象为一个无质量、无体积的粒子。因此，寻找最优解的过程即等同于模拟鸟类觅食的过程。

粒子群算法基于种群和进化的概念，通过个体间的竞争与协作来求解复杂空间中的最优解。与遗传算法不同，它不对个体执行“交叉”或“变异”操作，而是依据特定规则更新每个粒子的速度和位置，使粒子向自身历史最佳位置和全局历史最佳位置移动，从而推动整个种群向最优方向进化。

2.算法的具体描述：

2.1原理

在粒子群优化算法中，粒子通过信息共享机制获取其他粒子的发现与飞行经历。这种信息共享机制本质上是一种合作共生行为。在搜索最优解时，每个粒子既能记忆自身经过的最佳历史位置，又受群体中其他粒子的影响。因此，粒子的行为既受其他粒子的引导，也受自身经验的指导。

粒子群算法对鸟群的模拟遵循以下模式：假设一群鸟在空中搜索食物，每只鸟都知道当前距离食物的远近（由适应度值衡量）。最简单的搜索策略是飞向距离当前食物最近的鸟的周围空间。在粒子群算法中，每个粒子对应一只鸟，拥有适应度值和决定飞行距离与方向的速度。所有粒子追随当前最优粒子在解空间中搜索。每次搜索后，最优粒子可能变化，其他粒子随即追随新的最优粒子，如此反复迭代。

迭代开始时，每个粒子通过随机方式初始化其在空间中的速度和位置。在迭代过程中，粒子通过跟踪两个极值来更新自身在解空间中的位置和速度：一个是单个粒子在迭代过程中的最优位置（对应最优适应度值的空间解），称为个体极值；另一个是种群中所有粒子在迭代过程中找到的最优位置，称为全局极值。若粒子仅跟踪一个极值，则算法分别称为局部粒子群算法或全局粒子群算法。

PSO从该模型中获得启示并用于解决优化问题。在PSO中，每个优化问题的潜在解都是搜索空间中的一只鸟，称为粒子。所有粒子都有一个由被优化函数决定的适值（fitness value），每个粒子还有一个速度，决定其飞翔的方向和距离。随后，粒子们追随当前最优粒子在解空间中搜索。

PSO初始化为一群随机粒子（随机解），然后通过迭代找到最优解。在每次迭代中，粒子通过跟踪两个极值来更新自己：第一个是粒子本身找到的最优解，称为个体极值；另一个极值是整个种群目前找到的最优解，称为全局极值。此外，也可以不使用整个种群，而仅使用其中一部分作为粒子的邻居，此时所有邻居中的极值即为局部极值。

图解：

2.2标准粒子群算法流程

算法的流程如下：

Step1:种群初始化:可以进行随机初始化或根据被优化的问题设计特定的初始化方法,包括群体规模，每个粒子的位置 X i X_{i} Xi 和速度 V i V_i Vi ，然后计算每个粒子的适应度值，从而选择出个体的局部最优位置向量和种群的全局最优位置向量。

Step2:迭代设置:设置迭代次数 g m a x g_{max} gmax ，令当前迭代次数g=1。

Step3:根据公式更新每个粒子的速度向量V。

Step4:根据公式更新每个粒子的位置向量X。

Step5:局部位置向量和全局位置向量更新:更新每个粒子的Pbest,和种群的Gbest。

Step6:终止条件判断:判断迭代次数时都达到 g m a x g_{max} gmax 或误差是否足够小，如果满足则输出Gbest.否则继续进行迭代，跳转至步骤(3)。

对于粒子群优化算法的实际应用，因为主要是对速度和位置向量迭代算子的设计计,选代算子是否合理将决定整个PSO算法性能的优劣.，所以如何设计 t pso的迭代算子是算法应用的研究重点和难点。

3.代码案例

3.1问题

求解f(x,y)的最小值点

3.2绘图


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)

# 目标函数
Z=X**2+Y**2+X

# 绘图
fig=plt.figure()
ax=Axes3D(fig)
surf=ax.plot_surface(X,Y,Z,cmap=cm.coolwarm)
plt.show()

3.3计算适应度


# 速度
# 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维的值

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

mpl.rcParams['font.sans-serif'] = ['SimHei']  # 指定默认字体
mpl.rcParams['axes.unicode_minus'] = False  # 解决保存图像是负号'-'显示为方块的问题

def fitness_func(X):
    """计算粒子的的适应度值，也就是目标函数值，X 的维度是 size * 2 """
    A = 10
    pi = np.pi
    x = X[:, 0]
    y = X[:, 1]
    return x**2+y**2+x

3.4更新速度


def velocity_update(V, X, pbest, gbest, c1, c2, w, max_val):
    """
    根据速度更新公式更新每个粒子的速度
    :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))
    r2 = np.random.random((size, 1))
    V = w*V+c1*r1*(pbest-X)+c2*r2*(gbest-X)
    # 防止越界处理
    V[V < -max_val] = -max_val
    V[V > -max_val] = max_val
    return V

3.5更新粒子位置


def position_update(X, V):
    """
    根据公式更新粒子的位置
    :param X: 粒子当前的位置矩阵，维度是 20*2
    :param V: 粒子当前的速度举着，维度是 20*2
    """
    return X+V

3.6主要算法过程


def pos():
    w = 1
    c1 = 2
    c2 = 2
    r1 = None
    r2 = None
    dim = 2
    size = 20
    iter_num = 1000
    max_val = 0.5
    best_fitness = float(9e10)
    fitness_val_list = []
    # 初始化种群各个粒子的位置
    X = np.random.uniform(-5, 5, size=(size, dim))
    # 初始化各个粒子的速度
    V = np.random.uniform(-0.5, 0.5, size=(size, dim))
    # print(X)
    p_fitness = fitness_func(X)
    g_fitness = p_fitness.min()
    fitness_val_list.append(g_fitness)

    # 初始化的个体最优位置和种群最优位置
    pbest = X
    gbest = 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
            # 记录最优迭代记录
            fitness_val_list.append(g_fitness)
            i += 1

    # 输出迭代结果
    print("最优值是：%.5f" % fitness_val_list[-1])
    print("最优解是：x=%.5f,y=%.5f" % (gbest[0], gbest[1]))

    # 绘图
    plt.plot(fitness_val_list, color='r')
    plt.title('迭代过程')
    plt.show()

pos()