企业官网建设流程全解析

简介：面向无线传感器网络学习者的 LEACH 路由协议 MATLAB 实现代码，适合做课程实验、毕业设计或节能路由算法入门复现。LEACH 通过随机簇头选举、簇内数据聚合与轮换机制均衡节点负载，常用于观察网络寿命、能量消耗与数据传输效率。压缩包仅 2KB，内含 1 个 .m 文件，属于精简的 MATLAB 源码包，便于直接查看、运行和二次修改。已有 689 人学习。代码围绕节点随机分布、簇头概率选举、簇内与基站通信、能量消耗更新及多轮仿真展开，可帮助读者理解协议主流程，并在此基础上调整簇头比例、初始能量和网络规模，记录并绘制存活节点数、剩余能量等性能曲线。单文件结构也便于逐段调试、参数对照和结果复现，为算法对比、课程报告与论文实验提供可运行的基础脚本。

1. 先看结论：一份 LEACH 路由协议的 MATLAB 代码，能省掉你两周的踩坑时间

接手分簇路由这类课题时，最耗时间的通常不是公式推导，而是“把公式变成一条能跑出曲线的能量曲线”。这份 LEACH 路由协议的 MATLAB 代码解决的就是这件事：它把簇头选举、成簇、TDMA 时隙分配、稳态传输和一阶无线电能量模型串成一条可执行的链路，跑完一轮就能拿到节点剩余能量、簇头分布和死亡节点数。适合三类人：要快速复现 LEACH 基线做对比实验的研究生、要在自家协议前先搭一个可信参照系的工程师、以及被“代码跑得通但结果明显不对”折腾过的人。它不解决协议创新，但把基线做扎实，后面的对比才有意义。

2. 一阶无线电能量模型：把“耗电”换算成可对比的焦耳数

LEACH 之所以能画出那条经典的“轮次—存活节点”曲线，核心不在选举公式，而在能量模型足够诚实。换一套能量参数，同一份协议代码能给出完全相反的结论，所以动手写选举之前，先把这个模型固定下来。

2.1 为什么用一阶无线电模型而不是简单的 d² 衰减

一阶无线电模型（First Order Radio Model）把发送端能耗拆成两部分：射频电路消耗和功放消耗。功放部分按距离分成自由空间模型和多径衰落模型，分界点是阈值 d0。近距离下信号按 d² 衰减，远距离下按 d⁴ 衰减。

这个分界不是装饰。它直接决定了 LEACH 的一个关键行为：簇头到基站距离远时，簇头必须用 d⁴ 发送汇聚数据，这是簇头比普通节点死得快的物理原因，也是后续所有“能量均衡”改进（LEACH-C、SEP、DEEC）的出发点。

如果你图省事，全程用 d² 算发送能耗，会发现簇头寿命异常长、曲线漂亮得不像真的，然后你的改进算法在仿真里“提升有限”。这不是算法问题，是模型被你简化掉了关键非线性项。

常见参数取值如下，这份代码默认用的就是这一套：

d0 一定要用公式算，不要抄一个 87.7 进去。一旦你改 E_fs 或 E_mp 做敏感性分析，手填的 d0 会和公式脱节，出现“远距离反而比近距离省电”的诡异结果——这类翻车我在复现别人的对比实验时见过不止一次。

2.2 发送、接收、聚合三个函数的落地写法

把能量计算拆成三个独立函数，是整个工程里最值得的一个决定。后面做多跳、做中继、做异构节点，都只是换个地方调用它们，不用动主循环。

function E = tx_energy(k, d, prm) % TX_ENERGY 发送 k bit 到距离 d 处的能耗（J） % k : 数据包长度，单位 bit % d : 发送端到接收端的欧氏距离，单位 m % prm : 参数结构体，含 E_elec / E_fs / E_mp / d0 if d < prm.d0 % 自由空间模型：近距离，功放按 d^2 增长 E = k * prm.E_elec + k * prm.E_fs * d^2; else % 多径衰落模型：远距离，功放按 d^4 增长 E = k * prm.E_elec + k * prm.E_mp * d^4; end end function E = rx_energy(k, prm) % RX_ENERGY 接收 k bit 的能耗（J），与距离无关 E = k * prm.E_elec; end function E = da_energy(k, prm) % DA_ENERGY 簇头聚合 k bit 数据的额外能耗（J） E = k * prm.E_DA; end

