粒子群优化算法（Particle Swarm Optimization, PSO） 求解二维 Rastrigin 函数最小值问题

前言

提醒：
文章内容为方便作者自己后日复习与查阅而进行的书写与发布，其中引用内容都会使用链接表明出处（如有侵权问题，请及时联系）。
其中内容多为一次书写，缺少检查与订正，如有问题或其他拓展及意见建议，欢迎评论区讨论交流。

内容由AI辅助生成，仅经笔者审核整理，请甄别食用。

---

这段 MATLAB 代码实现的是 粒子群优化算法（Particle Swarm Optimization, PSO） 用于求解二维 Rastrigin 函数最小值问题。

我将从 整体框架、数学建模、公式推导和关键步骤 的角度详细分析此算法的原理和实现。

---

matlab代码

clc; clear; close all;% ------------------ PSO 参数设置 ------------------
num_particles = 30;        % 粒子数量
max_iter = 100;            % 最大迭代次数
dim = 2;                   % 问题维度
x_min = -5.12;             % 搜索空间下限
x_max = 5.12;              % 搜索空间上限
v_max = (x_max - x_min) / 2;w = 0.7;                   % 惯性权重
c1 = 1.5;                  % 个体学习因子
c2 = 1.5;                  % 社会学习因子% ------------------ 初始化 ------------------
% 初始化位置和速度
x = x_min + rand(num_particles, dim) * (x_max - x_min);
v = zeros(num_particles, dim);% 初始化个体最优
pbest = x;
pbest_val = arrayfun(@(i) rastrigin(x(i,:)), 1:num_particles)';% 初始化全局最优
[gbest_val, idx] = min(pbest_val);
gbest = pbest(idx, :);% ------------------ 迭代优化 ------------------
global_best_history = zeros(max_iter, 1);for iter = 1:max_iterfor i = 1:num_particles% 更新速度v(i,:) = w * v(i,:) ...+ c1 * rand() * (pbest(i,:) - x(i,:)) ...+ c2 * rand() * (gbest - x(i,:));% 限制速度v(i,:) = max(min(v(i,:), v_max), -v_max);% 更新位置x(i,:) = x(i,:) + v(i,:);x(i,:) = max(min(x(i,:), x_max), x_min);  % 限制位置在边界内% 计算适应度f_val = rastrigin(x(i,:));% 更新个体最优if f_val < pbest_val(i)pbest(i,:) = x(i,:);pbest_val(i) = f_val;end% 更新全局最优if f_val < gbest_valgbest = x(i,:);gbest_val = f_val;endendglobal_best_history(iter) = gbest_val;fprintf('迭代 %d：全局最优值 = %.6f\n', iter, gbest_val);
end% ------------------ 结果展示 ------------------
figure;
plot(global_best_history, 'LineWidth', 2);
xlabel('迭代次数'); ylabel('最优值');
title('PSO 优化 Rastrigin 函数过程');
grid on;fprintf('\n最终最优位置: (%.6f, %.6f)\n', gbest(1), gbest(2));
fprintf('函数值: %.6f\n', gbest_val);% ------------------ Rastrigin 函数 ------------------
function y = rastrigin(x)y = 20 + x(1)^2 + x(2)^2 - 10 * (cos(2 * pi * x(1)) + cos(2 * pi * x(2)));
end

运行结果

代码分析

🧠 一、问题建模

目标是优化 Rastrigin 函数，定义如下：

相关引用：Rastrigin函数简介

$$f(\mathbf{x})=10d+\sum _{i=1}^d\left[x_i^2-10\cos (2\pi x_i)\right]$$

对于本代码中 $d=2$，因此目标函数为：

$$f(x,y)=20+x^2+y^2-10\cos (2\pi x)-10\cos (2\pi y)$$

这是一个多峰函数，具有大量局部极小值，全局最小值在 $(0,0)$，其函数值为 0。

---

🚀 二、粒子群优化算法核心数学公式

粒子状态：

• 位置：${\mathbf{x}}_i^t\in {\mathbb{R}}^d$
• 速度：${\mathbf{v}}_i^t\in {\mathbb{R}}^d$

更新规则：

粒子速度与位置更新公式如下：

$${\mathbf{v}}_i^{t+1}=w\cdot {\mathbf{v}}_i^t+c_1\cdot r_1\cdot ({\mathbf{p}}_i-{\mathbf{x}}_i^t)+c_2\cdot r_2\cdot (\mathbf{g}-{\mathbf{x}}_i^t)$$

$${\mathbf{x}}_i^{t+1}={\mathbf{x}}_i^t+{\mathbf{v}}_i^{t+1}$$

其中：

