LEACH协议MATLAB仿真全解析：从原理到代码避坑指南

简介：面向无线传感器网络（WSN）研究者和相关课程学生，这份 LEACH 路由协议的 MATLAB 实现代码，可作为理解经典节能分簇算法和开展仿真实验的入门参考。LEACH（低能量自适应聚类层次）通过随机簇头选举与数据聚合来均衡节点能耗，资源包内的 leach.m 脚本覆盖节点分布、簇头选举、数据通信、能量模型等核心流程，便于快速运行并分析网络生命周期、能耗与数据成功率等指标。包体信息：压缩包含 1 个 m 文件，整体仅 2KB，轻量小巧，适合逐行阅读和二次改动。目前已有 687 人学习/下载，对初涉 WSN 仿真或准备毕业设计的读者具有一定参考价值。虽然源代码体量不大，但足以支撑跑通 LEACH 基线实验，借助脚本可调整节点数量、簇头概率和能量参数，观察不同策略对网络性能的影响，是课程设计、论文复现和协议优化的实用工具。 先说个真实经历：课程设计做到无线传感器网络，导师丢过来一句“你先把LEACH跑通”，你上网搜leach路由协议的MATLAB代码，能翻出几十个版本。挑了个星标最多、评论区一片“谢谢大佬”的下载下来，点运行，出图了——但那张存活节点曲线怎么看怎么不对劲。要么开局五十轮内节点成片暴毙，要么两百轮跑完还剩九十个节点“长生不老”，和论文里那种典型的阶梯状衰减曲线完全对不上。我当时就在这块折腾了将近一周，最后发现不是算法理解的问题，是代码里几个不起眼的坑在暗中作祟。这篇文章不打算再丢一份“能跑就行”的代码给你，而是把LEACH从协议机制到MATLAB实现逐层拆开，讲清楚每段关键代码在干什么、参数为什么这么定、网上流传版本里常见的坑到底在哪里。不管你是刚接触无线传感网的新手，还是正在做协议对比实验的研究生，这份拆解应该都能帮你省下不少时间。

1. LEACH协议到底干了什么，以及为什么非得用MATLAB

1.1 三件事：分簇、轮换、融合

LEACH全称Low Energy Adaptive Clustering Hierarchy，低功耗自适应分簇路由协议，它是无线传感器网络里最经典的分层路由协议之一。理解它的核心逻辑，一句话就能概括：别让每个节点都直接跟基站通信，选出一小部分节点当“簇头”，让普通节点把数据交给簇头，由簇头汇总后再转发给基站。

为什么这么做？因为无线传感器网络的节点靠电池供电，换电池不现实，而无线通信是最大的能耗来源。如果每个节点都直接往基站发数据，距离远的节点要不了多久就没电了，整个网络就会出现“感知空洞”。LEACH的思路用个生活类比特别容易懂：公司里如果所有人有事都直接找CEO汇报，CEO忙不过来，员工往返成本也高；更好的做法是每个部门选一个组长，统一收集意见再去跟CEO开会。组长这个活比较费精力，不能让同一个人一直干，所以要定期轮换，这回你当组长，下回换我来。LEACH的“轮”（round）就是干这个的。

具体机制分三步：

• 分簇：每一轮按照一个概率阈值选出一批簇头，其余节点根据距离远近加入某个簇头上报数据。
• 轮换：簇头不是固定的，每轮重新选举，把能耗负担分散到全网节点上。
• 融合：簇头收到簇内成员的数据后做数据融合，把多条消息压缩成一条再发给基站，大幅减少通信量。

这三个机制构成了LEACH协议的全部设计基础。后面你分析任何改进协议，比如LEACH-C、SEP、TEEN，本质上都是在回答同一个问题：我能不能把簇头选得更聪明一点、把能耗摊得更均匀一点。

1.2 选MATLAB做LEACH仿真，是有原因的