逻辑说明：发送能耗的第一项是电路开销，与距离无关，只要发就得付；第二项是功放开销，随距离呈平方或四次方增长。接收端只付电路开销，因为接收不需要功放。聚合开销只发生在簇头，且和“聚合了几个包”成正比，不是和比特数成正比——严格写应该是E_DA * 成员数，这里按比特折算是一种近似，误差小但要在论文里说清楚。

参数说明：prm结构体统一管理所有物理量，主脚本里一次性初始化，后面所有函数只传引用式的结构体，避免十几个全局变量满天飞。

prm.E_elec = 50e-9; % 50 nJ/bit prm.E_fs = 10e-12; % 10 pJ/bit/m^2 prm.E_mp = 0.0013e-12; % 0.0013 pJ/bit/m^4 prm.E_DA = 5e-9; % 5 nJ/bit/signal prm.d0 = sqrt(prm.E_fs / prm.E_mp); % 约 87.7 m，由公式导出 prm.E0 = 0.5; % 初始能量 J prm.k_data = 4000; % 数据包 bit prm.k_ctrl = 200; % 控制包 bit

2.3 用单点校验确认模型没写错

模型写完别急着跑全流程，先做一组手算对照。取 k=4000 bit、d=50 m：

• 自由空间项：4000 × 50e-9 = 2.0e-4 J；功放项：4000 × 10e-12 × 2500 = 1.0e-4 J；合计 3.0e-4 J。
• 若 d=100 m（超过 d0）：功放项变成 4000 × 0.0013e-12 × 1e8 = 5.2e-4 J，合计 7.2e-4 J。

距离翻倍，能耗从 3.0e-4 涨到 7.2e-4，翻了 2.4 倍而不止 4 倍——这就是模型跨界切换的直接证据。写个五行的脚本对一遍，比跑完 1000 轮发现曲线不像 LEACH 再去回溯要省事得多。

% 手算对照：验证 d0 两侧的能耗跳变 for d = [50, 87, 88, 100] fprintf('d=%5.1f m, E_tx=%.4e J\n', d, tx_energy(4000, d, prm)); end

输出的四个值应当呈现“88 m 处斜率明显变陡”的形态。看不到这个拐点，说明prm.d0被覆盖过，或者tx_energy里的判断写成了<=之外的别的比较。这个自检花两分钟，能拦住后面一连串说不清的现象。

3. 簇头选举：把阈值公式落地成不会选出死节点的代码

能量模型定死了，接下来是 LEACH 的灵魂——每轮按概率 p 选簇头。公式在纸上只有一行，落到代码里有三处容易出错：G 集合的维护方式、轮次的取模边界、以及“能量耗尽节点必须排除”这条硬约束。

3.1 阈值公式里每个符号的真实含义

标准写法是：节点 n 在第 r 轮被选为簇头的概率为

T(n) = p / (1 − p × (r mod 1/p))，当 n ∈ G；否则 T(n) = 0。

p 是期望簇头比例，比如 0.05 表示每轮大约 5% 的节点当簇头。G 是“最近 1/p 轮内还没当过簇头”的节点集合。分母那一项的意义是补偿：随着轮次推进，G 集合里剩下的节点越来越少，分母变小，单节点的入选概率被抬高，从而保证每 1/p 轮每个节点恰好当过一次簇头。

很多复现代码把 G 集合写成mod(r, round(1/p))的派生判断，能跑，但有两个隐患：一是当节点因能量耗尽退出网络时，它仍被算进 G 集合的分母，导致剩余节点被系统性低估概率；二是浮点轮次或非整数 1/p 时取模结果不稳。稳妥做法是给每个节点挂一个lastCH记录“上次当簇头的轮次”，用轮次差判断是否在冷却期内。

3.2 带能量约束的选举实现