• $w$：惯性权重，控制速度“保留性”
• $c_1,c_2$：个体学习因子、社会学习因子
• $r_1,r_2\sim \mathcal{U}(0,1)$：两个独立随机数
• ${\mathbf{p}}_i$：粒子 $i$ 的个体最优位置
• $\mathbf{g}$：所有粒子中的全局最优位置

---

⚙️ 三、代码分析

1. 参数设置

num_particles = 30;        % 粒子数量
max_iter = 100;            % 最大迭代次数
dim = 2;                   % 问题维度
x_min = -5.12;             % 搜索空间下限
x_max = 5.12;              % 搜索空间上限
v_max = (x_max - x_min) / 2;w = 0.7;                   % 惯性权重
c1 = 1.5;                  % 个体学习因子
c2 = 1.5;                  % 社会学习因子

设置粒子数量、最大迭代次数、搜索空间边界、速度最大值、惯性权重等。

2. 粒子初始化

x = x_min + rand(num_particles, dim) * (x_max - x_min); % 位置
v = zeros(num_particles, dim);                          % 速度

数学表示：

$${\mathbf{x}}_i^0=\mathcal{U}([x_{\min },x_{\max }{]}^d),{\mathbf{v}}_i^0=0$$

每个粒子的初始位置在边界范围内随机选取。

---

3. 适应度计算

相关引用：“适应度”简介

pbest_val = arrayfun(@(i) rastrigin(x(i,:)), 1:num_particles)';

pbest_val = arrayfun(@(i) rastrigin(x(i,:)), 1:num_particles)'

1. 1:num_particles 生成一个从1到num_particles的整数序列，表示粒子群的索引。例如，如果num_particles=30，则生成[1, 2, ..., 30]。
2. @(i) rastrigin(x(i,:)) 这是一个匿名函数（lambda函数），输入参数为i（粒子索引），输出为rastrigin(x(i,:))，即第i个粒子的位置向量x(i,:)在Rastrigin函数上的适应度值。
3. arrayfun arrayfun函数会对1:num_particles中的每个元素i执行匿名函数@(i) rastrigin(x(i,:))，并返回一个结果数组。例如：
   • 如果x是一个30×2的矩阵（30个粒子，每个粒子2维），则arrayfun会依次计算rastrigin(x(1,:)),
     rastrigin(x(2,:)), …, rastrigin(x(30,:))，并返回一个1×30的行向量。
4. '（转置符号） 由于arrayfun默认返回行向量，而pbest_val通常需要存储为列向量（便于后续操作），因此使用转置符号'将其转换为列向量。

即计算每个粒子的初始位置的函数值，并记录为个体最优。

数学表示：

$$f_i^0=f({\mathbf{x}}_i^0),{\mathbf{p}}_i={\mathbf{x}}_i^0$$

全局最优位置 $\mathbf{g}$ 是当前所有个体中函数值最小的位置。

---

4. 迭代优化主循环

for iter = 1:max_iter

- 更新速度（核心公式）：

v(i,:) = w * v(i,:) ...+ c1 * rand() * (pbest(i,:) - x(i,:)) ...+ c2 * rand() * (gbest - x(i,:));

数学形式（逐分量）：

$$v_{i,j}^{t+1}=w\cdot v_{i,j}^t+c_1\cdot r_1\cdot (p_{i,j}-x_{i,j}^t)+c_2\cdot r_2\cdot (g_j-x_{i,j}^t)$$

这是 PSO 的经典速度更新方程。

- 限制速度：

v(i,:) = max(min(v(i,:), v_max), -v_max);

防止粒子运动过快，跳出搜索空间。

- 更新位置：

x(i,:) = x(i,:) + v(i,:);
x(i,:) = max(min(x(i,:), x_max), x_min);

即：

$${\mathbf{x}}_i^{t+1}={\mathbf{x}}_i^t+{\mathbf{v}}_i^{t+1}$$

并保持在合法边界内。

- 更新个体最优和全局最优：

if f_val < pbest_val(i)pbest(i,:) = x(i,:);pbest_val(i) = f_val;
end
if f_val < gbest_valgbest = x(i,:);gbest_val = f_val;
end

反映出个体记忆和群体协作的机制。

---

📉 四、收敛性与可视化

global_best_history(iter) = gbest_val;
plot(global_best_history, ...);

记录并绘制每一代的全局最优值。

目标是找到函数的全局极小值点 $(x,y)=(0,0)$，函数值为 0。

---

✅ 五、总结与评价

组件	数学作用	代码实现
惯性项 $w\cdot v$	保持粒子动量	w * v(i,:)
认知项 $c_1\cdot r_1\cdot (p_i-x)$	自我学习	c1 * rand() * (pbest - x)
社会项 $c_2\cdot r_2\cdot (g-x)$	群体协作	c2 * rand() * (gbest - x)
边界约束	保证合法搜索	max(min(...))
适应度评估	目标函数值	rastrigin(x)

---