我记得有人问过：LEACH仿真不是有NS2、NS3、OMNeT++这些专用网络仿真器吗，为什么大家都用MATLAB？当然也可以用专业仿真器，但对绝大多数同学来说，MATLAB才是性价比最高的选择。NS3那套C++和TCL配置脚本的学习曲线，足够你再掉一层头发；而MATLAB的矩阵运算天然适合处理100个节点的网络状态数组，画图也方便，写完主循环加几行plot就能输出存活节点曲线、能量消耗曲线、数据吞吐量曲线。

更关键的一点是，MATLAB代码的逻辑是“透明”的。协议每一轮做了什么、每个节点当前什么状态、能量扣在哪一步，你都可以通过变量追踪看到。这对做课程设计、毕业设计的人来说价值极大——导师问起来，你能讲清楚每一行代码的理由，而不是给出一坨“仿真器黑盒”的运行结果。

2. 动手前先锁死能量模型和仿真参数

2.1 一阶无线模型：整个代码里的能量守恒基础

LEACH仿真里，所有能耗计算的依据是一个经典的一阶无线通信模型（first-order radio model），这个模型被绝大多数Leach论文沿用。它的规则非常直观：

• 节点发送k bit数据到距离为d的节点时，发射电路消耗的能量是 (E_{elec} \times k)，功率放大器的消耗取决于距离：当距离小于阈值 (d_0) 时，用自由空间模型，消耗 (E_{fs} \times k \times d^2)；当距离大于等于 (d_0) 时，改用多径衰落模型，消耗 (E_{mp} \times k \times d^4)。
• 节点接收k bit数据时，接收电路消耗的能量是 (E_{elec} \times k)。
• 簇头做数据融合时，每bit数据额外消耗 (E_{DA}) 的能量。

这个模型的直觉含义是：信号在短距离内传输衰减慢，能耗跟距离平方成正比；距离一旦超过某个临界值，信号衰减急剧恶化，能耗变为四次方关系。所以你必须先算出一个关键参数 (d_0)：

[ d_0 = \sqrt{\frac{E_{fs}}{E_{mp}}} ]

把下面参数表里的值代进去，(d_0) 大约是87.7米。在100m×100m的仿真区域内，有很大一部分节点与基站的距离会超过这个值，这直接决定了后续能耗计算的量级。这个细节特别重要——我在后面讲代码时会反复遇到它。

2.2 标准参数表：避免“复现不出来”的第一道坎

LEACH仿真里有一组约定俗成的默认参数。别看它们只是几行赋值，参数设得对不对，直接决定你的存活曲线跟别人的差多少。我把自己常用的参数整理成了一张表：

这里有个非常容易踩的坑：基站位置直接决定了你复现的结果是否跟原论文一致。基站放在(50, 175)时，绝大多数节点到基站的距离超过87.7米，走四次方能耗模型，网络整体寿命明显偏短；基站放在(50, 50)时，所有节点到基站的距离都在70米上下，基本走平方模型，网络能撑很长。不是说哪种设定“对”哪种“错”，而是你拿着别人的代码跑出来的曲线跟论文对不上时，第一件事不是怀疑代码有bug，而是先核对这两个参数。

3. LEACH核心算法代码逐行拆解

3.1 主循环：每轮到底做了哪几件事

LEACH的MATLAB仿真结构其实很清晰，主循环长这样：

