
简介面向MIMO无线通信系统物理层研究的MATLAB仿真源码包聚焦SVD、BD、ZF、MF、SLNR与MMSE六类线性预编码方案。资源针对多用户MIMO场景提供误码率、和速率随信噪比及天线数变化的完整对比实验适合通信工程专业学生、算法工程师及科研人员用于算法验证与性能分析。压缩包共34个文件其中14个m脚本覆盖预编码实现、QPSK映射、注水算法及接收机处理等核心模块20个fig结果为可直接复现的仿真曲线图整体仅207KB轻量易用。已有228人学习下载。通过阅读源码与运行程序可快速理解各预编码原理差异掌握多用户干扰抑制与速率折中分析方法也可基于现有代码扩展新算法或修改仿真参数是开展教学实践和课题预研的实用工具。1. 线性预编码性能对比为什么五种算法值得放在同一份仿真里比做多用户MIMO下行链路仿真时最耗时间的不是推导公式而是把不同预编码算法跑在同一套信道条件下还要保证曲线能对上别人的结果。下载这份svd_BD_ZF_MF_SLNR_MMSE线性预编码性能对比_源码你拿到的不只是五个算法的独立实现还有一个完整的对比框架两个主脚本分别负责 sum-rate 和误码率仿真ZF、MMSE、MF、SLNR、BD 共享同一套信道模型和调制方式。资源里还附带了多张.fig结果图能直接对照验证自己的仿真曲线是否跑偏。它适合正在做多用户 MIMO 方向课程设计或毕业论文的人尤其是想快速获得“天线数变化对性能影响”这组结果的研究者。2. 先看懂资源结构主脚本、子函数与结果图之间的调用关系拿到一个 MATLAB 仿真工程我一般不会急着跑main.m而是先把文件列表过一遍搞清楚谁是入口、谁是工具函数、谁负责画图。这份源码的文件命名比较直白但没有文档直接双击main1.m可能会因为缺失路径或函数名敲错而中断。2.1 线性预编码的适用边界为什么不做 THP 和 Vector Perturbation这份资源只对比线性预编码原因很实际非线性方法比如 THPTomlinson-Harashima Precoding或脏纸编码虽然在高信噪比下有更好的 sum-rate但实现涉及逐用户排序、取模运算和量化复杂度会随着用户数上升得非常快。在 4 个用户、每个用户 2 根接收天线的配置下线性预编码已经能看出算法之间的本质差异。从SNR5dB时4个2天线用户MIMO系统合速率随Nt变化.fig这个文件名可以看出仿真边界是固定的用户数 K4每个用户接收天线数 Nr2信噪比固定 5dB基站发射天线数 Nt 是变量。这是多用户 MIMO 里最经典的配置因为当 Nt 从 4 增加到 10 时自由度逐渐充足ZF 和 BD 的增益会被明显拉开。信道模型方面源码里大概率用的是平衰落瑞利信道也就是每个元素都服从标准复高斯分布。这个选择合理性很强如果是频率选择性信道就得引入 OFDM 子载波维度代码复杂度翻倍反而会模糊五种预编码之间差异的主线。2.2 关键文件清单与仿真流程的主线整理一下资源中的 MATLAB 文件可以按职责分成三类。类别文件名作用主脚本main1.m跑 sum-rate 随 Nt 变化生成sum-rate*.fig主脚本main2.m跑 BER 随 SNR 变化生成BER*.fig算法实现ZF.m、MMSE.m、MF.m、SLNR.m、BD.m各自返回预编码矩阵 W辅助工具svdprecoding.m、waterfilling.m、QPSK_mapper.m、receiver.m提供 SVD 分解、注水功率分配、调制映射、接收检测从函数命名看svdprecoding.m可能是把 SVD 分解封装了一次waterfilling.m负责给 BD/SLNR 后的等效信道分配功率。这两个文件很容易被忽略但它们直接决定 sum-rate 曲线是否合理。我见过很多人在复现时只改ZF.m却忘了注水功率分配会改变不同算法的功率基准结果 ZF 和 MF 的曲线相对位置发生平移结论整个反掉。调用关系可以用下面这段伪代码表示% main1.m 的核心调用链 H (randn(K*Nr, Nt) 1j*randn(K*Nr, Nt)) / sqrt(2); % 生成信道 W_ZF ZF(H, P); % 调用 ZF 预编码 W_MMSE MMSE(H, P); % 调用 MMSE 预编码 W_MF MF(H, P); % 调用 MF 预编码 W_SLNR SLNR(H, P, sigma2); % 注意需要噪声方差 W_BD BD(H, K, Nr, P); % BD 需要知道用户分组信息 rate calculate_sum_rate(H, W_ZF, sigma2); % 计算合速率H的维度是K*Nr行、Nt列因为所有用户接收天线是拼接在一起的。这也是最常见的翻车点如果某个函数内部把H转置了导致维度对不上MATLAB 会报错但更危险的是某些情况下隐式广播让计算能跑通结果却完全错误。后面第 5 章会专门讲这类问题。3. 五种预编码实现从公式到 MATLAB 核心代码这一章是资源的核心我会把每种算法在源码里的关键代码结构拆开讲。注意源码是.m文件我这里还原的是最常见的实现方式具体的归一化因子可能与你看到的版本有细微差异但算法主干是一致的。3.1 ZF 与 MF一对“零迫”与“匹配”的对称设计零迫预编码的目标是将干扰完全置零做法是把等效信道变为单位阵。对于下行链路接收信号可以写成y H W x n其中W就是预编码矩阵。ZF 的闭式解是W_ZF H^H (H H^H)^{-1}这里H^H是H的共轭转置H H^H的维度是K*Nr × K*Nr通常可逆。匹配滤波 MF 则简单得多它直接取W_MF H^H不考虑任何求逆所以低信噪比下噪声抑制最强高信噪比下会被用户间干扰限制住。MATLAB 实现里核心是这一行% ZF.m 核心代码 function W ZF(H, P) % H: K*Nr x Nt复数信道矩阵 % P: 总发射功率标量 W H / (H * H); % H 是共轭转置 % 功率归一化保证各流平均功率为 P/K W W ./ sqrt(mean(abs(W).^2, 1)) * sqrt(P / size(W, 2)); end逻辑说明H / (H * H)在 MATLAB 中等价于H * inv(H * H)但用右除可以避免显式求逆数值更稳定。归一化那一行是关键mean(abs(W).^2, 1)计算每一行对应发射天线的平均功率然后缩放让总功率等于P。如果注释掉归一化sum-rate 会随着 Nt 增大而虚高因为发射功率等于被放大了。参数说明P是总功率一般设置为K * Nr因为每个接收天线期望信噪比为 0dB 时功率维度要和噪声方差匹配。MF 的实现更简单把W H直接返回即可但归一化方式要和 ZF 保持一致否则两条曲线没有可比性。3.2 MMSE在干扰抑制与噪声放大之间取均衡点ZF 在高信噪比下性能很好但在低信噪比时求逆操作会把噪声分量放大。MMSE 的做法是在求逆矩阵里加一个对角加载项让解在“消除干扰”和“不放大噪声”之间取折中W_MMSE H^H (H H^H α I)^{-1}其中α σ_n^2 / σ_s^2也就是噪声方差除以信号方差。当 SNR 很高时α趋近于 0MMSE 会退化成为 ZF当 SNR 很低时α很大MMSE 会趋近于 MF 的方向。这个退化趋势经常被人当成代码 bug其实它是理论行为。% MMSE.m 核心代码 function W MMSE(H, P, sigma2) % sigma2: 噪声方差通常在 main 里按 10^(-SNR_dB/10) 计算 Nrx size(H, 1); % 总接收天线数 K*Nr alpha Nrx * sigma2 / P; % 对每根接收天线归一化的正则项 W H / (H * H alpha * eye(Nrx)); W W ./ sqrt(mean(abs(W).^2, 1)) * sqrt(P / size(W, 2)); end逻辑说明eye(Nrx)是单位阵alpha需要乘以接收天线数吗这是最容易混淆的地方。常见做法是alpha sigma2 / (P/Nrx)也就是每根接收天线上的信号功率。我一般直接写成alpha Nrx * sigma2 / P这样在总功率约束下高低 SNR 的过渡点会更准确。参数说明如果sigma2传错了比如传成10^(-SNR_dB/10)而忘了考虑符号和功率分配MMSE 曲线会明显偏离 ZF 和 MF 之间应有的位置。3.3 SLNR信漏噪比最大化与广义特征值分解SLNR 的出发点和 ZF/MMSE 不同它不追求完全消除对每个其他用户的干扰而是最大化“自己的信号功率”与“泄露给别人的干扰功率 噪声功率”之比。对用户 k定义干扰信道H_tilde_k为除用户 k 外所有用户的信道拼接SLNR ||H_k w_k||^2 / (||H_tilde_k w_k||^2 σ^2)最优的w_k是矩阵(σ^2 I H_tilde_k^H H_tilde_k)^{-1} H_k^H H_k的最大特征值对应的特征向量。% SLNR.m 核心代码 function W SLNR(H, K, Nr, sigma2) Nt size(H, 2); W zeros(Nt, K*Nr); for k 1:K idx (k-1)*Nr1 : k*Nr; % 当前用户的行范围 Hk H(idx, :); % 当前用户信道 idx_other setdiff(1:K*Nr, idx); Htilde H(idx_other, :); % 其他用户信道拼接 A Hk * Hk; % 信号空间 B Htilde * Htilde sigma2 * eye(Nt); % 干扰噪声空间 [V, D] eig(A, B); % 广义特征值分解 [~, pos] max(real(diag(D))); % 取最大广义特征值 w V(:, pos); % 对应最优方向 W(:, idx) w / norm(w); % 每流归一化 end W W * sqrt(P / sum(sum(abs(W).^2))); % 按总功率缩放 end逻辑说明eig(A,B)在 MATLAB 里解的是广义特征值问题A*v λ*B*v这正好是 SLNR 最优解的数学形式。注意B矩阵必须包含sigma2*eye(Nt)如果漏掉这一项SLNR 在低 SNR 时会失去噪声抑制能力曲线会和 ZF 重叠甚至更差。参数说明sigma2必须和 MMSE 里用的是同一个值不能一个按 dB 转换一个不转换。3.4 BD 与 SVD块对角化如何给多用户信道解耦BD块对角化的思路是先给每个用户构造一个干扰抑制矩阵使其能完全消除其他用户的信号。它依赖 SVD 分解。步骤是对用户 k先构造所有其他用户信道组成的矩阵H_tilde_k做 SVD 后取零空间基向量作为第一级预编码这样其他用户收到的来自用户 k 的干扰自然为零。然后再对等价信道做一次 SVD进行下一步功率分配。% BD.m 核心步骤以单个用户为例 function W_k BD_for_user(Hk, Htilde) % Hk: Nr x Nt目标用户信道 % Htilde: (K-1)*Nr x Nt其他用户信道 [~, ~, V_null] svd(Htilde); r rank(Htilde); B0 V_null(:, r1:end); % 零空间基向量 Hk_eff Hk * B0; % 等效信道 [U, S, V] svd(Hk_eff); W_k B0 * V(:, 1:size(S, 1)); % 完整预编码矩阵 end逻辑说明svd(Htilde)返回的三个矩阵中V_null的列是右奇异向量其中对应非零奇异值的列张成行空间剩余列张成零空间。rank(Htilde)计算的是有效行数如果 Htilde 不满秩r会小于行数零空间维度更大所以r1:end是常用写法。参数说明BD 需要知道每个用户的信道行范围所以BD.m的函数签名里必须有K和Nr。而svdprecoding.m可能是对这部分功能的封装方便 main 脚本统一调用。4. 跑通对比仿真参数设置与曲线读取方法主脚本不是上来就跑 10000 次蒙特卡洛而是先设定好固定的用户数、天线数和 SNR 范围。理解参数设置比看公式更重要因为这会直接影响曲线形状。4.1 天线数 Nt、用户数与 SNR 的默认配置怎么改从 fig 文件名可以确认资源默认是K4、Nr2、SNR5dBNt 从 4 变到 10。如果你关注误码率曲线那么main2.m里可能是固定Nt6或Nt8让 SNR 从 0 到 20dB 扫描。两套脚本的参数是独立的修改方式都是在文件开头找类似这样的配置段% main1.m 参数配置 K 4; % 用户数 Nr 2; % 每个用户的接收天线数 Nt_set 4:10; % 基站天线数扫描范围 SNR_dB 5; % 固定信噪比 5dB Monte_Carlo 2000; % 蒙特卡洛次数逻辑说明Nt_set里的每个值都要生成一批信道做平均。蒙特卡洛次数太少曲线会抖动太多则运行时间翻倍。针对 5 种算法、7 个 Nt 点的组合2000 次大概需要几分钟到几十分钟取决于 MATLAB 版本和 CPU。参数说明如果你想把 SNR 变成变量、Nt 固定就直接交换SNR_dB和Nt_set的作用域在循环外部固定Nt内部遍历 SNR 数组。4.2 sum-rate 与 BER 曲线分别暴露算法哪类问题sum-rate 反映的是平均可达速率公式为R sum_k log2(det(I H_k W_k W_k^H H_k^H / σ^2))实际计算时如果直接用W而忘记对每个用户的功率流归一化sum-rate 会高得离谱。源码里应该在主循环中先计算等效信道H_eff H * W再根据信号维度和噪声方差计算速率。BER 曲线则更直观。QPSK 调制下误码率下降斜率反映了预编码的抗干扰能力。MF 在低 SNR 表现不错但高 SNR 会出现地板效应ZF 和 BD 在高 SNR 时曲线斜率更陡SLNR 通常介于 MMSE 和 ZF 之间。如果你看到 ZF 的 BER 比 MMSE 还好那基本可以断定是 MMSE 里的噪声方差没算对。这里给出一个二选一的判断方法固定Nt8、K4、Nr2把 SNR 从 0 调到 20dB看 MMSE 和 ZF 曲线应该在 SNR 大于 10dB 后逐步汇合而 MF 的误码率在 15dB 后明显变平。如果你的结果没有这个趋势优先检查功率归一化是不是给每个求流都用了同样的平均功率。% 检查各算法等效信道增益的辅助脚本 H (randn(8, 8) 1j*randn(8, 8)) / sqrt(2); [Wz, Wm, Wmf, Ws, Wb] ... get_all_precoders(H, 4, 2, 1); % 示意函数 for W {Wz, Wm, Wmf, Ws, Wb} H_eff H * W{1}; gains diag(H_eff * H_eff); fprintf(每个流的增益范围: %.3f ~ %.3f\n, ... min(real(gains)), max(real(gains))); end逻辑很简单在对称信道下五个算法给出的等效信道增益应该都在接近 1 的量级。如果某个流增益只有 1e-3那说明这个流被功率归一化吃掉了要么是零空间方向取反要么是sqrt(P/trace(W*W))的缩放方式写错了。5. 常见问题与避坑五个最容易让仿真结果报废的细节以下问题都是我在复现同类预编码对比时踩过的坑每一条都足以让 sum-rate 或 BER 曲线形状发生肉眼可见的变化。5.1 功率归一化方式不一致导致曲线整体平移现象单独仿真 ZF 时 sum-rate 看起来没问题但把 MF 和 ZF 放在同一张图上时MF 整体比 ZF 低 3dB无论怎么调 SNR 都无法重合。原因MF 的实现里用了W W / norm(W, fro)而 ZF 用了W W / sqrt(mean(abs(W).^2))两种归一化导致的发射功率不同。解决强制所有算法的预编码矩阵在返回前都执行同一条归一化语句最好是W W ./ sqrt(sum(abs(W).^2, 1))先按流归一化再乘以sqrt(P / size(W, 2))。这样每个数据流的发射功率完全相同比较才有意义。5.2 SLNR 漏掉噪声方差曲线与 ZF 完全重合现象SLNR 的 sum-rate 曲线几乎和 ZF 重叠完全没有体现出它在低信噪比下的优势。原因在构造B Htilde * Htilde sigma2 * eye(Nt)时把后半部分注释掉了或者sigma2被赋值为 0。解决检查主脚本里sigma2 10^(-SNR_dB/10)是否正确并且在SLNR.m中确保sigma2 0。一个快速验证方法把SNR_dB设为 -5dBSLNR 应该比 ZF 高出不少如果还是重合说明噪声项确实没参与计算。5.3 BD 预编码里用户排序顺序影响最终曲线现象固定随机种子main1.m跑出来的 BD 曲线每次都有细微抖动尤其当 Nt 接近用户总天线数时。原因BD 算法要逐个用户处理如果信道状态较差的用户排在前面它的零空间基会占掉一些自由度削弱后面用户的性能。这属于数值排序导致的结果差异。解决在主循环里先按信道增益降序排序用户比如计算每个用户的norm(Hk, fro)从大到小排列后再做 BD。这样可复现性更好而且更贴近实际调度策略。5.4 注水算法传参错误导致 sum-rate 比平均功率分配还低现象waterfilling.m返回的功率分配向量出现负数甚至会报NaN。原因注水算法需要输入的是矩阵奇异值的平方也就是特征值λ_i s_i^2。如果直接把 SVD 返回的s传进去实际传的是奇异值导致水位线高于最大特征值部分通道得到负功率。解决在调用注水前执行lambda s.^2;。如果希望用更稳妥的写法可以直接用eig(chan * chan)拿到非负特征值。5.5 高 SNR 下 MMSE 和 ZF 曲线重合后误以为代码有 bug现象main2.m中 BER 曲线在 SNR 大于 20dB 时MMSE 和 ZF 的误码率完全相同曲线重叠。原因这是理论上的预期行为。alpha随着 SNR 增大而趋于 0MMSE 退化为 ZF。如果 SNR 达到 30dB 还看到明显差异反倒是alpha计算里多了个Nrx因子。解决不用解决但要写明注释避免复现者误删。在MMSE.m里加一行% Note: at high SNR, MMSE - ZF intentionally能帮你将来省下半天排查时间。6. 用现成 FIG 反向核验快速定位算法写对没有的土办法这份资源附带的.fig文件不只是用来展示结果的它们是最直接的调试参照。假设你只改了SLNR.m重新跑main1.m得到的新曲线到底可信不可信最简单的验证方式是直接在 MATLAB 里打开sum-rate8.fig用数据游标量一下特定 Nt 位置的数值再和你新跑的结果对比。如果 Nt8 时 SLNR 曲线的数值差超过 0.5bps/Hz那说明改动可能改变了功率分配而不是算法本身。更系统的方法是把 fig 中的数据提取出来做数值比较fig_list dir(sum-rate*.fig); for i 1:length(fig_list) openfig(fig_list(i).name, invisible); ax gca; lines findobj(ax, Type, Line); % 从第一条线提取 y 数据对比当前仿真的平均 sum-rate y_data get(lines(1), YData); fprintf(%s: Nt%d, sum-rate%.3f\n, ... fig_list(i).name, ax.XData(end), mean(y_data)); end逻辑说明fig文件保存的是画图后的图形对象findobj(ax, Type, Line)能拿到所有曲线的句柄。这里有两个需要注意的地方。第一gca在openfig(..., invisible)后可能被新图覆盖所以要用openfig(..., invisible)并手动指定图窗。第二XData可能是向量打印时最好取最后一个点或者整条曲线的平均值。我一直习惯在跑蒙特卡洛仿真前先用一个固定随机种子生成信道矩阵就能知道各算法的精确输出。比如下面这段rng(42); H (randn(8, 6) 1j*randn(8, 6)) / sqrt(2); W ZF(H, 1); % 手动计算 sum-rate和 fig 对照固定随机种子可以保证每次调试的信道完全一致配合 fig 中的现成曲线能快速锁定是算法问题还是随机信道波动。这个习惯帮我省了很多无意义的怀疑尤其是 SLNR 和 BD 这类依赖矩阵分解的实现小错误很容易被蒙特卡洛平均掩盖。从那以后我每次复现预编码对比都会先做三件事核对所有预编码矩阵的功率总和、跑一遍固定随机种子的单次仿真、再和附带的 fig 数值逐一对比。确认这三点没问题后才会放心地把 Nt 范围调大、SNR 调高去跑完整的曲线图。希望帮到你。本文还有配套的精品资源点击获取