function chFlag = select_ch(E, lastCH, p, r) % SELECT_CH 按 LEACH 阈值公式选出本轮簇头 % E : 1×n 节点剩余能量（J） % lastCH : 1×n 每个节点上次当选簇头的轮次，初始为 -inf % p : 期望簇头比例 % r : 当前轮次 % chFlag : 1×n 逻辑数组，true 表示本轮为簇头 n = numel(E); chFlag = false(1, n); cool = round(1 / p); % 冷却窗口长度，p=0.05 时为 20 轮 for i = 1:n if E(i) <= 0 continue; % 能量耗尽，永久退出选举 end % 冷却期内的节点不进 G 集合 if (r - lastCH(i)) <= cool continue; end T = p / (1 - p * mod(r, cool)); if rand() < T chFlag(i) = true; end end end

逻辑说明：三层过滤依次生效——死节点排除、冷却期排除、剩余节点按阈值随机判定。冷却窗口cool用round(1/p)而不是1/p本身，是为了兼容 p=0.07、0.1 这类不能被 1 整除的取值，避免mod作用在浮点上产生1e-16量级的残差把判断带偏。

参数说明：p不要超过 0.1。p 太大簇头过多，聚合开销和簇头到基站的远距离传输会同时上升，整网能耗反而变差；p 太小则单簇成员过多，簇头聚合压力大、TDMA 帧变长。0.05 是文献里最常用的起点，做 p 的敏感性扫描时建议取 0.03 / 0.05 / 0.08 / 0.1 四档。

还有一点值得单独强调：选举完必须做一次兜底检查。

if ~any(chFlag) && any(E > 0) % 极端情况：本轮一个簇头都没选出来（概率低但会发生） alive = find(E > 0); chFlag(alive(randi(numel(alive)))) = true; end

没有这段兜底，小规模网络（比如 20 个节点、p=0.05）偶尔会出现“全轮无簇头”，此时所有节点保持睡眠状态、本轮零能耗，曲线上表现为一个不该存在的平段。这种偶发 bug 最难查，因为重跑一次就没了。

3.3 轮次推进与“每轮重置什么”

主循环里最容易搞混的是“哪些状态跨轮保留、哪些每轮清零”。给个明确的分工：

这张表看似琐碎，但它对应的是无数人栽过的坑：把chFlag定义在循环外，第二轮开始所有节点都还是簇头，曲线直接躺平。我一般会把网络状态封装成一个 struct，主循环里只做“复制 → 修改 → 存回”三步，减少这种变量作用域错误。

net = struct('pos', pos, 'E', E0 * ones(1, n), ... 'lastCH', -inf(1, n), 'dead', zeros(1, n)); for r = 1:rounds state = net; % 本轮工作副本 state.chFlag = select_ch(state.E, state.lastCH, p, r); state.lastCH(state.chFlag) = r; % 当选的更新冷却计时 [state, stat] = run_round(state, prm, bs); net = state; % 写回，E 已扣减 record(r, :) = stat; end

lastCH初始化为-inf而不是 0，是为了让第 1 轮所有节点都满足r - lastCH > cool，直接进入候选池。初始化为 0 时，第 1 轮会算出1 - 0 = 1 <= 20，全网没人有资格当选，然后靠兜底逻辑强选一个——能跑，但不是 LEACH 的行为。

4. 成簇、TDMA 时隙与稳态传输：一次完整轮次的执行链路

选举只决定“谁当簇头”，真正的能量消耗发生在成簇和稳态传输里。这一章把一轮拆成四个阶段：簇头广播、成员入簇、时隙分配、数据上传，每个阶段都对应不同的包长和能量扣减。

4.1 簇头广播与成员入簇：距离矩阵怎么算才不吃亏

节点数在 100 左右时，直接建距离矩阵比在双重循环里反复调norm快得多，而且代码可读性好：