for r = 1:rmax % 1. 统计当前存活节点 alive = find(node.renergy > 0); if isempty(alive) break; % 全部节点死亡，提前结束 end % 2. 簇头选举阶段 for i = alive' if node(i).G == 0 % 本回合还没当过簇头 T = threshold_calc(r, p, node(i).G); if rand < T && node(i).renergy > 0 node(i).CH = 1; % 当选簇头 clusterHeads(end+1) = i; %#ok<SAGROW> node(i).G = 1; % 标记本回合已当过簇头 end end end % 3. 成簇阶段：普通节点就近入簇 for i = alive' if node(i).CH == 0 mindist = inf; for ch = clusterHeads d = sqrt((node(i).x - node(ch).x)^2 + (node(i).y - node(ch).y)^2); if d < mindist mindist = d; node(i).MCH = ch; end end end end % 4. 数据发送与能耗结算 for i = alive' if node(i).CH == 1 % 簇头：接收簇内数据 + 融合 + 发送到基站 E_rx = (n_cluster_members) * packetLength * Eelec; E_agg = n_cluster_members * packetLength * EDA; d_to_bs = sqrt((node(i).x - BS(1))^2 + (node(i).y - BS(2))^2); E_tx = tx_energy(packetLength, d_to_bs); node(i).renergy = node(i).renergy - E_rx - E_agg - E_tx; else % 普通节点：发送数据到簇头 d_to_ch = sqrt((node(i).x - node(node(i).MCH).x)^2 + ... (node(i).y - node(node(i).MCH).y)^2); node(i).renergy = node(i).renergy - tx_energy(packetLength, d_to_ch); end end % 5. 每轮记录存活、能量、簇头数等统计量 end

主循环一共有四个阶段：选簇头、普通节点入簇、数据收发、能耗结算。很多初学者一上来就去读那些几百行的完整代码，结果被数组下标绕晕了。其实你把这个主循环骨架画出来，整个仿真就是个“往复循环”：每轮开始时选领导，选完领导分组，分完组干活，干完活算账，算完账进入下一轮。

3.2 阈值公式：LEACH选举的“标尺”长什么样

LEACH最有代表性的就是簇头选举阈值公式：

[ T(n) = \frac{p}{1 - p \times (r \bmod \frac{1}{p})} ]

前提是节点 (n) 在本回合（最近 (1/p) 轮）内还没有当过簇头，即 (G=0)。这个公式的数学意义在于：每轮期望的簇头数正好近似是 (n \times p)。当节点当了簇头后，(G) 置为1，在本回合剩余轮次内失去候选资格，这样其他节点才有机会轮上。

MATLAB实现代码如下：

function T = threshold_calc(r, p, G) if G == 0 T = p / (1 - p * (mod(r, round(1/p)))); else T = 0; end end

注意这里用了round(1/p)，因为p=0.1时1/p=10正好是整数；如果p取别的值，比如0.05，20也是整数；一旦p取0.12，1/p≈8.33，就必须用round或floor明确取整，否则mod的行为会跟你预期的不一致。这个细节我在下面第四部分会展开讲，它牵扯到一个挺隐蔽的bug。

3.3 成簇阶段：普通节点凭什么选这个簇头

普通节点的入簇策略很简单：计算自己到每个簇头的欧氏距离，选最近的簇头加入。MATLAB里没有现成的“找最近簇头”函数，所以一般写成两层循环：

% 用距离矩阵一次性算完，避免嵌套循环过慢 distToCH = zeros(length(alive), length(clusterHeads)); for a = 1:length(alive) for c = 1:length(clusterHeads) distToCH(a, c) = sqrt((node(alive(a)).x - node(clusterHeads(c)).x)^2 + ... (node(alive(a)).y - node(clusterHeads(c)).y)^2); end end [minDist, chIdx] = min(distToCH, [], 2);

这里有个性能小技巧：节点数在100这个量级时，双循环完全没问题；但如果扩展到大网络仿真，比如500个节点、每轮50个簇头，双循环就有点吃亏了。可以用MATLAB的pdist2函数一次性算所有节点两两之间的距离，再从中抽取需要的列，速度会快很多。代码是给人读的，更是给机器跑的，别在小规模仿真里过早优化，但也要知道有更高效的工具可用。

3.4 能量结算：扣错一笔，整条曲线就废了

能耗结算是整个仿真里最容易出错、却又最体现细节的部分。先定义发送能耗函数：

function E = tx_energy(k, d) global Eelec Efs Emp d0 = sqrt(Efs / Emp); if d < d0 E = k * Eelec + k * Efs * d^2; else E = k * Eelec + k * Emp * d^4; end end

