ARTICLE · INTELLIGENCE

战地情报 · 详情页

来自尧图项目组的一线实战观察与深度解析

IRS辅助上行NOMA链路仿真代码详解:信道建模与SIC解码实践

IRS辅助上行NOMA链路仿真代码详解:信道建模与SIC解码实践

简介:面向IRS辅助NOMA上行链路研究者的论文复现代码包,围绕Sum Rate Maximization for IRS-Assisted Uplink NOMA图2与图3结果进行仿真,适合通信工程、信息与信号处理方向的硕博生、科研人员用于算法验证、曲线对比和课堂研讨。压缩包内共12个文件,其中8个Matlab脚本涵盖信道生成、IRS相位设计、NOMA功率分配、SDR优化求解及OMA基准对比等关键环节,4个fig图形文件可直接查看迭代收敛与性能变化趋势。整包大小约129KB,轻量易用,下载后配合主脚本运行即可快速产出论文对应结果。目前已有584人学习下载,属于典型的小型复现类学习资源。通过运行两个主脚本可分别得到图2与图3,帮助读者从代码层面理解IRS辅助NOMA和速率最大化的建模思路与实现细节,节省自行推导与编程调试的时间。

1. 这套IRS辅助上行NOMA链路代码,到底在仿真什么

仿真跑到第三天,我一度以为是 SIC 解码的索引写错了——同一套 IRS 辅助上行 NOMA 链路代码,换个相移初始化方式,两条用户速率曲线就交叉得完全没法解释。最后定位到的问题不在解码顺序,而在信道生成时把级联信道的方向搞反了。这类仿真,表面上是把 IRS 相移和 NOMA 功率分配叠起来算用户速率,实际上几乎所有的坑都藏在信道维度、线性/对数单位换算和解码顺序这三处。这套 IRS 辅助上行 NOMA 链路代码,就是用来量化智能反射面对上行非正交多址接入系统速率增益的蒙特卡洛仿真框架,适合做物理层算法预研、课程设计和毕业设计的开发者。它不碰硬件原型,但能在你投入真实验证平台之前,先把 IRS 的增益边界、NOMA 功率分配策略和相移量化的影响摸清楚。

2. 先理解仿真对象:IRS辅助上行NOMA的链路模型与三个核心模块

2.1 链路里到底有几条信道:直接链路与级联信道的建模差异

上行 NOMA 场景里,K 个单天线用户同时向基站发射信号,基站附近部署一个带有 N 个反射单元的 IRS。每个用户到基站有两条路径:一条是直接信道(用户→基站),另一条是反射信道(用户→IRS→基站)。反射信道是两跳信道的乘积,不能像直接信道那样一次性生成完事,代码里必须把「用户到 IRS」和「IRS 到基站」分开建模,再合成等效级联信道。

