简介一份面向无线通信研究者与学生的PDMAPosition Division Multiple Access仿真项目聚焦多用户利用位置信息进行空分多址传输的场景重点研究误块率BLER随信噪比SNR变化的规律可用于验证系统性能并支撑参数优化。压缩包仅含1个主程序文件采用Matlab源码.m形式编写体积仅2KB轻量易用可直接运行或按需修改。已有285人浏览学习适合对PDMA、大规模MIMO或无线物理层仿真感兴趣的中高级学习者。程序内置完整的AWGN信道模型与通信链路流程涵盖信号生成、调制编码、PDMA空间处理、接收检测与错误码块统计等环节通过循环调整噪声功率得到不同SNR下的BLER并绘制关系曲线。借助该仿真读者能直观理解PDMA的干扰分离机制观察高信噪比下误块率下降趋势还可进一步切换天线配置、编码方案等实验条件为实际系统设计提供定量参考。1. 拿到 PDMA 的 MATLAB 仿真包先搞清楚 BLER 曲线在仿什么搜到 PDMA.m.zip 这类压缩包的人基本是奔着两样东西去的一张能在 MATLAB 里复现的 BLER vs SNR 曲线以及能写进论文的仿真说明。PDMA模式分割多址是 5G 非正交多址NOMA家族的一员核心是让多个用户共用同一组时频资源靠模式矩阵拉开用户差异标题里的 BLER、SNR 已经把任务定死——链路级蒙特卡洛仿真不是系统级吞吐仿真。下面按信号模型、MATLAB 最小实现、参数口径、结果验证这条线把代码拆开讲适合刚解压不知道从哪看起的人也适合想把脚本改成自己信道模型的老手。全文对着最常见的 4×6 模式矩阵展开。2. 仿真骨架先行PDMA 模式矩阵、信号模型与 BLER 统计口径写代码之前先把三个概念钉死Pat 怎么读、接收信号长什么样、块错误怎么数。这三件事任何一个理解偏了后面跑出来的曲线都不可信。2.1 为什么 PDMA 的模式矩阵是 4×6 而不是方阵Pat 的行为资源组RE列为用户元素为 1 表示该用户占用该 RE。4 行 6 列的意思就是 6 个用户挤在 4 个 RE 上过载率 J/K 150%。一旦矩阵变成方阵系统退化为正交多址PDMA 相对 OFDMA 的增益也就不存在。压缩包里最常见的规则矩阵长这样Pat [1 1 1 0 0 0; 1 0 0 1 1 0; 0 1 0 1 0 1; 0 0 1 0 1 1]; % 6用户, 4RE; 每列列重2, 每行行重3这个矩阵的规则性体现在两处每个用户出现在两个 RE 上列重 d_j2提供二重分集每个 RE 上叠了三个用户行重 d_k3直接决定 MPA 的枚举复杂度。列重决定曲线高 SNR 段的斜率行重决定接收端计算量这两个数就是 PDMA 设计里最核心的权衡。还有一类论文会给 Pat 换成复加权值比如每列乘一个旋转相位那是增强版方案手头这份 zip 大概率是 0/1 矩阵下面的分析全部按 0/1 来。2.2 一个时隙的收发模型y (Pat .* H) s n记某个时隙 6 个用户的调制符号为 sJ×1信道矩阵 HK×J的元素 h_kj 是 RE k 到用户 j 的等效复信道。因为 Pat 把 s_j 复制到它占用的所有 RE 上接收向量能写成一条很紧凑的式子y (Pat .* H) * s n(Pat .* H) 是逐元素相乘把没连接的路径直接清零后面乘 s 完成叠加n 是复高斯白噪声每个元素方差 σ²。这个式子就是 PDMA 链路仿真的公共发送端模型不同代码的差别只在 Pat 的取法和 H 的分布。H 在块内固定、块间独立变化就是准静态 Rayleigh 信道直接把 H 设成全 1 矩阵模型退化成 AWGN 基线用来验证程序非常方便。这里容易和另一种写法混淆有人直接用 H * s 再乘 Pat 的转置把 Pat 放在信道后面当第二层预编码。语义上不对——PDMA 的模式映射发生在发送端符号先按列复制再进入各自的信道所以掩码必须作用在 H 上。代码里写成 (Pat .* h) * symPat 和 h 同尺寸逐元素乘正好表达只有占用路径参与叠加。2.3 BLER 的统计口径块怎么定义曲线就怎么解读块错误率的前提是块有明确边界。常见口径有两种按用户独立统计每个用户各成一个块出错单独计数或按系统统计任意用户出错就算整块错。压缩包里如果只有一个 bler 变量多半是系统口径混着平均读代码时要先确认清楚否则画出来的曲线对不上论文。仿真参数常见取值对曲线的影响每块时隙数 L2496块越短 BLER 越高抖动越大每 SNR 点块数10005000决定平滑度太少毛刺多MPA 迭代次数410太少出现平台太多只耗时间信道类型AWGN / 准静态 Rayleigh决定斜率与分集阶数块数不是拍脑袋定的。BLER 估计值的标准差大约是 sqrt(BLER·(1−BLER)/N)要让 10⁻² 附近的点不抖每个 SNR 点至少数到 100200 个错误块。反过来算所需的块数targetErr 100; % 希望数到的错误块数 estBler 1e-2; % 该点预估 BLER needN ceil(targetErr / estBler); % 需要 10000 块这里能看出仿真为什么越到高 SNR 越慢BLER 每降一个量级所需块数就升一个量级。第 4 章写退出条件时会基于这个公式。3. 用 MATLAB 自建 PDMA 的 BLER vs SNR 最小链路3.1 主脚本的参数区这套仿真不需要任何工具箱纯基础 MATLAB 就能跑这也是链路级仿真适合当毕业论文代码的原因。参数区把系统规模、调制、信道和统计量全部集中在一起方便整体替换。% pdma_bler.m —— PDMA 链路级 BLER vs SNR 仿真 clear; clc; rng(20240901); M 4; % QPSK, 每符号2比特 J 6; K 4; % 6个用户共享4个RE L 24; % 每块时隙数, 每用户每块48比特 nIter 6; % MPA 迭代次数 nBlock 3000; % 每个SNR点的块数 snrEbN0 0:2:16; % 横轴 Eb/N0, dB const [11i, -11i, 1-1i, -1-1i] / sqrt(2); % QPSK星座 Pat [1 1 1 0 0 0; 1 0 0 1 1 0; 0 1 0 1 0 1; 0 0 1 0 1 1];const 的四个星座点按格雷映射排列实部决定第 1 bit虚部决定第 2 bit判决时直接取符号比较不用查表。rng 固定种子是关键习惯换一台机器跑同样代码才能得到同一条曲线。提示先固定随机种子再谈结果复现这是论文和评审里最常被忽略的细节。把参数区、发收主循环、MPA 检测器拆成三个文件是更工程化的做法主循环里只留体系结构检测器单独写成函数方便替换成 SIC 或 MMSE 做对照。小规模仿真时放一个文件里反而好调试变量都在工作区可以直接看。3.2 发送端、信道与噪声功率QPSK 符号能量 Es1。若横轴是 Eb/N0(dB)每符号带 2 bit则 Es/N0 的线性值是 Eb/N0 线性值的 2 倍噪声方差 σ² 0.5·10^(−EbN0/10)。这个 0.5 是后面最容易错的地方。for si 1:numel(snrEbN0) sigma2 0.5 * 10^(-snrEbN0(si)/10); errBlk 0; for blk 1:nBlock bits randi([0 1], J, L*2); % 每用户48比特 sym (2*bits(:,1:2:end)-1 1i*(2*bits(:,2:2:end)-1))/sqrt(2); h (randn(K,J) 1i*randn(K,J)) / sqrt(2); % 准静态 y (Pat .* h) * sym; % K x L 用户叠加 y y sqrt(sigma2/2)*(randn(K,L) 1i*randn(K,L)); % 检测与判决见3.3节, 得到 bitsHat(J, L*2) % 若 any(bitsHat(:)~bits(:)), errBlk errBlk 1; end bler(si) errBlk / blk; endh 除以 √2 是让每个元素的实部虚部方差各为 1/2保证 |h|² 的期望为 1这样 Eb/N0 的定义在全平均意义下成立。噪声用 sqrt(sigma2/2) 拆到实部虚部生成的是方差 σ² 的复高斯。3.3 MPA 检测器消息在因子图上往返PDMA 接收端几乎都用消息传递算法MPA。每个 RE 是一个函数节点每个用户是一个变量节点Pat 里的 1 就是边。算法只有两个动作函数节点用接收值和高斯似然汇总相连用户的先验消息变量节点把相连函数节点的消息相乘再归一化。function Mf2v pdma_mpa(y1, h1, Pat, const, sigma2, nIter) % y1: Kx1 单时隙接收向量; h1: KxJ 信道 [K, J] size(Pat); P numel(const); Mv2f ones(J, K, P) / P; % 消息初值: 均匀分布 for it 1:nIter Mf2v ones(J, K, P); for k 1:K u find(Pat(k,:)); % 该RE上的用户集合 if isempty(u), continue; end d numel(u); g cell(1,d); [g{:}] ndgrid(1:P); idxG reshape(cat(d, g{:}), [], d); % 全部符号组合 symG const(idxG); residual y1(k) - sum(symG .* h1(k, u), 2); like exp(-abs(residual).^2 / (2*sigma2)); % 高斯似然 for t 1:d % 对RE上每个用户 j u(t); pm ones(size(idxG,1),1); for q 1:d if q ~ t % 其余用户的先验 col reshape(Mv2f(u(q), k, :), P, 1); pm pm .* col(idxG(:,q)); end end for p 1:P mask idxG(:,t) p; Mf2v(j, k, p) sum(like(mask) .* pm(mask)); end Mf2v(j, k, :) Mf2v(j,k,:) / sum(Mf2v(j,k,:)); % 归一化 end end for j 1:J % 变量节点更新 rsc find(Pat(:,j)); for k rsc pr ones(P,1); for k2 rsc if k2 ~ k pr pr .* reshape(Mf2v(j,k2,:), P, 1); end end Mv2f(j,k,:) pr / sum(pr); end end end end复杂度由行重决定每个 RE 每轮要枚举 M^d 个符号组合d3 时是 64 个列重越大分集越好但行重跟着涨d4 就跳到 256 个。判决时把每个用户相连 RE 的 Mf2v 相乘取最大概率符号按实虚部映射回比特。一帧 24 个时隙、6 轮迭代在普通笔记本上跑完整个 SNR 区间也不到一分钟。大规模跑批时在 SNR 循环外加 tic/toc或打印进度否则高 SNR 点跑多久心里没数。3.4 曲线该长什么样AWGN 下 BLER 曲线应有明显瀑布区准静态 Rayleigh 下斜率变缓高频段斜率约等于列重 2也就是二重分集。如果曲线在高 SNR 变成水平线优先查迭代次数和噪声方差而不是怀疑 Pat 配错。观测现象最可能原因排查顺序高 SNR 平台nIter 太小或消息未归一化先加迭代数再看发散整条曲线左右平移σ² 少了 0.5 倍核对 Eb/N0 口径曲线比预期陡H 用了实高斯确认复数信道4. 三个必调参数与两个经典坑迭代次数、块数和 SNR 口径4.1 MPA 迭代次数消息不动了就停迭代次数少了必然出平台多了只是白耗时间。判断标准不是拍脑袋而是看消息变化量连续两轮每个用户后验分布的最大变化小于 1e-3就可以认为收敛。4×6 矩阵在 QPSK 下 nIter6 通常够用保守一点取 8。真正容易翻车的是消息没归一化小概率连乘会快速下溢大约 10 轮之后就出现 NaN整条曲线在高 SNR 段断裂这就是常说的仿真发散。上面代码里每个消息算完都归一化就是为了规避这个问题。若想在更大矩阵上跑建议直接转到 log 域用加法代替乘法。如果熵监控发现长期震荡给更新加阻尼系数消息新旧值按 0.7/0.3 混合通常能压住这是工程上处理因子图环的通用手段。4.2 每 SNR 点的块数100 个错误块是底线固定块数跑完再算 BLER最大的问题是高 SNR 点可能一个错误都没数到曲线尾部全是噪声。工程上我会按错误块数控制循环至少数到 100 个错误块才停再用实际执行的块数做分母。maxBlk 50000; % 上限, 防止低SNR点跑太久 minErr 100; % 至少数到100个错误块 for blk 1:maxBlk % ... 一帧仿真 ... if any(bitsHat(:) ~ bits(:)) errBlk errBlk 1; end if errBlk minErr blk 1000 break; end end bler(si) errBlk / blk; % 分母是实际块数, 不是maxBlk按错误块数退出的好处是自适应低 SNR 点不会过早截断高 SNR 点不会白跑。配合第 2 章的公式一个点最多也就几万块运行时间完全可控。4.3 Eb/N0 与 Es/N0 的 3 dB 之谜这是 PDMA 仿真返工率最高的原因。横轴写 Eb/N0 时σ² 0.5·10^(−EbN0/10)横轴写 Es/N0 时σ² 10^(−EsN0/10)。两者差一个 0.5对应 3 dB 平移。很多人把 Es/N0 的噪声公式配上 Eb/N0 的横轴整条曲线右移 3 dB还以为是算法实现错了。横轴名称QPSK 下 σ² 表达式 (Es1)常见错误Eb/N0 (dB)0.5·10^(−ebdB/10)漏乘 0.5Es/N0 (dB)10^(−esdB/10)误乘 0.5注意σ² 是复噪声的方差不是每维方差。生成噪声时用 sqrt(sigma2/2) 乘在实部和虚部上这样才是总方差 σ²。另外还有一个归一化口径列重 d_j1 时同一符号在多个 RE 上重复发送。本代码定义每用户符号能量 Es1重复发送不额外摊功率若改成每个 RE 上每符号能量为 1曲线会整体平移两种口径都自洽但论文里必须写明用了哪一种。4.4 列重与过载率曲线斜率的来源列重直接等于分集阶数。把 Pat 换成 4×4 单位阵PDMA 退化成正交 OFDMA每条曲线斜率明显变缓这是验证PDMA 有没有带来增益最直观的对照实验。要提高分集就加列重代价是行重上升、MPA 组合数指数增长d_k4 时每 RE 每轮要枚举 256 个组合跑一帧的时间直接翻倍。过载率选 1.5 还是 2.0 取决于场景过载越高谱效率越好但 MPA 组合数随行重 d_k 指数膨胀仿真时间会先于性能到达极限。除了改 0/1 位置还可以给 Pat 引入复权值例如给第 i 个用户乘相位旋转 e^{jφ_i}在不增加行重的前提下改善分集这就是论文里常说加权模式矩阵的来源。5. 让 PDMA 仿真结果可信三条验证手段5.1 用单用户退化特例核对整条链路把 Pat 改成全 1 列J1、K2 时 PDMA 退化为单用户两路重复发送接收端等价于最大比合并。QPSK 在 Rayleigh 下 MRC 的 BLER 有解析表达式仿真曲线和解析曲线对上说明发送端、信道、噪声、判决链路全部正确问题只能出在 MPA 内部。这步十分钟能做完比反复检查索引高效得多。Pat1 ones(2,1); % 1用户, 2个RE重复 % 等效接收SNR (|h1|^2|h2|^2)*Es/N0, 服从2阶Chi-square分布 % 用解析BLER与仿真曲线对照, 偏差超过0.5dB就要查代码5.2 固定信道做 MAP 对照为确认 MPA 收敛性每个时隙把所有 4^64096 个用户符号组合枚举一次做最大后验判决再与 MPA 判决结果比对。两者在小 L 下应完全一致若某个 SNR 段出现分歧基本是迭代轮数不足或消息归一化位置写错。这个对照只跑几十帧即可但能一次性排除因子图实现的索引错误。5.3 用后验熵判断要不要加迭代BLER 曲线对收敛不敏感后验熵却非常灵敏。每个用户后验分布的熵降到接近 0说明判决已确信再加迭代只是空转熵长期不降则说明因子图有环导致消息震荡要靠 Pat 设计或消息阻尼处理。循环里顺手打印H -sum(post .* log(post eps), 2); % 逐用户后验熵 if max(H) 1e-3, break; end % 已收敛这套退化特例 MAP 对照 熵监控三板斧下来BLER vs SNR 曲线基本不会再被挑出硬伤。最后一个技巧AWGN 场景下把 MPA 里的 exp 计算展开成查表能省约三成运行时间SNR 点一多差别就很明显。本文还有配套的精品资源点击获取 PDMA链路仿真：MATLAB实现BLER曲线与MPA检测全解析
简介一份面向无线通信研究者与学生的PDMAPosition Division Multiple Access仿真项目聚焦多用户利用位置信息进行空分多址传输的场景重点研究误块率BLER随信噪比SNR变化的规律可用于验证系统性能并支撑参数优化。压缩包仅含1个主程序文件采用Matlab源码.m形式编写体积仅2KB轻量易用可直接运行或按需修改。已有285人浏览学习适合对PDMA、大规模MIMO或无线物理层仿真感兴趣的中高级学习者。程序内置完整的AWGN信道模型与通信链路流程涵盖信号生成、调制编码、PDMA空间处理、接收检测与错误码块统计等环节通过循环调整噪声功率得到不同SNR下的BLER并绘制关系曲线。借助该仿真读者能直观理解PDMA的干扰分离机制观察高信噪比下误块率下降趋势还可进一步切换天线配置、编码方案等实验条件为实际系统设计提供定量参考。1. 拿到 PDMA 的 MATLAB 仿真包先搞清楚 BLER 曲线在仿什么搜到 PDMA.m.zip 这类压缩包的人基本是奔着两样东西去的一张能在 MATLAB 里复现的 BLER vs SNR 曲线以及能写进论文的仿真说明。PDMA模式分割多址是 5G 非正交多址NOMA家族的一员核心是让多个用户共用同一组时频资源靠模式矩阵拉开用户差异标题里的 BLER、SNR 已经把任务定死——链路级蒙特卡洛仿真不是系统级吞吐仿真。下面按信号模型、MATLAB 最小实现、参数口径、结果验证这条线把代码拆开讲适合刚解压不知道从哪看起的人也适合想把脚本改成自己信道模型的老手。全文对着最常见的 4×6 模式矩阵展开。2. 仿真骨架先行PDMA 模式矩阵、信号模型与 BLER 统计口径写代码之前先把三个概念钉死Pat 怎么读、接收信号长什么样、块错误怎么数。这三件事任何一个理解偏了后面跑出来的曲线都不可信。2.1 为什么 PDMA 的模式矩阵是 4×6 而不是方阵Pat 的行为资源组RE列为用户元素为 1 表示该用户占用该 RE。4 行 6 列的意思就是 6 个用户挤在 4 个 RE 上过载率 J/K 150%。一旦矩阵变成方阵系统退化为正交多址PDMA 相对 OFDMA 的增益也就不存在。压缩包里最常见的规则矩阵长这样Pat [1 1 1 0 0 0; 1 0 0 1 1 0; 0 1 0 1 0 1; 0 0 1 0 1 1]; % 6用户, 4RE; 每列列重2, 每行行重3这个矩阵的规则性体现在两处每个用户出现在两个 RE 上列重 d_j2提供二重分集每个 RE 上叠了三个用户行重 d_k3直接决定 MPA 的枚举复杂度。列重决定曲线高 SNR 段的斜率行重决定接收端计算量这两个数就是 PDMA 设计里最核心的权衡。还有一类论文会给 Pat 换成复加权值比如每列乘一个旋转相位那是增强版方案手头这份 zip 大概率是 0/1 矩阵下面的分析全部按 0/1 来。2.2 一个时隙的收发模型y (Pat .* H) s n记某个时隙 6 个用户的调制符号为 sJ×1信道矩阵 HK×J的元素 h_kj 是 RE k 到用户 j 的等效复信道。因为 Pat 把 s_j 复制到它占用的所有 RE 上接收向量能写成一条很紧凑的式子y (Pat .* H) * s n(Pat .* H) 是逐元素相乘把没连接的路径直接清零后面乘 s 完成叠加n 是复高斯白噪声每个元素方差 σ²。这个式子就是 PDMA 链路仿真的公共发送端模型不同代码的差别只在 Pat 的取法和 H 的分布。H 在块内固定、块间独立变化就是准静态 Rayleigh 信道直接把 H 设成全 1 矩阵模型退化成 AWGN 基线用来验证程序非常方便。这里容易和另一种写法混淆有人直接用 H * s 再乘 Pat 的转置把 Pat 放在信道后面当第二层预编码。语义上不对——PDMA 的模式映射发生在发送端符号先按列复制再进入各自的信道所以掩码必须作用在 H 上。代码里写成 (Pat .* h) * symPat 和 h 同尺寸逐元素乘正好表达只有占用路径参与叠加。2.3 BLER 的统计口径块怎么定义曲线就怎么解读块错误率的前提是块有明确边界。常见口径有两种按用户独立统计每个用户各成一个块出错单独计数或按系统统计任意用户出错就算整块错。压缩包里如果只有一个 bler 变量多半是系统口径混着平均读代码时要先确认清楚否则画出来的曲线对不上论文。仿真参数常见取值对曲线的影响每块时隙数 L2496块越短 BLER 越高抖动越大每 SNR 点块数10005000决定平滑度太少毛刺多MPA 迭代次数410太少出现平台太多只耗时间信道类型AWGN / 准静态 Rayleigh决定斜率与分集阶数块数不是拍脑袋定的。BLER 估计值的标准差大约是 sqrt(BLER·(1−BLER)/N)要让 10⁻² 附近的点不抖每个 SNR 点至少数到 100200 个错误块。反过来算所需的块数targetErr 100; % 希望数到的错误块数 estBler 1e-2; % 该点预估 BLER needN ceil(targetErr / estBler); % 需要 10000 块这里能看出仿真为什么越到高 SNR 越慢BLER 每降一个量级所需块数就升一个量级。第 4 章写退出条件时会基于这个公式。3. 用 MATLAB 自建 PDMA 的 BLER vs SNR 最小链路3.1 主脚本的参数区这套仿真不需要任何工具箱纯基础 MATLAB 就能跑这也是链路级仿真适合当毕业论文代码的原因。参数区把系统规模、调制、信道和统计量全部集中在一起方便整体替换。% pdma_bler.m —— PDMA 链路级 BLER vs SNR 仿真 clear; clc; rng(20240901); M 4; % QPSK, 每符号2比特 J 6; K 4; % 6个用户共享4个RE L 24; % 每块时隙数, 每用户每块48比特 nIter 6; % MPA 迭代次数 nBlock 3000; % 每个SNR点的块数 snrEbN0 0:2:16; % 横轴 Eb/N0, dB const [11i, -11i, 1-1i, -1-1i] / sqrt(2); % QPSK星座 Pat [1 1 1 0 0 0; 1 0 0 1 1 0; 0 1 0 1 0 1; 0 0 1 0 1 1];const 的四个星座点按格雷映射排列实部决定第 1 bit虚部决定第 2 bit判决时直接取符号比较不用查表。rng 固定种子是关键习惯换一台机器跑同样代码才能得到同一条曲线。提示先固定随机种子再谈结果复现这是论文和评审里最常被忽略的细节。把参数区、发收主循环、MPA 检测器拆成三个文件是更工程化的做法主循环里只留体系结构检测器单独写成函数方便替换成 SIC 或 MMSE 做对照。小规模仿真时放一个文件里反而好调试变量都在工作区可以直接看。3.2 发送端、信道与噪声功率QPSK 符号能量 Es1。若横轴是 Eb/N0(dB)每符号带 2 bit则 Es/N0 的线性值是 Eb/N0 线性值的 2 倍噪声方差 σ² 0.5·10^(−EbN0/10)。这个 0.5 是后面最容易错的地方。for si 1:numel(snrEbN0) sigma2 0.5 * 10^(-snrEbN0(si)/10); errBlk 0; for blk 1:nBlock bits randi([0 1], J, L*2); % 每用户48比特 sym (2*bits(:,1:2:end)-1 1i*(2*bits(:,2:2:end)-1))/sqrt(2); h (randn(K,J) 1i*randn(K,J)) / sqrt(2); % 准静态 y (Pat .* h) * sym; % K x L 用户叠加 y y sqrt(sigma2/2)*(randn(K,L) 1i*randn(K,L)); % 检测与判决见3.3节, 得到 bitsHat(J, L*2) % 若 any(bitsHat(:)~bits(:)), errBlk errBlk 1; end bler(si) errBlk / blk; endh 除以 √2 是让每个元素的实部虚部方差各为 1/2保证 |h|² 的期望为 1这样 Eb/N0 的定义在全平均意义下成立。噪声用 sqrt(sigma2/2) 拆到实部虚部生成的是方差 σ² 的复高斯。3.3 MPA 检测器消息在因子图上往返PDMA 接收端几乎都用消息传递算法MPA。每个 RE 是一个函数节点每个用户是一个变量节点Pat 里的 1 就是边。算法只有两个动作函数节点用接收值和高斯似然汇总相连用户的先验消息变量节点把相连函数节点的消息相乘再归一化。function Mf2v pdma_mpa(y1, h1, Pat, const, sigma2, nIter) % y1: Kx1 单时隙接收向量; h1: KxJ 信道 [K, J] size(Pat); P numel(const); Mv2f ones(J, K, P) / P; % 消息初值: 均匀分布 for it 1:nIter Mf2v ones(J, K, P); for k 1:K u find(Pat(k,:)); % 该RE上的用户集合 if isempty(u), continue; end d numel(u); g cell(1,d); [g{:}] ndgrid(1:P); idxG reshape(cat(d, g{:}), [], d); % 全部符号组合 symG const(idxG); residual y1(k) - sum(symG .* h1(k, u), 2); like exp(-abs(residual).^2 / (2*sigma2)); % 高斯似然 for t 1:d % 对RE上每个用户 j u(t); pm ones(size(idxG,1),1); for q 1:d if q ~ t % 其余用户的先验 col reshape(Mv2f(u(q), k, :), P, 1); pm pm .* col(idxG(:,q)); end end for p 1:P mask idxG(:,t) p; Mf2v(j, k, p) sum(like(mask) .* pm(mask)); end Mf2v(j, k, :) Mf2v(j,k,:) / sum(Mf2v(j,k,:)); % 归一化 end end for j 1:J % 变量节点更新 rsc find(Pat(:,j)); for k rsc pr ones(P,1); for k2 rsc if k2 ~ k pr pr .* reshape(Mf2v(j,k2,:), P, 1); end end Mv2f(j,k,:) pr / sum(pr); end end end end复杂度由行重决定每个 RE 每轮要枚举 M^d 个符号组合d3 时是 64 个列重越大分集越好但行重跟着涨d4 就跳到 256 个。判决时把每个用户相连 RE 的 Mf2v 相乘取最大概率符号按实虚部映射回比特。一帧 24 个时隙、6 轮迭代在普通笔记本上跑完整个 SNR 区间也不到一分钟。大规模跑批时在 SNR 循环外加 tic/toc或打印进度否则高 SNR 点跑多久心里没数。3.4 曲线该长什么样AWGN 下 BLER 曲线应有明显瀑布区准静态 Rayleigh 下斜率变缓高频段斜率约等于列重 2也就是二重分集。如果曲线在高 SNR 变成水平线优先查迭代次数和噪声方差而不是怀疑 Pat 配错。观测现象最可能原因排查顺序高 SNR 平台nIter 太小或消息未归一化先加迭代数再看发散整条曲线左右平移σ² 少了 0.5 倍核对 Eb/N0 口径曲线比预期陡H 用了实高斯确认复数信道4. 三个必调参数与两个经典坑迭代次数、块数和 SNR 口径4.1 MPA 迭代次数消息不动了就停迭代次数少了必然出平台多了只是白耗时间。判断标准不是拍脑袋而是看消息变化量连续两轮每个用户后验分布的最大变化小于 1e-3就可以认为收敛。4×6 矩阵在 QPSK 下 nIter6 通常够用保守一点取 8。真正容易翻车的是消息没归一化小概率连乘会快速下溢大约 10 轮之后就出现 NaN整条曲线在高 SNR 段断裂这就是常说的仿真发散。上面代码里每个消息算完都归一化就是为了规避这个问题。若想在更大矩阵上跑建议直接转到 log 域用加法代替乘法。如果熵监控发现长期震荡给更新加阻尼系数消息新旧值按 0.7/0.3 混合通常能压住这是工程上处理因子图环的通用手段。4.2 每 SNR 点的块数100 个错误块是底线固定块数跑完再算 BLER最大的问题是高 SNR 点可能一个错误都没数到曲线尾部全是噪声。工程上我会按错误块数控制循环至少数到 100 个错误块才停再用实际执行的块数做分母。maxBlk 50000; % 上限, 防止低SNR点跑太久 minErr 100; % 至少数到100个错误块 for blk 1:maxBlk % ... 一帧仿真 ... if any(bitsHat(:) ~ bits(:)) errBlk errBlk 1; end if errBlk minErr blk 1000 break; end end bler(si) errBlk / blk; % 分母是实际块数, 不是maxBlk按错误块数退出的好处是自适应低 SNR 点不会过早截断高 SNR 点不会白跑。配合第 2 章的公式一个点最多也就几万块运行时间完全可控。4.3 Eb/N0 与 Es/N0 的 3 dB 之谜这是 PDMA 仿真返工率最高的原因。横轴写 Eb/N0 时σ² 0.5·10^(−EbN0/10)横轴写 Es/N0 时σ² 10^(−EsN0/10)。两者差一个 0.5对应 3 dB 平移。很多人把 Es/N0 的噪声公式配上 Eb/N0 的横轴整条曲线右移 3 dB还以为是算法实现错了。横轴名称QPSK 下 σ² 表达式 (Es1)常见错误Eb/N0 (dB)0.5·10^(−ebdB/10)漏乘 0.5Es/N0 (dB)10^(−esdB/10)误乘 0.5注意σ² 是复噪声的方差不是每维方差。生成噪声时用 sqrt(sigma2/2) 乘在实部和虚部上这样才是总方差 σ²。另外还有一个归一化口径列重 d_j1 时同一符号在多个 RE 上重复发送。本代码定义每用户符号能量 Es1重复发送不额外摊功率若改成每个 RE 上每符号能量为 1曲线会整体平移两种口径都自洽但论文里必须写明用了哪一种。4.4 列重与过载率曲线斜率的来源列重直接等于分集阶数。把 Pat 换成 4×4 单位阵PDMA 退化成正交 OFDMA每条曲线斜率明显变缓这是验证PDMA 有没有带来增益最直观的对照实验。要提高分集就加列重代价是行重上升、MPA 组合数指数增长d_k4 时每 RE 每轮要枚举 256 个组合跑一帧的时间直接翻倍。过载率选 1.5 还是 2.0 取决于场景过载越高谱效率越好但 MPA 组合数随行重 d_k 指数膨胀仿真时间会先于性能到达极限。除了改 0/1 位置还可以给 Pat 引入复权值例如给第 i 个用户乘相位旋转 e^{jφ_i}在不增加行重的前提下改善分集这就是论文里常说加权模式矩阵的来源。5. 让 PDMA 仿真结果可信三条验证手段5.1 用单用户退化特例核对整条链路把 Pat 改成全 1 列J1、K2 时 PDMA 退化为单用户两路重复发送接收端等价于最大比合并。QPSK 在 Rayleigh 下 MRC 的 BLER 有解析表达式仿真曲线和解析曲线对上说明发送端、信道、噪声、判决链路全部正确问题只能出在 MPA 内部。这步十分钟能做完比反复检查索引高效得多。Pat1 ones(2,1); % 1用户, 2个RE重复 % 等效接收SNR (|h1|^2|h2|^2)*Es/N0, 服从2阶Chi-square分布 % 用解析BLER与仿真曲线对照, 偏差超过0.5dB就要查代码5.2 固定信道做 MAP 对照为确认 MPA 收敛性每个时隙把所有 4^64096 个用户符号组合枚举一次做最大后验判决再与 MPA 判决结果比对。两者在小 L 下应完全一致若某个 SNR 段出现分歧基本是迭代轮数不足或消息归一化位置写错。这个对照只跑几十帧即可但能一次性排除因子图实现的索引错误。5.3 用后验熵判断要不要加迭代BLER 曲线对收敛不敏感后验熵却非常灵敏。每个用户后验分布的熵降到接近 0说明判决已确信再加迭代只是空转熵长期不降则说明因子图有环导致消息震荡要靠 Pat 设计或消息阻尼处理。循环里顺手打印H -sum(post .* log(post eps), 2); % 逐用户后验熵 if max(H) 1e-3, break; end % 已收敛这套退化特例 MAP 对照 熵监控三板斧下来BLER vs SNR 曲线基本不会再被挑出硬伤。最后一个技巧AWGN 场景下把 MPA 里的 exp 计算展开成查表能省约三成运行时间SNR 点一多差别就很明显。本文还有配套的精品资源点击获取