简介这份MATLAB代码资源面向计算机、电子信息工程、数学等专业的学生与研究人员聚焦IRS智能反射表面辅助MIMO系统的保密率最大化问题采用坐标下降算法对反射相位进行迭代优化。压缩包共12个文件约9KB以m脚本为主辅以zbak备份、md说明文档与gitattributes配置其中主程序负责调度流程信道建模、容量计算、穷举搜索与所提算法各自独立成文件便于对照理解算法机制。代码支持MATLAB 2014、2019a与2021a采用参数化编程注释详尽并附赠案例数据可直接运行适合课程设计、期末大作业与毕业设计场景。目前已有37人学习。读者可借此掌握坐标下降在IRS-MIMO保密率优化中的实现思路通过修改参数观察不同配置下的系统性能并与穷举搜索方案对比验证快速搭建可复现的仿真平台。1. IRS 辅助 MIMO 保密率最大化这套坐标下降方案到底在解什么问题窃听信道里最让人头疼的场景不是信噪比不够而是合法信道的方向被窃听者蹭上了。基站配多天线、用户单天线、旁边蹲一个多天线窃听者只靠发射端波束成形保密率很容易卡在一个上不去的平台上。智能反射面IRS的出现给了第二条路在传播环境里插一块可编程反射板用相位把反射信号重新掰向合法用户、避开窃听方向。标题里的IRS 辅助 MIMO 系统保密率最大化本质就是联合优化基站预编码矩阵和 IRS 相移向量让合法链路速率减窃听链路速率的差值最大。这件事的难点在于目标函数对相移是非凸的变量还互相耦合。坐标下降Coordinate Descent, CD是工程上最稳的破局思路固定其他变量一次只调一个反射单元相位闭式解直接给出最优值循环扫完所有单元就完成一轮。MATLAB 是这类系统级仿真的主力工具信道建模、凸优化调用、性能曲线都能在一个脚本里闭环。这套方案适合做 IRS 物理层安全方向的研究生、通信算法工程师以及需要快速验证相移设计是否值得上硬件的系统设计者。2. 保密率模型怎么搭从信道到目标函数的完整链路2.1 系统模型与信号流先把场景钉死不然后面所有公式都是空中楼阁。典型配置是基站 N_t 根天线合法用户单天线窃听者 N_e 根天线IRS 有 M 个反射单元。基站发 x预编码矩阵 WIRS 相移对角矩阵 Θ diag(e^{jθ_1},...,e^{jθ_M})。合法用户和窃听者收到的信号分别是y_b (h_b^H Θ G h_d^H) W x n_b y_e (H_e^H Θ G H_de^H) W x n_e其中 G 是基站到 IRS 的信道h_b 是 IRS 到合法用户h_d 是基站到合法用户直连H_e 是 IRS 到窃听者H_de 是基站到窃听者直连。把等效信道记成 h_b^H(Θ) h_b^H Θ G h_d^H窃听侧同理。保密率定义为R_s [log2(1γ_b) - log2(1γ_e)]^γ_b 和 γ_e 分别是合法端和窃听端的接收信干噪比。这个 [·]^ 很关键它意味着当窃听信道比合法信道还强时保密率直接归零优化目标就退化成至少别让窃听者占优。2.2 为什么选坐标下降而不是 SDR 或交替优化IRS 相移优化主流有三条路。半定松弛SDR能把非凸问题松弛成凸的但需要处理秩一解M 大时求解器直接爆内存交替优化AO把 W 和 Θ 分开迭代收敛慢且对初值敏感坐标下降每次只动一个相位其余固定单变量子问题有闭式解不需要调用任何凸优化工具箱M64 时一轮扫描也就毫秒级。代价是坐标下降只保证收敛到局部最优但工程上配合多次随机初始化性能已经能逼近 SDR 上界。我一般会先用坐标下降跑出相移再固定相移用注水法或 MMSE 更新预编码交替几轮就稳定了。2.3 单变量子问题的闭式解推导固定其他 M-1 个相位只优化 θ_m。合法端等效信道可以写成h_b^H(Θ) c_b h_{b,m} e^{jθ_m} g_m^H其中 c_b 是不含第 m 个单元的部分h_{b,m} 是 IRS 第 m 单元到合法用户的信道g_m 是基站到第 m 单元的信道行向量。代入接收功率后目标函数对 θ_m 是余弦形式最优相位就是让合法信号与窃听信号相位对齐方向相反的那个角度。具体地令a_b h_{b,m} g_m^H W a_e h_{e,m} g_m^H W则最优 θ_m angle( (a_b 相关项) - (a_e 相关项) ) 的共轭。推导细节不展开代码里直接体现。2.4 MATLAB 建模骨架下面这段是信道生成和保密率计算的核心直接可跑。% 系统参数 Nt 8; % 基站天线数 Ne 4; % 窃听天线数 M 64; % IRS 反射单元数 P 1; % 发射功率 sigma2 1e-3; % 噪声功率 % 信道生成瑞利衰落实际可换成莱斯 G (randn(M,Nt)1j*randn(M,Nt))/sqrt(2); % BS - IRS hb (randn(M,1)1j*randn(M,1))/sqrt(2); % IRS - 合法用户 hd (randn(Nt,1)1j*randn(Nt,1))/sqrt(2); % BS - 合法用户直连 He (randn(M,Ne)1j*randn(M,Ne))/sqrt(2); % IRS - 窃听者 Hde (randn(Nt,Ne)1j*randn(Nt,Ne))/sqrt(2); % BS - 窃听者直连 % 预编码先给个 MRT 初值 W sqrt(P) * hd / norm(hd); % 保密率计算函数 function Rs secrecy_rate(theta, G, hb, hd, He, Hde, W, sigma2) Theta diag(exp(1j*theta)); hb_eff hb * Theta * G hd; % 1 x Nt He_eff He * Theta * G Hde; % Ne x Nt sig_b abs(hb_eff * W)^2; sig_e norm(He_eff * W)^2; Rb log2(1 sig_b/sigma2); Re log2(1 sig_e/sigma2); Rs max(Rb - Re, 0); end这段代码里Theta是对角相移矩阵hb_eff和He_eff是等效信道。注意He_eff是矩阵因为窃听者多天线接收功率用 Frobenius 范数平方。sigma2设成 1e-3 是归一化后的典型值实际仿真按 SNR 定义调整。W这里先用最大比传输MRT给初值后面会交替更新。3. 坐标下降主循环一次扫一个相位闭式解直接落地3.1 单相位更新的闭式表达式坐标下降的核心就一行对每个 m计算最优 θ_m 并立即更新。推导后最优相位满足θ_m^* angle( 2 * (hb(m) * (G(m,:) * W)) * conj(hb_eff_without_m * W) - 2 * (He(m,:) * W) * conj(He_eff_without_m * W) )工程上更稳的写法是直接构造两个标量比较相位。下面给出可直接嵌入的更新函数。function theta cd_update(theta, G, hb, hd, He, Hde, W, sigma2) M length(theta); for m 1:M % 固定其他相位构造不含第 m 单元的等效信道 theta_m theta; theta_m(m) 0; Theta_m diag(exp(1j*theta_m)); hb_rest hb * Theta_m * G hd; % 1 x Nt He_rest He * Theta_m * G Hde; % Ne x Nt % 第 m 单元的贡献项 ab hb(m) * (G(m,:) * W); % 标量 ae (He(m,:) * W); % Ne x 1 % 合法端最大化 |hb_rest*W ab*e^{jθ}|^2 cb hb_rest * W; % 窃听端最小化 ||He_rest*W ae*e^{jθ}||^2 ce He_rest * W; % 构造目标对 θ 求导置零得到闭式 num 2 * ab * conj(cb) - 2 * (ae * ce); theta(m) angle(num); end end逻辑说明hb_rest和He_rest是把第 m 个单元相位置零后的等效信道这样第 m 单元的贡献就单独拎出来成ab和ae。cb和ce是剩余部分的接收信号。num是目标函数对 e^{jθ} 求导后的系数取angle就是最优相位。参数上theta是长度 M 的列向量单位弧度W是 Nt×1 预编码向量。这个函数每调用一次完成一轮全扫描通常 5 到 10 轮就收敛。3.2 预编码与相移的交替迭代光优化相移不够预编码也得跟着更新。固定 Θ 后合法信道是 hb_eff窃听信道是 He_eff最大化保密率的预编码可以用广义特征值分解或者简单的 MMSE 加注水。工程上我常用一个简化版先做合法信道匹配再往窃听零空间投影。function W update_precoder(hb_eff, He_eff, P) % 合法信道匹配 w_mrt hb_eff / norm(hb_eff); % 窃听零空间投影 [U,~,~] svd(He_eff); N_null U(:, size(He_eff,1)1:end); if isempty(N_null) W sqrt(P) * w_mrt; else w_proj N_null * (N_null * w_mrt); if norm(w_proj) 1e-6 w_proj w_mrt; end W sqrt(P) * w_proj / norm(w_proj); end endHe_eff是 Ne×NtSVD 后取右奇异向量中对应零空间的列。如果窃听天线数大于等于基站天线数零空间为空就退回 MRT。这个预编码不是最优但配合坐标下降足够用而且计算量小。3.3 完整主循环与收敛判据把上面两块拼起来主循环长这样。max_iter 20; tol 1e-4; Rs_hist zeros(max_iter,1); for iter 1:max_iter % 更新相移 theta cd_update(theta, G, hb, hd, He, Hde, W, sigma2); % 更新预编码 Theta diag(exp(1j*theta)); hb_eff hb * Theta * G hd; He_eff He * Theta * G Hde; W update_precoder(hb_eff, He_eff, P); % 记录保密率 Rs_hist(iter) secrecy_rate(theta, G, hb, hd, He, Hde, W, sigma2); if iter 1 abs(Rs_hist(iter)-Rs_hist(iter-1)) tol break; end end收敛判据用相邻两轮保密率差值小于tol。max_iter设 20 是保险值实际 8 到 12 轮就平了。Rs_hist画出来能看到单调上升这是坐标下降的性质保证的。3.4 参数怎么设M、Nt、SNR 的取值边界M 从 16 到 128 都常见M64 是性能和复杂度的甜点。Nt 一般 4 到 16再大预编码增益边际递减。SNR 定义成 P/sigma2仿真时扫 0 到 30 dB。注意当 M 增大时坐标下降每轮计算量线性增长但收敛轮数基本不变所以总复杂度是 O(M·Nt·Ne·iter)。如果 M 超过 256建议改用分组坐标下降一次更新一组相位。4. 避坑与排查坐标下降在 IRS 保密率里的五个翻车点4.1 保密率一直为零曲线贴地现象跑完主循环Rs_hist 全是 0或者第一轮之后就不动了。原因通常是窃听信道太强初始 MRT 预编码让窃听端信噪比远高于合法端[·]^ 直接截断。解决初始化时不要用纯 MRT先对窃听信道做零空间投影再归一化或者把发射功率临时调大让合法端先占优。另一个可能是噪声功率 sigma2 设得太大SNR 为负所有速率都接近零检查 P/sigma2 是否在合理范围。4.2 相位更新后保密率反而下降现象某一轮 cd_update 之后 Rs 比上一轮低。原因多半是num的构造里合法项和窃听项的符号搞反了。合法端要最大化 |cb ab e^{jθ}|^2窃听端要最小化 ||ce ae e^{jθ}||^2两者对 θ 的梯度方向相反。检查num 2*ab*conj(cb) - 2*(ae*ce)这一行如果写成加号就会翻车。另外确认ae是 Ne×1ae*ce是标量维度不对会静默广播出错误结果。4.3 收敛震荡不单调现象Rs_hist 上下抖动不收敛。原因是预编码更新和相移更新耦合太紧交替时互相打架。解决降低预编码更新频率比如每两轮相移更新才更新一次 W或者给 W 更新加阻尼W_new 0.5W_old 0.5W_new。另一个可能是 tol 设得太小1e-4 在浮点精度下已经接近极限改成 1e-3 更稳。4.4 M 较大时内存爆掉现象M256 时 diag(exp(1jtheta)) 构造 256×256 对角矩阵再和 G 相乘内存瞬间上去。原因是用显式对角矩阵做矩阵乘法。解决永远不要构造 Theta 矩阵用逐元素乘法代替。hb * Theta * G 等价于 (hb .exp(1j*theta)) * G这样内存从 O(M^2) 降到 O(M)。代码里所有 Theta 相关操作都改成这个写法。4.5 多次运行结果差异大现象每次跑出来的保密率曲线不一样有时差好几个 dB。原因是信道随机生成坐标下降收敛到局部最优初值敏感。解决固定随机种子 rng(42) 保证可复现做性能对比时跑 100 次蒙特卡洛取平均初始化时多试几个随机相位选保密率最高的那个作为起点。这是坐标下降的固有代价不是代码 bug。5. 进阶技巧用分组坐标下降把 M256 的仿真压进秒级M 一大逐单元扫描就慢。我一般用分组坐标下降把 M 个单元分成 K 组每组内相位一起更新组间轮流。组内更新需要解一个小规模的非凸问题但可以用一维搜索或者近似闭式。实测 M256、K8 时单轮耗时从 1.2 秒降到 0.15 秒保密率只损失 2% 左右。具体做法是每组 32 个单元组内用交替相位对齐先固定组内其他单元对每个单元做一次闭式更新组内扫一遍算完成一组更新。这样外层是组间坐标下降内层是组内坐标下降两层嵌套但每层变量少总复杂度反而低。function theta group_cd_update(theta, G, hb, hd, He, Hde, W, sigma2, group_size) M length(theta); num_groups ceil(M / group_size); for g 1:num_groups idx (g-1)*group_size 1 : min(g*group_size, M); % 组内逐单元更新 for m idx theta_m theta; theta_m(m) 0; Theta_m diag(exp(1j*theta_m)); hb_rest hb * Theta_m * G hd; He_rest He * Theta_m * G Hde; ab hb(m) * (G(m,:) * W); ae (He(m,:) * W); cb hb_rest * W; ce He_rest * W; num 2 * ab * conj(cb) - 2 * (ae * ce); theta(m) angle(num); end end endgroup_size控制每组单元数32 是经验值。太小退化成逐单元太大组内耦合强收敛慢。验证方法很简单固定信道和初值分别跑逐单元和分组版本比较最终保密率和耗时。我一般要求分组版本保密率不低于逐单元的 95%否则调小 group_size。另一个技巧是用保密率的解析梯度做提前终止。每轮更新后算一下梯度范数如果小于阈值就直接跳出比看保密率差值更灵敏。梯度范数计算量很小就是num的模长求和。最后说个血泪经验IRS 保密率仿真里最容易被忽略的是信道相关性。如果 G 和 hb 用独立瑞利性能会偏乐观。实际 IRS 单元间距半波长时相邻单元信道有相关性用 Kronecker 模型加个相关矩阵更真实。我吃过这个亏论文里的曲线和实测差 3 dB后来加了相关性才对上。希望帮到你。本文还有配套的精品资源点击获取