这个函数判断距离是否小于临界值 (d_0)，用不同模型计算发送能耗。接收和融合能耗相对简单：

E_rx = packetLength * Eelec; E_agg = packetLength * EDA;

簇头的能耗模型是三重负担：接收所有成员的数据、融合数据、把融合结果发送到基站。普通节点只需要把自己的数据传输给簇头即可。每一轮算完，把消耗量从节点的剩余能量里减掉，剩余能量一旦小于等于0，就标记为死亡，后续轮次不再参与计算。

4. 网上流传代码里最致命的三个Bug

这部分是重头戏。我在帮别人调试LEACH代码的过程中发现，网上的版本虽然多，但问题高度集中，基本上是三个地方在反复出错。

4.1 Bug 1：阈值公式里的取模运算写错，导致所有节点疯狂竞选簇头

不少流传代码的阈值公式是这么写的：

T = p / (1 - p * (r * mod(1/p)));

看到问题了吗？r * mod(1/p)和mod(r, 1/p)完全是两回事。前者是 r 乘以一个常数，随着 r 不断增大，分母 (1 - p \times (r \times mod(1/p))) 会迅速变成负值，T 就变成负数或者超过1。当 T >= 1 时，rand < T永远成立，意味着每一轮每个节点都在竞选簇头。一轮下来上百个簇头，所有节点的能量在成簇和数据传输阶段就消耗殆尽，最直接的现象就是存活节点曲线在前二十轮内断崖式下跌。

正确的写法是用mod(r, round(1/p))，让分母保持在一个周期内循环，T 始终落在区间 (0, p] 附近。

4.2 Bug 2：G集合没维护好，簇头轮换形同虚设

G集合是LEACH实现能耗均衡的核心机制。一个节点当选簇头后，在本回合内应当从候选集合中移除，不再参与下一轮选举，直到 (r \bmod (1/p)) 重新归零。如果G集合不维护好，会出现什么现象？每一轮都是同一个高能量节点反复当选簇头，其他节点永远没机会，最终这个簇头因为能量耗尽而阵亡，网络很快就出现大面积覆盖空洞。

我见过一个版本的代码，它虽然判断了G==0，但每轮开头都把G重置为0，等于根本没限制：

for r = 1:rmax for i = 1:n node(i).G = 0; % 每轮强行重置，错！ end % 后面照常选举 end

这样的写法，逻辑上G形同虚设，簇头轮换完全乱套。正确做法是：在上一轮当选簇头时把G置1，并且只在 (mod(r, round(1/p)) == 0) 的那一轮，才把所有节点的G重置为0。

4.3 Bug 3：随机数种子处理不当，复现结果全靠“缘分”

MATLAB的rand函数每次运行都生成不同的随机序列，这本身不是问题，问题是做协议对比实验时需要可复现结果。很多开源代码把rand直接用得毫无控制，你跑三次，出三张完全不同的存活曲线，你根本没法判断算法的改进到底是有效的还是随机波动。

一个简单的做法是在仿真开头固定种子：

rng(2024); % 固定随机种子，保证结果可复现

但还有一个更隐蔽的问题：有些代码里，在轮询所有节点时对已死亡节点也调用了rand和距离计算，导致死亡的节点也在消耗随机序列。这样即使你设置了种子，不同版本或不同平台上得到的选举结果也可能不一致，因为死亡节点的布尔判断在数值上可能有微小差异。最好在节点初始化后就把所有随机事件统一规划好，或者在每次迭代之前都先筛选存活节点，让死节点不再参与任何随机计算。

这三点是我总结的“LEACH代码三大致命伤”。对照排查一遍，你会发现很多网上代码跑不出论文曲线的谜底，其实都在这里。

4.4 一个小问题的排查链路记录