function D = dist_matrix(pos, bs) % DIST_MATRIX 返回 (n+1)×(n+1) 距离矩阵，最后一行为基站 P = [pos, bs]; % 把基站拼在末尾 D = squareform(pdist(P')); % pdist 要求按列排，注意转置 end

pdist要求输入是“每列一个观测”，所以P'是必须的。少了这个转置，矩阵会变成把每个坐标分量当成一个节点，节点数直接翻倍，得到的距离矩阵维度对不上，后面索引全乱。这类错误不会报错，只会让结果莫名其妙。

簇头广播阶段每个簇头发一个控制包，扣tx_energy(k_ctrl, d)；非簇头收到所有簇头广播后，选择距离最近的那个入簇，扣若干次rx_energy(k_ctrl)。严格说，非簇头应该先收到全部广播再决定，所以接收次数等于簇头数：

chList = find(state.chFlag); memberOf = zeros(1, n); % 0 表示自己是簇头 for i = find(~state.chFlag & state.E > 0) [~, idx] = min(D(i, chList)); memberOf(i) = chList(idx); state.E(i) = state.E(i) - numel(chList) * rx_energy(prm.k_ctrl, prm); end

成员入簇时向簇头发一个控制包，扣一次发送能耗。这一步常被省略，因为占比小；但如果你要做“控制开销占比”这类分析，漏掉它就是错的。我的习惯是把所有控制包都算进去，宁可多算，不要事后被问“你这广播开销去哪了”。

4.2 TDMA 时隙分配与帧长计算

LEACH 的稳态阶段用 TDMA：簇头给自己簇内每个成员分配一个时隙，成员只能在自己的时隙发数据，其余时间休眠。时隙表本身不耗能（计算开销忽略），但它决定了稳态阶段的时长，进而影响“每轮时间”这个指标。

function slot = build_tdma(memberOf, chList) % BUILD_TDMA 为每个簇头生成成员时隙表 % 返回 cell 数组，第 k 个元素是第 k 个簇头的成员索引向量 slot = cell(1, numel(chList)); for k = 1:numel(chList) slot{k} = find(memberOf == chList(k)); end end

稳态阶段每个成员在各自时隙里发k_data比特给簇头，簇头接收并聚合。簇头自己也要把聚合后的数据发给基站——注意这一步的距离是“簇头到基站”，不是簇内距离，通常是整个网络里最长的一跳，d⁴ 项在这里发力。

for k = 1:numel(chList) c = chList(k); mem = slot{k}; if isempty(mem) && state.E(c) <= 0 continue; % 空簇且已死亡，跳过 end % 成员上行 for m = mem d = D(m, c); state.E(m) = state.E(m) - tx_energy(prm.k_data, d, prm); state.E(c) = state.E(c) - rx_energy(prm.k_data, prm); end % 簇头聚合 state.E(c) = state.E(c) - prm.E_DA * prm.k_data * max(numel(mem), 1); % 簇头向基站发送 dBS = D(c, end); state.E(c) = state.E(c) - tx_energy(prm.k_data, dBS, prm); end

这里有个细节值得注意：max(numel(mem), 1)保证空簇（成员全死）时仍有 1 个包的聚合开销，而不是 0。空簇意味着簇头自己也要上传自己的数据，聚合开销不该消失。

4.3 死亡判定与每轮统计量

能量扣完之后必须立刻判死，且判死要放在所有节点的能量都更新完之后，否则先判死的节点会影响后面对它的能量扣减：

function state = check_death(state) newlyDead = state.E <= 0 & state.dead == 0; state.dead(newlyDead) = 1; state.E(state.dead == 1) = 0; % 防止负能量累积 state.nAlive = sum(state.dead == 0); end

state.E(state.dead == 1) = 0这一行看着多余，实则是必须的。节点死亡后，如果后续某轮它仍然被当成簇成员参与能量扣减，能量会变成 -0.3 J 这样的负值；累加总剩余能量时这些负值会把总量抵消掉，画出来的“全网剩余能量”曲线出现反常抬升。把负能量截断为 0 是标准做法，代价是可忽略的数值误差。

每轮统计我一般记这几个量，后面画图和做对比全靠它们：

stat = [r, state.nAlive, sum(state.E), sum(state.chFlag), numel(chList)];

依次是轮次、存活节点数、全网剩余能量、簇头数、实际有成员的簇数。第三个和第五个常被忽略——簇头数和“有效簇数”不是一回事，当出现空簇时两者会分叉，这个差值本身就是诊断信息。

4.4 让整轮函数保持无副作用

把上面四段拼成run_round时，坚持一个原则：只通过入参和返回值通信，不读全局变量。好处是你可以随时拿某轮的中间状态单独重跑、对比、打日志，而不用重跑前面几百轮。

function [state, stat] = run_round(state, prm, bs) D = dist_matrix(state.pos, bs); chList = find(state.chFlag); % ……成簇、TDMA、稳态、判死 state = check_death(state); stat = [sum(state.nAlive), sum(state.E), numel(chList)]; end

这种写法在调试阶段的回报极高。有一次我发现第 340 轮之后曲线斜率突变，直接把该轮的状态存成 mat 文件，单独喂给run_round重跑十几遍，十秒定位到是空簇的聚合开销写成了 0。要是每次都得从第 1 轮跑起，这个 bug 能吃掉一整天。

5. 避坑与排查：五条踩过的血泪记录

仿真代码的特点是“不报错但结果不对”，下面五条都是我在复现和给别人看代码时反复遇到的，按出现频率排序。

5.1 曲线前几轮就掉节点

现象：第 3 到第 5 轮开始出现死亡节点，和文献里“前 200 轮全部存活”对不上。

原因：初始能量 E0 单位搞混。文献里常写 0.5 J 或 2 J，也有人写 500 mJ、写成 0.5 却按焦耳以外的量级在算；更隐蔽的是能量模型里的 nJ 没转成 J，E_elec = 50而不是50e-9，导致单次发送能耗是实际的一亿倍。

解决：跑之前先打印一次单轮总能耗，按 100 节点、5 个簇头估算应在 0.1～0.5 J 量级；如果打出来是 1e7 数量级，一定是单位没转。把prm初始化单独抽成一个函数，所有量都带e-9、e-12后缀，杜绝手写漏乘。

5.2 簇头数每轮都远大于预期

现象：设了 p=0.05、100 个节点，每轮却选出十几个甚至几十个簇头。

原因：阈值公式的分母写错，或者mod(r, cool)里 cool 用了1/p的浮点值导致分母接近 0，单节点概率被放大。另一种情况是冷却期判断写反（r - lastCH < cool当成>），每年节点每轮都重新参与竞争。

解决：在选举函数里加一行断言，跑完检查sum(chFlag)是否落在[0.5×n×p, 2×n×p]区间；超出就打日志。同时把cool强制取整，mod的两个入参都保证是整数。

5.3 存活曲线是一条直线下坠

现象：存活节点数从第一轮开始匀速下降，画出来是标准斜线，而不是文献里“先平后陡”的形态。

原因：通常是没有做能量耗尽节点退出选举的处理，或者判死后没有把E截断为 0。前者让死节点继续被选为簇头，后者让负能量在总和里抵消，两种情况都会把曲线压成直线。

解决：确认select_ch第一行有if E(i) <= 0, continue, end，确认check_death里有能量截断。再检查簇成员分配时是否还包含死节点——发现死节点仍在接收广播扣能量，说明find(~chFlag & E > 0)漏了E > 0这个条件。

5.4 同样的参数两次运行曲线不一样

现象：不改任何参数，连跑两次，首节点死亡轮次差了二三十轮。

原因：rand()没有固定种子。这不是 bug，但对比实验必须有可复现性，尤其当你要证明“改进算法比 LEACH 好 8%”时，随机波动可能就有 8%。

解决：每次实验开头rng(seed)，且不同算法用同一组种子。做多组实验取平均时，建议 20 组以上，种子固定成1:20，这样别人拿你的代码能跑出完全一致的图。

seeds = 1:20; res = zeros(numel(seeds), 3); % FND / HND / LND for s = 1:numel(seeds) rng(seeds(s)); res(s, :) = run_experiment(prm, p, rounds); end fprintf('FND = %.1f ± %.1f\n', mean(res(:,1)), std(res(:,1)));

5.5 簇头到基站的能耗看起来不合理地低

现象：簇头明明离基站 150 m，按 d⁴ 算应该有明显开销，但日志里簇头每轮只掉一点点能量。

原因：d0被手填成了 200 或者忘了用公式算，导致 150 m 仍走自由空间分支；或者距离矩阵里基站坐标拼错了位置，D(c, end)取到的其实是某个普通节点。

解决：打印prm.d0确认约等于 87.7；单测一次tx_energy(4000, 150, prm)，手算应为 4000×50e-9 + 4000×0.0013e-12×150⁴ ≈ 2.0e-4 + 2.63e-3 ≈ 2.83e-3 J，如果打出来是 1e-4 量级，就是走了自由空间分支。基站坐标建议在dist_matrix里最后一项拼入，并用size(D)断言等于n+1。

提示：这五条里有四条都能靠“单轮手算 + 断言”提前拦住，别等跑完 1000 轮再回头查。

6. 从单轮到多轮批量：指标落盘与三条值得试的改法

单轮跑通只是热身，真正拿去做对比要的是批量实验和可复现的指标。我通常固定三件事：种子、轮数、指标定义。指标就取三个——首节点死亡轮次（FND）、半数节点死亡轮次（HND）、最后节点死亡轮次（LND），再补一条全网剩余能量曲线。

function [fnd, hnd, lnd, eCurve] = run_experiment(prm, p, rounds, n) % 初始化网络并跑满 rounds 轮，返回三个死亡轮次与能量曲线 pos = 100 * rand(2, n); % 100m×100m 区域 bs = [50; 150]; % 基站位于区域外 E = prm.E0 * ones(1, n); lastCH = -inf(1, n); eCurve = zeros(1, rounds); fnd = rounds; hnd = rounds; lnd = rounds; for r = 1:rounds chFlag = select_ch(E, lastCH, p, r); if ~any(chFlag) && any(E > 0) alive = find(E > 0); chFlag(alive(randi(numel(alive)))) = true; end lastCH(chFlag) = r; state = struct('pos', pos, 'E', E, 'chFlag', chFlag, ... 'lastCH', lastCH, 'dead', zeros(1, n)); [state, ~] = run_round(state, prm, bs); E = state.E; nAlive = sum(state.dead == 0); eCurve(r) = sum(E); if fnd == rounds && nAlive < n, fnd = r; end if hnd == rounds && nAlive <= n/2, hnd = r; end if lnd == rounds && nAlive == 0, lnd = r; break; end end end

逻辑上要注意fnd的初值设为rounds而不是 0，这样“没死过节点”和“第 1 轮就死”能区分开，避免把未死亡误判成第 0 轮。hnd用<= n/2，节点数为奇数时向下取整，比== n/2稳。

参数上，rounds至少给到 LND 的两倍，否则最后一段曲线被截断，看不出全网死亡后的收敛形态。基站位置放在区域外（比如 y=150）是常见设定，它制造的远距离上行正是簇头负载集中的来源，也是后面所有改进型算法的着力点。

三条值得动手改的方向，按性价比排序：一是把簇头到基站改成多跳中继，选离基站更近的簇头转发，对 d⁴ 项的削减立竿见影；二是引入剩余能量因子，把阈值里的 p 换成 p×(E/E0)，让高能量节点多当簇头，这条改动只有一行但效果明显；三是给每轮加上真实的时隙时长，把“每轮时间”也统计出来，因为很多场景比的是“单位时间内的存活轮数”而不是纯粹的轮数。

最后说个习惯上的教训。早些年我做对比实验，图一出来发现改进算法只在某些种子下领先，换个种子就落后，当时第一反应是“改进不够显著”，差点把方向推翻。后来把种子固定成 20 组、把方差打进表格，才发现原方案在 20 组里只有 3 组领先，均值差被两组异常值撑起来的。从那以后我每次跑对比实验都强制走一遍“多种子 + 标准差 + 原始曲线截图留档”，不再看单次结果下结论。希望帮到你。

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