维度约定是第一个容易翻车的地方。我一般统一用这套约定:直接信道 h_d 存成 1×K 行向量,用户到 IRS 的信道 h_ru 存成 K×N 矩阵,IRS 到基站的信道 h_br 存成 N×1 列向量,相移向量 phi 也是 N×1 列向量。这样等效反射信道可以写成(h_ru .* phi.') * h_br,得到 K×1 再转置成 1×K,和 h_d 直接相加就是用户 k 的完整等效信道。

符号维度含义
h_d1×K用户直接到基站的复信道
h_ruK×N各用户到 IRS 各反射单元的复信道
h_brN×1IRS 到基站的复信道
phiN×1IRS 各单元相移,单位模约束

衰落模型的选择也要区分两条链路。直接链路通常用瑞利衰落就够,路径损耗指数按城市宏蜂窝取 3.5 左右;级联链路因为 IRS 一般部署在建筑高处,视距概率高,建议用莱斯衰落,莱斯因子取 5~10 dB。代码里如果统一用瑞利,会把 IRS 的增益低估不少,因为少了视距分量的稳定贡献。

路径损耗方面,级联链路要算两段损耗的和:用户到 IRS 一段,IRS 到基站一段。初版代码可以合并成一个等效路径损耗 PL_r,这样参数少、跑得快;等后面要对比 IRS 部署位置时,再把两段拆开。小尺度衰落生成时注意每个复数元素的实部和虚部要独立同分布、方差各为 1/2,这样信道模平方的期望才是 1,路径损耗才能单独叠加。

2.2 无源波束成形与功率分配:代码里两个可调旋钮

IRS 的本质是一个无源反射阵列,它不能放大信号,只能对每个反射单元的信号做相位旋转,这就是「无源波束成形」。代码里每个单元的相移是一个单位模复数exp(1j*theta),整个 IRS 的相移向量 phi 必须满足逐元素模为 1 的约束。如果代码里出现幅度不为 1 的相移,仿真结果就会凭空多出增益或损耗,这是物理上不允许的。

最常见的启发式相移设计是相位对齐:取某个用户到基站的级联信道,对每个反射单元的相位取负,让所有反射信号在基站处同相叠加。这个做法不是最优的,但它是所有对比的基线。真正的最优相移需要联合优化功率分配,属于非凸问题,仿真初期先跑启发式,后面再上优化算法对比,否则你很难判断优化带来的增益是真实的还是 SIC 写错了。

NOMA 部分的旋钮是用户发射功率。上行 NOMA 和下行 NOMA 有本质区别:下行是基站总功率约束,基站给不同用户分配功率;上行是每个用户各自有功率约束,用户之间不能转移功率。代码里把所有用户的发射功率存成一个 1×K 数组,初版直接全部设成等功率,先把链路跑通,再改成功率分配策略。

2.3 代码组织方式:一个主脚本加三个函数文件

这套仿真不要写成一个 800 行的脚本,我习惯按功能拆成四个文件,放在同一个目录下:

  • main_irs_noma_ul.m:主脚本,负责参数设置、蒙特卡洛循环、结果累计和画图
  • gen_channels.m:生成直接信道和级联信道,返回 h_d、h_ru、h_br
  • design_phase.m:设计 IRS 相移向量,返回 phi
  • compute_rate_ul_noma.m:给定等效信道、功率和噪声,按 SIC 顺序计算各用户速率

函数化的理由很简单:信道生成、相移设计、速率计算这三个模块在调试阶段会被反复调用,你改参数、换相移策略、加信道估计误差时,只需要动对应函数的内部实现,主脚本不用改。每个函数头部我都习惯写三行注释:输入维度、输出维度、单位。IRS 辅助 NOMA 这种多维复数运算,维度一致性和单位一致性决定了仿真代码的命。

3. 从零搭好IRS辅助上行NOMA链路代码:参数表、信道生成与相移初始化

3.1 仿真参数先定死:一张可抄的参数表

参数不统一的代码最容易翻车。我建议先把下面这张表当作初始配置写进主脚本开头,跑通之后再按需调整。K 先取 2,双用户场景的 SIC 索引关系最直观,排查问题容易;N 取 32 起步,既能看到 IRS 增益,又不会因为复杂度卡住。

参数符号建议取值说明
用户数K2先跑双用户,SIC 索引不容易错
IRS 单元数N32单元数翻倍时理想对齐增益约 3 dB
用户发射功率P10 dBm每用户独立,代码内用线性功率存
带宽B1 MHz决定噪声底
噪声系数NF5 dB接收机底噪抬升
直接链路损耗PL_d100 dB用户到基站的平均路径损耗
级联链路损耗(合并)PL_r80 dB用户-IRS-基站两段合计
莱斯因子Rice10 dB级联信道视距分量占比
蒙特卡洛次数MC10000少于 5000 曲线抖动会很严重
相移量化q_bitinf(连续)先连续相移,后加 2 bit 量化对比

噪声功率的换算这里最容易出错。标准做法是用noise_power = 10^((-174 + 10*log10(B) + NF)/10)算出瓦特数,再和用户功率一起参与 SINR 计算。更省事的做法是把噪声功率归一化为 1,用信噪比直接反推用户线性功率,但无论哪种,都要在参数区注释里写清楚用的哪种约定,不然过两周回来看代码一定会懵。

3.2 生成信道:瑞利衰落与莱斯信道的代码写法

信道生成函数返回三个量:直接信道、用户到 IRS 的信道、IRS 到基站的信道。路径损耗在线性域计算,小尺度衰落部分要区分瑞利和莱斯两种写法。

function [h_d, h_ru, h_br] = gen_channels(params) K = params.K; N = params.N; % 路径损耗 dB -> 线性 PL_d = 10^(-params.PL_d_dB/10); PL_r = 10^(-params.PL_r_dB/10); rice_lin = 10^(params.Rice_dB/10); % 莱斯因子转线性 % 直接链路:瑞利衰落,实部虚部分别是方差1/2的高斯 h_d = sqrt(PL_d/2) * (randn(1,K) + 1j*randn(1,K)); % 用户->IRS:莱斯衰落,含视距分量 los_ru = exp(1j*2*pi*rand(K,N)); % 简化视距相位 nlos_ru = sqrt(PL_r/(2*(rice_lin+1))) * (randn(K,N) + 1j*randn(K,N)); los_ru = sqrt(PL_r*rice_lin/(rice_lin+1)) * los_ru; h_ru = nlos_ru + los_ru; % IRS->基站:同样的莱斯结构 los_br = exp(1j*2*pi*rand(N,1)); nlos_br = sqrt(PL_r/(2*(rice_lin+1))) * (randn(N,1) + 1j*randn(N,1)); los_br = sqrt(PL_r*rice_lin/(rice_lin+1)) * los_br; h_br = nlos_br + los_br; end

这段代码的关键在莱斯衰落的功率分配:nlos 分量除以(rice_lin+1)再开方,los 分量按rice_lin/(rice_lin+1)缩放,这样合成后信道平均功率恰好是PL_r,路径损耗和小尺度衰落互不干扰。如果这里写反,IRS 的平均增益会被低估一半以上。

参数说明:randn生成的实部虚部各有方差 1,所以除以sqrt(2)让每个分量方差变成 1/2,合成信道模平方期望为 1。莱斯因子建议先传 dB 值再在函数内转线性,参数区看着更直观。视距分量的相位这里用均匀随机量简化,严格仿真应该用阵列响应向量,但初版这个精度足够看趋势。

3.3 相移初始化:先对齐相位再看NOMA增益

相移设计函数只做一件事:把某个用户的级联信道相位取负,让各反射单元的反射信号在基站处同相叠加。

function phi = design_phase(h_ru, h_br, target_user) % 输入: h_ru (KxN), h_br (Nx1), target_user 对齐目标用户 % 输出: phi (Nx1) 单位模相移向量 cascade = h_ru(target_user,:).' .* h_br; % Nx1 各单元级联信道 phi = exp(-1j * angle(cascade)); % 相位取负实现对齐 % 单位模约束:abs(phi) 全为 1 end

cascade是逐元素乘积,不是矩阵乘法,含义是每个反射单元从目标用户到基站各自形成的级联信道。angle取出相位,exp(-1j * angle)把它旋转到零度,这样每个单元的反射信号到达基站时相位一致、直接相加。

参数说明:target_user选谁直接决定多用户公平性——对齐强用户,弱用户反射增益会很小;对齐弱用户,弱用户速率提升明显但和速率不一定最高。初版先固定对齐用户 1,跟随机相移做对比。随机相移的基准也很简单,phi = exp(1j*2*pi*rand(N,1))就行,做一个函数或直接在循环里生成都可以。

4. 把NOMA解码和速率统计写进链路:叠加编码、SIC顺序与蒙特卡洛循环

4.1 上行NOMA的叠加信号与SIC解码顺序

基站接收信号是所有用户发射信号经过各自等效信道后的叠加,表达式为 y = Σ sqrt(P_k) h_k s_k + n。每个用户的等效信道是直接信道加反射信道:h_eff(k) = h_d(k) + sum(h_ru(k,:) .* phi.' .* h_br.'),代码里用矩阵乘法一次算完。

上行 NOMA 的 SIC 解码顺序和下行相反。下行基站给弱用户分配更多功率,所以基站先解强用户信号;上行每个用户功率独立,常见做法是按信道增益升序解码,也就是信道最弱的用户最先解。原因是弱用户不先解的话,会被强用户信号完全压住,第一步判决就出错,后面误差传播更严重。

解码次序解码对象此时干扰来源
第 1 个信道增益最小的用户其余所有用户
第 2 个信道增益次小的用户除去已解码用户的其余用户
第 K 个信道增益最大的用户无(只剩噪声)

这个顺序不是唯一选择,有些方案会按速率公平性动态调整,但初版代码用「弱先解」最容易验证。SINR 计算时要注意:对当前解码用户,干扰只来自排在它后面的未解码用户,已经解掉的用户信号要从干扰项里剔除。

4.2 速率计算的完整函数:一个蒙特卡洛循环的骨架

速率计算函数接收等效信道、线性功率、噪声功率和解码顺序,输出各用户速率。这个函数把 SIC 的干扰逻辑封装好,主脚本只需要调用它。

function rate = compute_rate_ul_noma(h_eff, Ptx_lin, noise_power, order) K = length(h_eff); rate = zeros(1,K); for idx = 1:K m = order(idx); % 当前要解码的用户 inter = 0; % 干扰只来自还没解码的用户 for j = idx+1:K l = order(j); inter = inter + Ptx_lin(l) * abs(h_eff(l))^2; end sinr = Ptx_lin(m) * abs(h_eff(m))^2 / (inter + noise_power); rate(m) = log2(1 + sinr); % bps/Hz end end

这里用「累计剩余干扰」的写法,而不是先算完整信号再逐次相减,避免了完美 SIC 下信号重建的数值误差。如果有人把 SIC 写成了每次减去Ptx_lin(m)*h_eff(m)再继续,碰到信道条件差、SINR 为负的情况,速率会变成 NaN。

主脚本的蒙特卡洛骨架如下:

%% main_irs_noma_ul.m clear; clc; rng(20240618); % 固定种子,结果可复现 params.K = 2; params.N = 32; params.Ptx_dBm = 10; params.noise_dBm = -109; params.PL_d_dB = 100; params.PL_r_dB = 80; params.Rice_dB = 10; Ptx_lin = 10^(params.Ptx_dBm/10) / 1000 * ones(1, params.K); noise_power = 10^(params.noise_dBm/10) / 1000; MC = 10000; rate_acc = zeros(1, params.K); for mc = 1:MC [h_d, h_ru, h_br] = gen_channels(params); phi = design_phase(h_ru, h_br, 1); % 对齐用户1 h_eff = h_d + ((h_ru .* phi.') * h_br).'; % 1xK [~, order] = sort(abs(h_eff).^2, 'ascend'); rate = compute_rate_ul_noma(h_eff, Ptx_lin, noise_power, order); rate_acc = rate_acc + rate; end rate_avg = rate_acc / MC;

(h_ru .* phi.') * h_br这一步做了三件事:phi.'把 N×1 转成 1×N,点乘让每个用户的信道逐单元乘上相移,再乘h_br完成反射链路叠加,最后转置成 1×K 和h_d对齐。MATLAB 的隐式扩展需要 R2016b 之后的版本,如果用的是老版本,要手动repmat(phi.', K, 1)。

参数说明:Ptx_lin的单位是瓦特,noise_power也是瓦特,两者量纲一致才能算 SINR。sort(..., 'ascend')对应弱用户先解。速率累计除以 MC 就是遍历速率,横轴扫发射功率时,这个循环结构完全不用改,只动参数区。

4.3 对比基准怎么选:OMA、随机相移、全用户等功率

没有基准的仿真曲线没有说服力。最少要跑四个场景放在一张图里:对齐相移加 NOMA、随机相移加 NOMA、对齐相移加 OMA、无 IRS 直接链路。每个场景换一下design_phase的返回值和接入方式即可。

对比场景设置看什么
对齐相移 + 等功率 NOMAphi 对齐用户 1,Ptx 等功率本文主推场景
随机相移 + NOMAphi 随机生成IRS 无源波束成形的增益
对齐相移 + OMA每个用户独占总资源NOMA 相对 OMA 的频谱效率增益
无 IRS 直接链路去掉反射项,h_eff = h_dIRS 部署本身的价值

OMA 基准有个新手必踩的坑:每个用户只能分到 1/K 的资源,速率要乘 1/K,即rate_oma = (1/K) * log2(1 + Ptx*abs(h_d)^2 / noise_power)。很多人的 OMA 曲线比 NOMA 还高,就是忘了除资源份额。

横轴一般扫用户发射功率,从 0 dBm 到 20 dBm,纵轴画和速率或最弱用户速率。这时候能清楚看到两个问题:IRS 对齐带来的增益随功率怎么变化,NOMA 相对 OMA 的增益在哪个功率区间被放大或吞掉。

5. IRS辅助上行NOMA链路代码避坑:5 个高频翻车点与修正方法

5.1 翻车点1:相移向量维度写错,级联信道拼不起来

现象:运行时报矩阵维度不匹配,或者h_eff变成 N×K 的矩阵,画出来的速率曲线乱成一团。这事我遇到过不止一次,基本都是phi用了行向量、h_br是列向量,点乘时隐式广播把维度搞爆炸。

原因:h_ru是 K×N,h_br是 N×1,phi如果是 1×N,那么h_ru .* phi得到 K×N 没问题,但(h_ru .* phi) * h_br维度是 K×1,如果你忘了转置,后面和h_d的 1×K 对不上;如果phi用了 N×1,直接h_ru .* phi就会报错。

解决:统一约定phi永远存 N×1 列向量,h_br也是 N×1。在design_phase函数头部加断言,让问题在源头暴露,而不是等到主脚本里报错:

assert(size(phi,1) == N && size(phi,2) == 1, 'phi必须为Nx1列向量'); assert(size(h_br,1) == N && size(h_br,2) == 1, 'h_br必须为Nx1列向量');

5.2 翻车点2:功率用dBm直接参与线性叠加,SINR变成负数

现象:计算出来的 SINR 是负值,速率全是 NaN,曲线根本没法看。第一次跑 NOMA 仿真的人十个有八个栽在这。

原因:dBm 是对数量纲,两个 dBm 值相加相当于线性值相乘,不是叠加。把 Ptx_dBm 和 noise_dBm 直接带进Ptx*abs(h)^2 / noise,分母分子量纲混乱,结果完全失真。

解决:所有功率在进入公式前统一转成瓦特。dBm 转瓦特是P_watt = 10^(P_dBm/10) / 1000,只在最后画图、标纵轴刻度时才转回 dBm。主脚本里换算完加一行注释,避免下次手滑。

5.3 翻车点3:SIC把想要的信号一起消掉了

现象:强用户的速率算出来比弱用户还低,甚至低于 OMA 基准,曲线完全违背直觉。

原因:解码顺序写反了,从强用户开始解,且重建信号时把当前用户也减掉了。更隐蔽的原因是干扰项里包含了已解码用户,导致强用户 SINR 被恶意压低。

解决:算 SINR 用「剩余干扰」写法,当前用户的信号功率永不出现在干扰项里。解码顺序用排序索引取,不要手写order = [weak_user, strong_user]这种硬编码,因为信道排序每次蒙特卡洛实现都不一样。

% 正确:当前用户m不在干扰项里 inter = 0; for j = idx+1:K l = order(j); inter = inter + Ptx_lin(l) * abs(h_eff(l))^2; end

5.4 翻车点4:蒙特卡洛次数不够,相移增益看不出来

现象:对齐相移的曲线和随机相移交叉,甚至低于随机相移,看起来像是 IRS 没有用。这种问题最容易让人怀疑「无源波束成形是不是玄学」。

原因:相位对齐本身就是针对随机信道实现的,单次信道里对齐带来的增益方差很大,MC 次数少时随机信道没有充分平均,曲线抖动超过实际增益差。

解决:MC 至少 10000 起步,曲线不平滑就加到 50000。固定rng种子,确保每次跑同参数得到同一条曲线,否则你无法判断曲线抖动是算法问题还是随机性。调参阶段先用 N=16、MC=2000 快速验证趋势,确认无误再放大参数。

5.5 翻车点5:IRS单元数N一大就跑不动,复杂度爆炸

现象:N=64 还能忍,N=256 时一个蒙特卡洛循环跑几分钟,内存也不断攀升。

原因:有人在相移设计里写了for n = 1:N逐单元循环,或者把每次信道实现都存进一个大矩阵,MC 还没跑完内存先爆了。

解决:信道生成和相移设计全部矢量化,不要对 IRS 单元写循环。主循环里只累加速率,不存中间信道矩阵,内存占用和 MC 次数无关。N 很大时可以把 MC 分块,每块 1000 次,块结束清空临时变量,防止 MATLAB 内存碎片化。

6. 让仿真结果更可信:验证曲线、基准对比与进阶调试习惯

拿到第一版曲线先别急着加复杂算法,先做退化验证。把 K 设成 1,NOMA 退化成单用户上行链路,此时没有 SIC,信道只剩h_eff = h_d + sum(|cascade|)。固定一个信道实现,手算Ptx * abs(h_eff)^2 / noise_power对应的速率,和代码打印值对比。对不上就逐个检查信道生成和相移函数,对得上再去跑多用户。

第二个习惯是固定信道快照调试。把蒙特卡洛循环关掉,只跑一次gen_channels→design_phase→compute_rate_ul_noma,打印每个中间变量的disp。重点核对三处:h_eff 的实部虚部量级是否和路径损耗匹配,排序后的 order 是不是符合信道强弱直觉,SINR 公式里干扰项拉出来手算一遍。这套操作能过滤掉八成以上的逻辑错误。

验证通过后,进阶方向有三个:相移量化、信道估计误差、用户配对。相移量化是把连续phi改成 2 bit 或 1 bit 量化,观察 IRS 增益损失,这个改动很小但非常值得做;信道估计误差是在 h_eff 上叠加一个百分比误差,用来评估 IRS 辅助系统的鲁棒性;K 大于 2 时还要加用户配对策略,不同用户组合对 NOMA 增益影响很大。

我自己的习惯是每个函数文件头部固定写三行注释:输入维度、输出维度、单位。IRS 辅助 NOMA 这种到处是复数矩阵、功率换算的仿真,维度一致性和单位一致性就是代码的生命线。这个习惯帮我少走了很多弯路,也希望帮到你。

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

RELATED READING

延伸阅读

更多一线实战笔记与深度复盘,助您持续精进