我记得有一次调研一个改进协议，作者声称网络生命周期比LEACH延长了80%。我拿到代码之后直接跑，发现LEACH基准算法的首节点死亡轮数居然不到30轮，而作者论文里画的是150轮。起初我怀疑是能量参数设置问题，检查了一遍发现参数是对的。后来又怀疑是基站坐标不对，从一个版本换到另一个版本还是对不上。

最后把所有节点每轮选的簇头数打印出来，发现每轮簇头数量高达几十个，远远超过理论值 (n \times p = 10)。再追到阈值函数，果然又见到r * mod(1/p)这个写法。改回mod(r, round(1/p))之后，簇头数量基本稳定在10上下，首节点死亡轮数也恢复到了120轮左右。整个过程其实没有太多“灵光一现”，就是沿着数据合理性质疑、逐步定位的过程。

5. 怎么用数据和曲线判断协议好坏

5.1 三大生命周期指标：FND、HND、LND

仿真跑完之后，不能光看一张图。工程上通常用三个指标来量化网络生命周期：

LEACH的设计目标就是尽量推迟FND，因为第一个节点死亡往往意味着某个区域失去感知覆盖。如果你的改进算法只在LND上提升明显，FND却提前了，那这个改进方向可能把能耗负担集中到少数节点上去了，未必是好事。

用MATLAB统计这些指标其实很简单：

fnd = find(aliveHistory <= n-1, 1); % 第一个节点死亡的轮数 hnd = find(aliveHistory <= n/2, 1); % 一半节点死亡的轮数 lnd = find(aliveHistory == 0, 1); % 全部节点死亡的轮数

5.2 绘图：让数据自己说话

画图是MATLAB的主场。至少要把下面几张图画出来：

• 存活节点数随轮数的变化曲线：核心指标图，x轴是轮数，y轴是存活节点数。
• 网络总剩余能量随轮数的变化曲线：反映每个节点的平均能耗速度。
• 每轮簇头数柱状图：验证选举算法的稳定性，正常应该在p×n附近波动。
• 基站接收数据总量：体现网络的实际吞吐能力。

代码大概长这样：

figure; plot(1:rmax, aliveHistory, 'LineWidth', 1.5); xlabel('轮数 (round)'); ylabel('存活节点数'); grid on; figure; plot(1:rmax, totalEnergyHistory, 'LineWidth', 1.5); xlabel('轮数 (round)'); ylabel('网络总剩余能量 (J)'); grid on;

我看到太多人在做课程设计时，只给一张存活节点图，然后写一大堆文字描述。说实话，图不够，数据也不够。把能量曲线和簇头数波动都画出来，整个工作的说服力会高一个档次。

5.3 进阶方向：还想继续做的话，套路都在这里

如果你做完基础LEACH仿真，还想在这个方向继续深入，常见且省力的路线有这么几个：

• LEACH-C（集中式LEACH）：改成基站集中规划分簇，每轮基站根据节点剩余能量和位置算出最优簇头组合。对比指标就是FND和HND是否有提升。
• 节点能量异构场景：假设一部分节点初始能量高于其他节点，这时候LEACH的原始阈值公式是不是还合理？SEP协议就是干这个的。
• 多跳传输改进：簇头不再直接发数据到基站，而是在簇头之间选多跳路径，主要解决基站距离网络很远的场景。

做这些扩展时，基础仿真框架完全不需要重写，只需要改动能量计算和簇头选举两个函数即可。这也是当初用MATLAB写协议仿真的最大好处：模块化非常清晰，改一两个函数就能验证一个想法。

最后再分享一句经验：这类协议仿真代码，拿到手先别急着跑，花十分钟把阈值公式、G集合维护、随机数种子这三个地方看一遍，能帮你避开我在这个领域踩过的大半的坑。改完代码，再把参数跟原论文对齐一遍，这时候你跑出来的曲线，才真正具有可比性。这份“避坑手记”如果能让你的无线传感网仿真之路少绕几次弯，那今天这些内容就没白写。

本文还有配套的精品资源，点击获取