
简介面向无线传感器网络WSN研究者和相关课程学生这份 LEACH 路由协议的 MATLAB 实现代码可作为理解经典节能分簇算法和开展仿真实验的入门参考。LEACH低能量自适应聚类层次通过随机簇头选举与数据聚合来均衡节点能耗资源包内的 leach.m 脚本覆盖节点分布、簇头选举、数据通信、能量模型等核心流程便于快速运行并分析网络生命周期、能耗与数据成功率等指标。包体信息压缩包含 1 个 m 文件整体仅 2KB轻量小巧适合逐行阅读和二次改动。目前已有 687 人学习/下载对初涉 WSN 仿真或准备毕业设计的读者具有一定参考价值。虽然源代码体量不大但足以支撑跑通 LEACH 基线实验借助脚本可调整节点数量、簇头概率和能量参数观察不同策略对网络性能的影响是课程设计、论文复现和协议优化的实用工具。 先说个真实经历课程设计做到无线传感器网络导师丢过来一句“你先把LEACH跑通”你上网搜leach路由协议的MATLAB代码能翻出几十个版本。挑了个星标最多、评论区一片“谢谢大佬”的下载下来点运行出图了——但那张存活节点曲线怎么看怎么不对劲。要么开局五十轮内节点成片暴毙要么两百轮跑完还剩九十个节点“长生不老”和论文里那种典型的阶梯状衰减曲线完全对不上。我当时就在这块折腾了将近一周最后发现不是算法理解的问题是代码里几个不起眼的坑在暗中作祟。这篇文章不打算再丢一份“能跑就行”的代码给你而是把LEACH从协议机制到MATLAB实现逐层拆开讲清楚每段关键代码在干什么、参数为什么这么定、网上流传版本里常见的坑到底在哪里。不管你是刚接触无线传感网的新手还是正在做协议对比实验的研究生这份拆解应该都能帮你省下不少时间。1. LEACH协议到底干了什么以及为什么非得用MATLAB1.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仿真里有一组约定俗成的默认参数。别看它们只是几行赋值参数设得对不对直接决定你的存活曲线跟别人的差多少。我把自己常用的参数整理成了一张表参数符号取值含义节点总数n100部署在监测区域的传感器节点网络区域—100m × 100m节点随机均匀分布基站坐标BS(50, 175) 或 (50, 50)前者在区域外后者在区域内初始能量E00.5 J每个节点初始电池能量电路能耗Eelec50 nJ/bit发送/接收电路每bit能耗自由空间放大系数Efs10 pJ/bit/m²短距离传输放大能耗多径放大系数Emp0.0013 pJ/bit/m⁴长距离传输放大能耗数据融合能耗EDA5 nJ/bit簇头融合每bit数据能耗数据包大小packetLength4000 bit传感数据包控制包大小controlLength100 bit控制消息包簇头比例p0.1每轮期望成为簇头的节点比例这里有个非常容易踩的坑基站位置直接决定了你复现的结果是否跟原论文一致。基站放在(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(end1) i; %#okSAGROW 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) 轮内还没有当过簇头即 (G0)。这个公式的数学意义在于每轮期望的簇头数正好近似是 (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)因为p0.1时1/p10正好是整数如果p取别的值比如0.0520也是整数一旦p取0.121/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 2G集合没维护好簇头轮换形同虚设G集合是LEACH实现能耗均衡的核心机制。一个节点当选簇头后在本回合内应当从候选集合中移除不再参与下一轮选举直到 (r \bmod (1/p)) 重新归零。如果G集合不维护好会出现什么现象每一轮都是同一个高能量节点反复当选簇头其他节点永远没机会最终这个簇头因为能量耗尽而阵亡网络很快就出现大面积覆盖空洞。我见过一个版本的代码它虽然判断了G0但每轮开头都把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仿真跑完之后不能光看一张图。工程上通常用三个指标来量化网络生命周期指标全称含义为什么重要FNDFirst Node Dies第一个节点死亡的轮数网络是否还能提供完整覆盖的关键点HNDHalf Nodes Die一半节点死亡的轮数网络容量衰减到警戒线的时刻LNDLast Node Dies最后一个节点死亡的轮数网络整体寿命极限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集合维护、随机数种子这三个地方看一遍能帮你避开我在这个领域踩过的大半的坑。改完代码再把参数跟原论文对齐一遍这时候你跑出来的曲线才真正具有可比性。这份“避坑手记”如果能让你的无线传感网仿真之路少绕几次弯那今天这些内容就没白写。本文还有配套的精品资源点击获取