做通信系统仿真的朋友应该都清楚LDPC码的性能很大程度上押在稀疏校验矩阵上。最近我完成了一个用大衍数构造稀疏校验矩阵的LDPC误码率Matlab仿真工程对比了不同译码迭代次数、码率和码长对误码率曲线的影响。整套代码能直接跑改参数就能出图对研究校验矩阵构造、学习LDPC译码的同学很实用。这篇博客把从构造原理到仿真落地的过程完整梳理一遍尽量少讲虚的全是实操。1. 项目需求拆解与整体仿真方案设计1.1 这个仿真到底在做什么LDPC仿真很容易陷入“跑出两条BER曲线就算了”的误区。这个项目的要求拆开看其实只有三个层次构造出确定性的稀疏校验矩阵不能靠随机碰运气。实现BP译码并让迭代次数、码率、码长都成为可配置参数。用蒙特卡洛仿真统计误码率用曲线支撑结论。我一开始也犹豫过是否直接用dv、dc随机构造但随机矩阵容易出现短环尤其是在中短码长下误码率平台很突出。后来从数论里借鉴大衍数的方法用同余方程生成非零元位置矩阵结构确定四环可控跑出来的曲线稳定很多。这个仿真的目标不是追求极致性能而是把“矩阵构造—译码—性能评估”这条链路打通并定量观察三个参数的影响。1.2 为什么选择大衍数构造方法在通信领域校验矩阵构造大体分随机、代数和结构三类。Gallager的随机构造最简单但每次生成的矩阵环分布不同有限几何构造性能好但参数不灵活。大衍数构造本质上是准循环LDPC的一个变种核心思想是用模逆元确定循环移位量。好处是H矩阵由基础矩阵和扩展因子唯一确定不需要存储大矩阵。结构确定利于硬件实现码长码率扩展也方便。使用大衍求一术求逆元参数选取有数论依据不是拍脑袋。关于大衍数很多朋友第一反应是《周易》的“大衍之数五十其用四十有九”。但在工程里我主要用到的是秦九韶《数书九章》中发展出的大衍求一术也就是求解一次同余式组的系统性方法。为构造LDPC矩阵我把扩展因子z作为模数用大衍求一术求解同余方程得到一组与z互素的乘性逆元再按行、列坐标生成循环移位偏移量。这样生成的矩阵具有准循环结构同时可以控制环长。1.3 仿真对比维度与指标设计项目要求同时观察译码迭代次数、码率、码长那就要把对比维度设计清楚否则蒙特卡洛仿真时间会爆炸。我最后确定的设计如下表对比维度取值固定条件译码迭代次数5、10、20、30码率1/2码长648码率1/2、2/3、3/4迭代次数20码长648码长324、648、1296码率1/2迭代次数20误码率统计采用误比特率BER同时记录误帧率FER因为LDPC译码失败往往整帧错。每组信噪比点至少统计100帧错误低误码率时限制最大帧数10万帧保证曲线平滑。这个设计看似简单但跑起来就知道码长1296、码率3/4、10万帧在普通电脑上要几个小时。后面我会讲怎么用Matlab的并行工具加快扫描。2. LDPC编译码与大衍数构造原理2.1 LDPC码稀疏校验矩阵的直观理解LDPC码的核心是用一个稀疏矩阵H描述码字约束关系。所谓稀疏是非零元素远少于零元素。比如648码长、1/2码率的规则LDPC码如果列重是3那么H矩阵非零元素只有648×31944个而总元素是324×648209952个稀疏度只有0.9%。译码时信息通过Tanner图的边在变量节点和校验节点间传递短环会破坏消息独立性所以构造时最怕出现长度为4的环。行列重也很关键。列重太小最小距离低误码平台高列重太大译码门限变差。1/2码率规则码常用列重3、行重6。用大衍数构造时通过基础矩阵的设计可以控制行列重。以码率1/2、码长648为例取扩展因子z54则基础矩阵维度mb6、nb12列重设为3行重自然就是6。2.2 大衍数构造稀疏校验矩阵的算法步骤我实现的构造方法可以拆成四步。第一步确定基础矩阵维度。假设最终H矩阵是m行n列扩展因子为z那么基础矩阵维度是mbm/z行、nbn/z列。为了便于编码通常要求mb远小于nb而且m取z的整数倍。第二步构造基础矩阵的非零图案。对于规则码每个基础矩阵行有dc个1列有dv个1。我用的是确定性循环放置方式而不是随机放置。第j列的非零行索引由下式确定row(k) mod(j - 1 s * (k - 1), mb) 1, k 1,2,...,dv其中s是一个与mb互素的步长。这个步长也由大衍数候选序列中选取从而避免重边问题。第三步用大衍求一术生成循环移位系数。对基础矩阵中每个非零位置(i,j)需要计算一个0到z-1的偏移p(i,j)。具体做法是取一组与z互素的整数序列a然后对每个非零位置求解同余方程a(i) * p(i,j) ≡ (j - c(i)) (mod z)其中的c(i)是行方向偏移由行号i和“其用四十有九”的49取模确定。由于a(i)与z互素这个同余方程有唯一解解出的p(i,j)就用做大衍数偏移量。大衍求一术本质上就是扩展欧几里得算法Matlab里可以直接用gcd函数实现。第四步扩展。把基础矩阵中每个“1”替换成z×z单位阵循环右移p位的矩阵每个“0”替换成全零矩阵得到最终的H矩阵。这个过程在Matlab里可以用稀疏索引构造避免生成完整稠密矩阵。这里还要提一句“问数化定数”。秦九韶在大衍术里把非两两互素的问数转化为两两互素的定数我在构造时对扩展因子z做了类似处理。如果z本身是合数就先把z拆成互素因子的乘积再分别对每个因子求逆最后用中国剩余定理合并偏移量。为了简化首版仿真我直接选择了与所有候选a都互素的z代码更简单后续再扩展这个处理也不难。2.3 校验矩阵质量检查构造完H矩阵后我在代码里强制做了三项检查稀疏度检查非零元占比是否在预期范围。行重列重统计是否满足规则码设定。四环检查遍历任意两列计算公共非零行数量如果大于1就存在四环。四环检查的Matlab实现非常简单function has4 checkCycle4(H) colIdx cell(size(H,2), 1); for j 1:size(H,2) colIdx{j} find(H(:,j)); end has4 false; for j1 1:size(H,2)-1 for j2 j11:size(H,2) if length(intersect(colIdx{j1}, colIdx{j2})) 1 has4 true; return; end end end end如果检查出四环我会返回第一步调整基础矩阵的步长s或换一组a序列。实测下来用大衍数方法比纯随机生成命中无四环矩阵的概率高很多。3. Matlab代码实现与参数配置3.1 主程序结构整个工程的主脚本流程是清空环境、加载参数。调用constructHDayan函数生成H矩阵。求出生成矩阵G完成编码。对每个信噪比点做蒙特卡洛仿真生成随机信息位编码BPSK映射加高斯白噪声LLR-BP译码统计错误。保存BER/FER结果并画图。主脚本的核心参数区块我习惯集中放在文件头部% 仿真参数 EbN0dB 0:0.5:8; z 54; % 扩展因子 mb 6; nb 12; % 基础矩阵维度码率约1/2 R 1 - mb/nb; maxIter 20; % 最大译码迭代次数 numFrames 1e5;这里码率是1/2。如果要做码率对比就换mb和nb组合例如mb4、nb12得到2/3码率mb3、nb12得到3/4码率。码长通过z同步调整z27、54、108分别对应码长324、648、1296。参数之间的换算关系我建议写进脚本注释改起来不容易乱。3.2 构造H矩阵的Matlab实现下面给出constructHDayan函数的核心部分细节做了注释function [H, Hbase] constructHDayan(mb, nb, z, dv) % 第一步确定性基础矩阵非零图案 Hbase zeros(mb, nb); s 1; % 如果与mb不互素换下一个候选 while gcd(s, mb) ~ 1 s s 1; end for col 1:nb for k 1:dv row mod(col - 1 s * (k - 1), mb) 1; Hbase(row, col) 1; end end % 第二步准备与z互素的大衍数候选基数 candidates mod(50 (0:48), z); candidates unique(candidates); candidates(candidates 0) []; aSeq candidates(gcd(candidates, z) 1); if isempty(aSeq) error(扩展因子与候选数不互素请调整z); end % 第三步大衍求一术求模逆元并扩展 H sparse(mb*z, nb*z); cnt 0; for i 1:mb for j 1:nb if Hbase(i,j) 1 cnt cnt 1; a aSeq(mod(cnt-1, length(aSeq)) 1); c mod(49 * i j, z); [g, x, ~] gcd(a, z); if g ~ 1 error(求逆失败); end p mod(x * (j - c), z); % 将单位阵循环右移p位后的非零位置填入H for r 1:z colPos mod(r - 1 p, z) 1; H((i-1)*z r, (j-1)*z colPos) 1; end end end end end这个函数有几个关键点。循环移位方向要统一否则译码端需要额外调整。使用sparse预分配大矩阵能避免内存爆炸。基础矩阵如果有重边需要额外处理我通过步长s与mb互素来避免重边。3.3 译码器实现与迭代次数参数接口译码器用的是对数似然比置信传播也就是SPA。输入是信道软信息LLR输出是硬判决码字。核心迭代过程包括校验节点更新和变量节点更新。为了方便对比迭代次数我把最大迭代次数作为函数参数传入function [xHat, iterUsed] decodeLDPC_BP(H, rxLLR, maxIter) [M, N] size(H); [rowIdx, colIdx] find(H); Lv rxLLR(:).; Lvc zeros(length(rowIdx), 1); % 变量到校验的消息 Lcv zeros(length(rowIdx), 1); % 校验到变量的消息 for iter 1:maxIter % 变量节点更新 % 校验节点更新这里使用tanh规则 % 计算伪后验概率 % 硬判决 iterUsed iter; % 如果校验方程全零提前退出 end endSPA的校验节点更新是最耗时的部分我推荐用tanh规则或者Min-Sum近似。本项目为了精确保真使用tanh但要对数值稳定性做保护比如限制LLR绝对值不超过30。迭代次数接口就是maxIter在扫描脚本里直接赋值。3.4 参数扫描脚本一次性跑三个维度脚本边界要清楚。我专门写了runSweep.m用struct数组存结果results struct(iter, {}, rate, {}, len, {}, EbN0dB, {}, BER, {}); % 迭代次数扫描 for it [5 10 20 30] ber runBER(H, R, z, EbN0dB, maxIterit); results(end1) struct(...); end这里用到了Matlab较新版本的命名参数语法如果版本比较老就改成普通参数传入。每次仿真结束把BER数据保存成mat文件防止电脑崩溃后白跑。4. 仿真结果对比与影响分析4.1 译码迭代次数对误码率的影响先看固定码率1/2、码长648条件下迭代次数从5次增加到30次的BER曲线。结果符合直觉迭代次数从5次提升到10次时性能提升非常明显高信噪比区域误码率能下降一个数量级以上从10次到20次还有可感知的改善再往上到30次曲线几乎没有变化。这说明在中等信噪比下BP译码的收敛速度存在饱和。迭代次数增加带来的是译码时延和功耗线性上升性能收益却递减。实际系统选取迭代次数不能只看误码率曲线还要看吞吐率约束。我习惯在仿真报告里同时给出平均实际迭代次数的曲线它会帮你在“性能天花板”和“复杂度”之间找到平衡点。4.2 码率对系统性能与吞吐的影响码率直接决定冗余度。仿真的三组码率1/2、2/3、3/4中1/2码率曲线最好3/4码率曲线最差。以误码率1e-4为参考2/3码率大概比1/2码率差0.8 dB左右3/4码率差更多。但码率越高有用信息占比越多频谱效率越高。因此实际通信系统选码率本质是误码性能和吞吐性能的折中。自适应调制编码中信道条件差时降到1/2码率信道条件好时切到3/4靠的就是这套仿真数据支撑门限配置。如果把仿真结果换算成频谱效率能更直观地看到码率选择的意义。4.3 码长对误码率性能与复杂度的双向影响码长324、648、1296对比码长越长性能越好这是因为长码的随机化程度好更接近香农限。1296码长在误码率1e-5处比324码长有大约1.2dB的编码增益差异。但长码的代价不只是仿真时间。译码器校验矩阵尺寸变大BP的每一轮迭代复杂度都随之上升存储LLR消息的内存也线性增长。从工程影响范围看LDPC码在高速传输场景中采用较大码长做高速率传输而在物联网短包场景中码长受限就需要用短码结合适当迭代次数来控制时延。这套仿真正好提供了三个关键维度对性能影响的量化结果。5. 常见问题与调优经验5.1 校验矩阵构造中的短环与重边处理用大衍数方法构造时即便有数论基础基础矩阵设计不好仍然会出现四环。我踩过的一个坑是扩展因子z如果与步长s不互素会导致同一列中多个非零行冲突矩阵出现重边。解决办法是检查gcd(a,z)1并且对每个扩展块做非零行列检测。如果真的检测到四环不建议直接重开随机矩阵更快的方法是保持既有非零图案只对冲突的偏移量p(i,j)加一个固定的模增量再做一次求逆。我试过通常调整两三个位置就能消除四环。5.2 编码端G矩阵过慢或内存爆炸用H矩阵直接高斯消元求G矩阵在码长648以后会非常慢而且G矩阵是稠密的1/2码率648码长就要648×324个double内存还好但1296码长就很吃紧。我的处理是用稀疏LU分解解H*x0的零空间或者直接把H构造成近似下三角形式做迭代编码。这样虽然编码有少量开销但整体仿真时间下降非常明显。如果你对编码时延不敏感更简单的方式是用Matlab内置的ldpcEncoderCfg等函数但那样会把矩阵构造过程隔离掉不利于研究校验矩阵本身。所以我这里还是保留了自己求G的流程。5.3 Matlab仿真性能优化与并行扫描蒙特卡洛仿真最忌讳的是在循环里反复找H的非零索引。我的经验是进入信噪比循环之前把所有非零边的行列索引提取出来存成数组BP译码时直接索引。另外LLR初始化用向量化运算不要用for逐比特。信噪比点之间相互独立用parfor并行跑四核机器能省2/3时间。记得先对每个worker预先加载H矩阵否则内存缓存反复传递性能反而下降。还要注意随机数流控制。并行仿真时如果不设置每个worker的随机种子不同信噪比点可能产生相同的随机序列BER曲线不独立。我用的是parfor循环内根据信噪比索引手动设置RandStream确保每个worker的数据独立。5.4 这套仿真工程的扩展方向这个工程不仅限于对比三个参数。把译码器换成Min-Sum或者NMS可以研究低复杂度算法损失把H矩阵构造部分改成5G标准里定义的QC-LDPC可以直接评估标准码的性能加入BICM与高阶调制也能分析编码调制联合方案。大衍数构造核心价值在于参数化与可复现性给大家提供了一个在Matlab里快速验证新想法的底座。最后讲一点我做完这个工程后的体会。开始时我也迷信“迭代次数越多越好”真正把迭代次数从5扫到30才发现中短码长下20次以后几乎是原地踏步。还有码率选择不要只看编码增益要结合频谱效率和使用场景。仿真代码的价值不在于把曲线画得多漂亮而在于每根曲线背后都能回答一个工程问题。这套程序我大概测了一周中间被矩阵奇异、内存爆掉各种折磨但调完之后再去看5G的LDPC参数就顺畅多了。如果你也在做类似对比建议先小码长小迭代数跑通流程再放大参数准能少走弯路。
基于大衍数构造稀疏校验矩阵的LDPC码误码率仿真实现
做通信系统仿真的朋友应该都清楚LDPC码的性能很大程度上押在稀疏校验矩阵上。最近我完成了一个用大衍数构造稀疏校验矩阵的LDPC误码率Matlab仿真工程对比了不同译码迭代次数、码率和码长对误码率曲线的影响。整套代码能直接跑改参数就能出图对研究校验矩阵构造、学习LDPC译码的同学很实用。这篇博客把从构造原理到仿真落地的过程完整梳理一遍尽量少讲虚的全是实操。1. 项目需求拆解与整体仿真方案设计1.1 这个仿真到底在做什么LDPC仿真很容易陷入“跑出两条BER曲线就算了”的误区。这个项目的要求拆开看其实只有三个层次构造出确定性的稀疏校验矩阵不能靠随机碰运气。实现BP译码并让迭代次数、码率、码长都成为可配置参数。用蒙特卡洛仿真统计误码率用曲线支撑结论。我一开始也犹豫过是否直接用dv、dc随机构造但随机矩阵容易出现短环尤其是在中短码长下误码率平台很突出。后来从数论里借鉴大衍数的方法用同余方程生成非零元位置矩阵结构确定四环可控跑出来的曲线稳定很多。这个仿真的目标不是追求极致性能而是把“矩阵构造—译码—性能评估”这条链路打通并定量观察三个参数的影响。1.2 为什么选择大衍数构造方法在通信领域校验矩阵构造大体分随机、代数和结构三类。Gallager的随机构造最简单但每次生成的矩阵环分布不同有限几何构造性能好但参数不灵活。大衍数构造本质上是准循环LDPC的一个变种核心思想是用模逆元确定循环移位量。好处是H矩阵由基础矩阵和扩展因子唯一确定不需要存储大矩阵。结构确定利于硬件实现码长码率扩展也方便。使用大衍求一术求逆元参数选取有数论依据不是拍脑袋。关于大衍数很多朋友第一反应是《周易》的“大衍之数五十其用四十有九”。但在工程里我主要用到的是秦九韶《数书九章》中发展出的大衍求一术也就是求解一次同余式组的系统性方法。为构造LDPC矩阵我把扩展因子z作为模数用大衍求一术求解同余方程得到一组与z互素的乘性逆元再按行、列坐标生成循环移位偏移量。这样生成的矩阵具有准循环结构同时可以控制环长。1.3 仿真对比维度与指标设计项目要求同时观察译码迭代次数、码率、码长那就要把对比维度设计清楚否则蒙特卡洛仿真时间会爆炸。我最后确定的设计如下表对比维度取值固定条件译码迭代次数5、10、20、30码率1/2码长648码率1/2、2/3、3/4迭代次数20码长648码长324、648、1296码率1/2迭代次数20误码率统计采用误比特率BER同时记录误帧率FER因为LDPC译码失败往往整帧错。每组信噪比点至少统计100帧错误低误码率时限制最大帧数10万帧保证曲线平滑。这个设计看似简单但跑起来就知道码长1296、码率3/4、10万帧在普通电脑上要几个小时。后面我会讲怎么用Matlab的并行工具加快扫描。2. LDPC编译码与大衍数构造原理2.1 LDPC码稀疏校验矩阵的直观理解LDPC码的核心是用一个稀疏矩阵H描述码字约束关系。所谓稀疏是非零元素远少于零元素。比如648码长、1/2码率的规则LDPC码如果列重是3那么H矩阵非零元素只有648×31944个而总元素是324×648209952个稀疏度只有0.9%。译码时信息通过Tanner图的边在变量节点和校验节点间传递短环会破坏消息独立性所以构造时最怕出现长度为4的环。行列重也很关键。列重太小最小距离低误码平台高列重太大译码门限变差。1/2码率规则码常用列重3、行重6。用大衍数构造时通过基础矩阵的设计可以控制行列重。以码率1/2、码长648为例取扩展因子z54则基础矩阵维度mb6、nb12列重设为3行重自然就是6。2.2 大衍数构造稀疏校验矩阵的算法步骤我实现的构造方法可以拆成四步。第一步确定基础矩阵维度。假设最终H矩阵是m行n列扩展因子为z那么基础矩阵维度是mbm/z行、nbn/z列。为了便于编码通常要求mb远小于nb而且m取z的整数倍。第二步构造基础矩阵的非零图案。对于规则码每个基础矩阵行有dc个1列有dv个1。我用的是确定性循环放置方式而不是随机放置。第j列的非零行索引由下式确定row(k) mod(j - 1 s * (k - 1), mb) 1, k 1,2,...,dv其中s是一个与mb互素的步长。这个步长也由大衍数候选序列中选取从而避免重边问题。第三步用大衍求一术生成循环移位系数。对基础矩阵中每个非零位置(i,j)需要计算一个0到z-1的偏移p(i,j)。具体做法是取一组与z互素的整数序列a然后对每个非零位置求解同余方程a(i) * p(i,j) ≡ (j - c(i)) (mod z)其中的c(i)是行方向偏移由行号i和“其用四十有九”的49取模确定。由于a(i)与z互素这个同余方程有唯一解解出的p(i,j)就用做大衍数偏移量。大衍求一术本质上就是扩展欧几里得算法Matlab里可以直接用gcd函数实现。第四步扩展。把基础矩阵中每个“1”替换成z×z单位阵循环右移p位的矩阵每个“0”替换成全零矩阵得到最终的H矩阵。这个过程在Matlab里可以用稀疏索引构造避免生成完整稠密矩阵。这里还要提一句“问数化定数”。秦九韶在大衍术里把非两两互素的问数转化为两两互素的定数我在构造时对扩展因子z做了类似处理。如果z本身是合数就先把z拆成互素因子的乘积再分别对每个因子求逆最后用中国剩余定理合并偏移量。为了简化首版仿真我直接选择了与所有候选a都互素的z代码更简单后续再扩展这个处理也不难。2.3 校验矩阵质量检查构造完H矩阵后我在代码里强制做了三项检查稀疏度检查非零元占比是否在预期范围。行重列重统计是否满足规则码设定。四环检查遍历任意两列计算公共非零行数量如果大于1就存在四环。四环检查的Matlab实现非常简单function has4 checkCycle4(H) colIdx cell(size(H,2), 1); for j 1:size(H,2) colIdx{j} find(H(:,j)); end has4 false; for j1 1:size(H,2)-1 for j2 j11:size(H,2) if length(intersect(colIdx{j1}, colIdx{j2})) 1 has4 true; return; end end end end如果检查出四环我会返回第一步调整基础矩阵的步长s或换一组a序列。实测下来用大衍数方法比纯随机生成命中无四环矩阵的概率高很多。3. Matlab代码实现与参数配置3.1 主程序结构整个工程的主脚本流程是清空环境、加载参数。调用constructHDayan函数生成H矩阵。求出生成矩阵G完成编码。对每个信噪比点做蒙特卡洛仿真生成随机信息位编码BPSK映射加高斯白噪声LLR-BP译码统计错误。保存BER/FER结果并画图。主脚本的核心参数区块我习惯集中放在文件头部% 仿真参数 EbN0dB 0:0.5:8; z 54; % 扩展因子 mb 6; nb 12; % 基础矩阵维度码率约1/2 R 1 - mb/nb; maxIter 20; % 最大译码迭代次数 numFrames 1e5;这里码率是1/2。如果要做码率对比就换mb和nb组合例如mb4、nb12得到2/3码率mb3、nb12得到3/4码率。码长通过z同步调整z27、54、108分别对应码长324、648、1296。参数之间的换算关系我建议写进脚本注释改起来不容易乱。3.2 构造H矩阵的Matlab实现下面给出constructHDayan函数的核心部分细节做了注释function [H, Hbase] constructHDayan(mb, nb, z, dv) % 第一步确定性基础矩阵非零图案 Hbase zeros(mb, nb); s 1; % 如果与mb不互素换下一个候选 while gcd(s, mb) ~ 1 s s 1; end for col 1:nb for k 1:dv row mod(col - 1 s * (k - 1), mb) 1; Hbase(row, col) 1; end end % 第二步准备与z互素的大衍数候选基数 candidates mod(50 (0:48), z); candidates unique(candidates); candidates(candidates 0) []; aSeq candidates(gcd(candidates, z) 1); if isempty(aSeq) error(扩展因子与候选数不互素请调整z); end % 第三步大衍求一术求模逆元并扩展 H sparse(mb*z, nb*z); cnt 0; for i 1:mb for j 1:nb if Hbase(i,j) 1 cnt cnt 1; a aSeq(mod(cnt-1, length(aSeq)) 1); c mod(49 * i j, z); [g, x, ~] gcd(a, z); if g ~ 1 error(求逆失败); end p mod(x * (j - c), z); % 将单位阵循环右移p位后的非零位置填入H for r 1:z colPos mod(r - 1 p, z) 1; H((i-1)*z r, (j-1)*z colPos) 1; end end end end end这个函数有几个关键点。循环移位方向要统一否则译码端需要额外调整。使用sparse预分配大矩阵能避免内存爆炸。基础矩阵如果有重边需要额外处理我通过步长s与mb互素来避免重边。3.3 译码器实现与迭代次数参数接口译码器用的是对数似然比置信传播也就是SPA。输入是信道软信息LLR输出是硬判决码字。核心迭代过程包括校验节点更新和变量节点更新。为了方便对比迭代次数我把最大迭代次数作为函数参数传入function [xHat, iterUsed] decodeLDPC_BP(H, rxLLR, maxIter) [M, N] size(H); [rowIdx, colIdx] find(H); Lv rxLLR(:).; Lvc zeros(length(rowIdx), 1); % 变量到校验的消息 Lcv zeros(length(rowIdx), 1); % 校验到变量的消息 for iter 1:maxIter % 变量节点更新 % 校验节点更新这里使用tanh规则 % 计算伪后验概率 % 硬判决 iterUsed iter; % 如果校验方程全零提前退出 end endSPA的校验节点更新是最耗时的部分我推荐用tanh规则或者Min-Sum近似。本项目为了精确保真使用tanh但要对数值稳定性做保护比如限制LLR绝对值不超过30。迭代次数接口就是maxIter在扫描脚本里直接赋值。3.4 参数扫描脚本一次性跑三个维度脚本边界要清楚。我专门写了runSweep.m用struct数组存结果results struct(iter, {}, rate, {}, len, {}, EbN0dB, {}, BER, {}); % 迭代次数扫描 for it [5 10 20 30] ber runBER(H, R, z, EbN0dB, maxIterit); results(end1) struct(...); end这里用到了Matlab较新版本的命名参数语法如果版本比较老就改成普通参数传入。每次仿真结束把BER数据保存成mat文件防止电脑崩溃后白跑。4. 仿真结果对比与影响分析4.1 译码迭代次数对误码率的影响先看固定码率1/2、码长648条件下迭代次数从5次增加到30次的BER曲线。结果符合直觉迭代次数从5次提升到10次时性能提升非常明显高信噪比区域误码率能下降一个数量级以上从10次到20次还有可感知的改善再往上到30次曲线几乎没有变化。这说明在中等信噪比下BP译码的收敛速度存在饱和。迭代次数增加带来的是译码时延和功耗线性上升性能收益却递减。实际系统选取迭代次数不能只看误码率曲线还要看吞吐率约束。我习惯在仿真报告里同时给出平均实际迭代次数的曲线它会帮你在“性能天花板”和“复杂度”之间找到平衡点。4.2 码率对系统性能与吞吐的影响码率直接决定冗余度。仿真的三组码率1/2、2/3、3/4中1/2码率曲线最好3/4码率曲线最差。以误码率1e-4为参考2/3码率大概比1/2码率差0.8 dB左右3/4码率差更多。但码率越高有用信息占比越多频谱效率越高。因此实际通信系统选码率本质是误码性能和吞吐性能的折中。自适应调制编码中信道条件差时降到1/2码率信道条件好时切到3/4靠的就是这套仿真数据支撑门限配置。如果把仿真结果换算成频谱效率能更直观地看到码率选择的意义。4.3 码长对误码率性能与复杂度的双向影响码长324、648、1296对比码长越长性能越好这是因为长码的随机化程度好更接近香农限。1296码长在误码率1e-5处比324码长有大约1.2dB的编码增益差异。但长码的代价不只是仿真时间。译码器校验矩阵尺寸变大BP的每一轮迭代复杂度都随之上升存储LLR消息的内存也线性增长。从工程影响范围看LDPC码在高速传输场景中采用较大码长做高速率传输而在物联网短包场景中码长受限就需要用短码结合适当迭代次数来控制时延。这套仿真正好提供了三个关键维度对性能影响的量化结果。5. 常见问题与调优经验5.1 校验矩阵构造中的短环与重边处理用大衍数方法构造时即便有数论基础基础矩阵设计不好仍然会出现四环。我踩过的一个坑是扩展因子z如果与步长s不互素会导致同一列中多个非零行冲突矩阵出现重边。解决办法是检查gcd(a,z)1并且对每个扩展块做非零行列检测。如果真的检测到四环不建议直接重开随机矩阵更快的方法是保持既有非零图案只对冲突的偏移量p(i,j)加一个固定的模增量再做一次求逆。我试过通常调整两三个位置就能消除四环。5.2 编码端G矩阵过慢或内存爆炸用H矩阵直接高斯消元求G矩阵在码长648以后会非常慢而且G矩阵是稠密的1/2码率648码长就要648×324个double内存还好但1296码长就很吃紧。我的处理是用稀疏LU分解解H*x0的零空间或者直接把H构造成近似下三角形式做迭代编码。这样虽然编码有少量开销但整体仿真时间下降非常明显。如果你对编码时延不敏感更简单的方式是用Matlab内置的ldpcEncoderCfg等函数但那样会把矩阵构造过程隔离掉不利于研究校验矩阵本身。所以我这里还是保留了自己求G的流程。5.3 Matlab仿真性能优化与并行扫描蒙特卡洛仿真最忌讳的是在循环里反复找H的非零索引。我的经验是进入信噪比循环之前把所有非零边的行列索引提取出来存成数组BP译码时直接索引。另外LLR初始化用向量化运算不要用for逐比特。信噪比点之间相互独立用parfor并行跑四核机器能省2/3时间。记得先对每个worker预先加载H矩阵否则内存缓存反复传递性能反而下降。还要注意随机数流控制。并行仿真时如果不设置每个worker的随机种子不同信噪比点可能产生相同的随机序列BER曲线不独立。我用的是parfor循环内根据信噪比索引手动设置RandStream确保每个worker的数据独立。5.4 这套仿真工程的扩展方向这个工程不仅限于对比三个参数。把译码器换成Min-Sum或者NMS可以研究低复杂度算法损失把H矩阵构造部分改成5G标准里定义的QC-LDPC可以直接评估标准码的性能加入BICM与高阶调制也能分析编码调制联合方案。大衍数构造核心价值在于参数化与可复现性给大家提供了一个在Matlab里快速验证新想法的底座。最后讲一点我做完这个工程后的体会。开始时我也迷信“迭代次数越多越好”真正把迭代次数从5扫到30才发现中短码长下20次以后几乎是原地踏步。还有码率选择不要只看编码增益要结合频谱效率和使用场景。仿真代码的价值不在于把曲线画得多漂亮而在于每根曲线背后都能回答一个工程问题。这套程序我大概测了一周中间被矩阵奇异、内存爆掉各种折磨但调完之后再去看5G的LDPC参数就顺畅多了。如果你也在做类似对比建议先小码长小迭代数跑通流程再放大参数准能少走弯路。
