智能反射面安全传输复现：N变化曲线代码与避坑指南

简介：面向无线通信安全与智能反射面（IRS）研究的论文复现资料包，聚焦“随N变化”部分，重点解决IRS单元数量N对安全速率、信号聚焦与抗窃听性能影响的仿真分析问题，适合毕业设计、科研复现及课程实验场景。压缩包共7个文件，大小约87KB，主要包含4个mlx实时脚本、1个m函数文件、1个fig结果图与1个png示意图，覆盖均匀线性阵/均匀矩形阵的阵列响应、SDR安全传输率计算、N值扫描等核心模块。这份资料可直接在MATLAB中运行，通过调整N观察系统安全性能的变化规律，帮助理解相位优化与波束成形机制；同时m函数可复用，fig和png便于快速查看关键结果。已有2147人学习下载，对需要快速复现IRS安全通信论文结果、梳理仿真代码结构的读者有较高参考价值。

1. 智能反射面安全传输复现：N 的变化是调试入口，不是画图参数

复现智能反射面（IRS）在无线通信安全领域应用论文时，我习惯跳过花哨的算法部分，先把“保密速率随反射单元数 N 变化”那一节跑出来。原因很简单：N 是 IRS 所有性能红利的源头，这根曲线如果不成立，后面的相移优化、功率分配、用户调度全都建立不起来。很多同学把 N 当作横坐标，随便生成信道就画图，结果曲线不增反降，或者抖动到没法看。

这篇笔记会从系统模型、信道生成、保密速率计算到 N 扫描的完整代码，以及我踩过的 5 个坑，按可复现的步骤讲清楚。适合正在复现类似工作、卡在曲线对不上的研究生和无线通信工程师。如果你只想要一段能跑通的 N 扫描代码，可以直接看第 3 章；如果你想让曲线和论文一样平滑，请重点读第 4 章和第 5 章的避坑记录。

2. IRS 安全传输的系统模型：N 是通过什么路径进入保密速率的

开始写代码之前，我建议先花半小时把 N 在保密速率表达式里的位置摆清楚。IRS 可以看成一块由 N 个无源反射单元组成的面阵，每个单元只改变入射信号的相位和幅度，不发射信号。发射端到接收端之间除了直达链路，还有一条“发射端→IRS→接收端”的反射链路。因为 IRS 没有射频链路，所以它带来的增益来自相位对齐后的相干叠加，而这个叠加的次数正好是 N。这一章我会沿着系统模型、保密速率公式、信道模型三条线展开，最终解释为什么 N 越大，保密速率近似按 log2 N 增长。

2.1 从发射端到接收端：IRS 的信道作用链路

考虑最常见的单天线安全传输场景：发射端 Alice 单天线，合法接收端 Bob 单天线，窃听端 Eve 单天线，IRS 部署在 Alice 附近或 Bob 附近。令 h1 表示 Alice 到 IRS 的 N 维信道向量，h2 表示 IRS 到 Bob 的 N 维信道向量，g2 表示 IRS 到 Eve 的 N 维信道向量。IRS 的反射系数矩阵是 Θ = diag(e^{jθ1}, e^{jθ2}, ..., e^{jθN})，这里先假设每个反射单元幅度为 1，只调相位。合法接收信号可以写成：

y_b = ( h_d + h1^H Θ h2 ) s + n_b

其中 h_d 是直达链路的复标量，h1^H Θ h2 = Σ_{n=1}^N conj(h1[n]) e^{jθn} h2[n]。窃听者的接收信号结构完全一样，只要把 h_d、h2、n_b 换成 g_d、g2、n_e。这里 h1^H 是 h1 的共轭转置，所以 IRS 分量的每一项都包含一个共轭项，这个共轭项是后面相位公式最容易出错的地方。

IRS 分量的每一项都是三个复数值相乘：conj(h1[n])、e^{jθn}、h2[n]。三个值乘在一起后的相位是 -arg(h1[n]) + θn + arg(h2[n])。如果 IRS 每个单元都用随机相位，这 N 项就是一次随机游走，合法接收功率期望与 N 成正比；如果 θn 能抵消掉那个相位，让每一项都是正实数，那么 N 项同相相加，幅度期望与 N 成正比，功率就与 N 成正比。这就是 IRS 增强安全传输的核心：合法链路可以被波束成形，窃听链路因为不知道合法链路的相位信息，仍然停留在随机游走水平。

实际复现时，还需要确认论文是否包含直达链路。很多 IRS 安全论文为了突出 IRS 的作用，假设合法链路和窃听链路都被阻挡，即 h_d = g_d = 0。这种情况下曲线最干净，也最适合先跑通。如果论文包含直达链路，最优相位不能简单逐项对齐，因为还得分一部分相位去贴合 h_d 的方向，通常要把目标函数拆成两部分，或者用迭代算法搜索。我的建议是先跑无直达版本，确认 N 的曲线形状正确后再加直达链路，否则问题混在一起，很难定位。

2.2 保密速率公式：N 在容量表达式中的两个位置

安全传输里最常用的评价指标是保密速率，不是误码率。因为安全问题关心的是窃听者能否解码消息，而不是本端是否解码错误。对于加性高斯窃听信道，瞬时保密速率定义为合法信道容量与窃听信道容量之差的正值部分：

R_s = max(0, log2(1 + γ_b) - log2(1 + γ_e))

其中 γ_b 和 γ_e 分别是 Bob 和 Eve 的接收信噪比。如果发射功率为 P，Bob 处噪声功率为 σ_b²，Eve 处噪声功率为 σ_e²，那么：

γ_b = P |h_d + h1^H Θ h2|² / σ_b²
γ_e = P |g_d + h1^H Θ g2|² / σ_e²

N 在这里出现了两次：一次在合法链路的绝对值平方里，一次在窃听链路的绝对值平方里。假设没有直达链路，并且所有信道元素独立同分布复高斯 CN(0,1)，当 IRS 相位对齐 Bob 时，统计意义上 |h1^H Θ h2|² 的均值与 N² 成正比，而 |h1^H Θ g2|² 的均值只与 N 成正比，因为相位是参照 h2 设计的，g2 无法享受对齐收益。

把这两个量代入保密速率，高信噪比下平均保密速率大约等于 log2(1 + SNR·N²) 减去 log2(1 + SNR·N)。当 N 足够大时，这个式子趋近于 log2 N。这条结论非常重要，它就是你第 4 章里校验仿真曲线是否正确的理论基准。如果跑出来的曲线斜率不符合 log 增长，而是接近线性增长，说明合法链路没有获得 N² 增益；如果完全不增长，说明相位公式直接写错了。

还有一点要提醒：有些论文用安全中断概率替代平均保密速率，同样画出随 N 变化的曲线。安全中断概率定义是瞬时保密速率小于某门限 R_th 的概率。这种情况下 N 会出现在中断概率的表达式里，曲线形态和平均保密速率不同，但背后的信道模型和相位设计方法完全一样。复现前先看清论文的指标，不要拿平均速率的理解去套中断概率的图。

2.3 为什么不能直接用随机生成的 N 个独立信道？信道去相关假设

前面推导依赖“各个反射单元信道独立同分布”这个假设。实际 IRS 单元间距如果不够大，单元之间的信道是相关的。比如一个 32 单元的线阵，单元间距只有半波长时，相邻单元相关系数可能到 0.5 以上，这时 N 项叠加不会产生理想的 N² 增益，曲线斜率会变缓。复现论文时，必须先确认论文用的是理想独立信道模型，还是相关信道模型。大多数论文为了突出 IRS 的增益，会采用独立瑞利衰落，并且写明“N independent reflecting elements”。如果论文没有特别说明，默认是独立同分布，仿真时直接生成随机信道就行，但心里要清楚这是理想假设。

如果你想验证相关性对 N 曲线的影响，可以构造一个空间相关矩阵。常见做法是用相邻单元相关系数 ρ 构造托普利兹矩阵 R，R 的 (i,j) 元素为 ρ^{|i-j|}。给定独立复高斯向量 h_iid，相关信道可以表示为 h_corr = R^{1/2} h_iid。在 Python 里可以用 scipy.linalg.sqrtm 计算矩阵平方根。但我建议不要把相关信道放在第一次复现里，因为相关信道会压低 N² 增益，曲线抖动也更明显，你反而分不清是模型差异还是代码 bug。

还要注意信道功率归一化。我生成 h1 用 (randn + j·randn)/sqrt(2)，这样每个元素的实部和虚部方差各为 1/2，总功率为 1，等价于 CN(0,1)。有些论文用 CN(0,1/2) 表示同相和正交分量方差各为 1/2，总功率也是 1，这两种写法都行。但如果你生成复高斯时没有除 sqrt(2)，总功率就是 2，相当于整条链路 SNR 抬高了 3dB。后面遇到与他人数值对不上时，回来检查这一步，往往能立刻发现问题。

3. 复现“随 N 变化”曲线：最小可运行代码与逐步解释

有了第 2 章的公式，写代码就很顺了。我用 Python 而不是 MATLAB，主要是因为核函数只有几十行，Python 免费、矩阵运算和绘图都顺手。如果你想移植到 MATLAB，逻辑完全一样，把 randn 的写法换一下即可。这一章给的是完整的最小脚本，目标不是代码量，而是让你能在十几分钟里跑出趋势正确的 N 曲线。我会把每一段代码都拆开讲，避免你复制完也不知道为什么这样写。

3.1 环境与脚本结构：先定 4 个常量

脚本依赖 numpy 和 matplotlib，Python 3.8 以上即可。第一步定义四个常量：信噪比、N 列表、蒙特卡洛次数、发射功率。我建议把发射功率固定为 1 W，噪声功率由 SNR 反推，这样代码最直观，也不会出现功率单位混乱。先看代码：

import numpy as np import matplotlib.pyplot as plt # 固定参数：先跑出趋势，再按论文调整 SNR_dB = 20 # 平均信噪比，单位 dB P = 1.0 # 发射功率，归一化为 1 W sigma2 = 10 ** (-SNR_dB / 10) # 噪声功率 W，由 SNR 反推 N_list = [16, 32, 64, 128, 256] # 反射单元数扫描点 mc_iters = 5000 # 每个 N 点的蒙特卡洛次数

SNR_dB 是平均信噪比，定义为 10log10(P/sigma2)。固定 P=1 后，sigma2 由 SNR_dB 直接反推。如果你需要模拟噪声功率随带宽变化，可以把 sigma2 改成带内噪声功率计算值，但这里保持简单。N_list 建议从 16 开始，因为 N 太小时 IRS 增益很弱，保密速率经常被截断成 0，画出来是一条平线，看不出规律。mc_iters 最开始可以用 1000 跑通流程，确认曲线趋势后再升到 5000 或更高。

3.2 生成信道并计算保密速率：核心函数

核心函数做三件事：生成 N 组独立的复高斯信道、计算 IRS 最优相位、计算单次信道实现下的保密速率。下面给的是连续相位版本，离散相位后面小节扩展。代码：

def secrecy_rate_single(N, P, sigma2, phase_quant_bits=None): """ 计算一次信道实现下的保密速率 phase_quant_bits: None 表示连续相位；1 表示 1bit，2 表示 2bit """ # 生成 Alice->IRS, IRS->Bob, IRS->Eve 信道 # 除以 sqrt(2) 使每个元素方差为 1，即 E[|h|^2]=1 h1 = (np.random.randn(N) + 1j * np.random.randn(N)) / np.sqrt(2) h2 = (np.random.randn(N) + 1j * np.random.randn(N)) / np.sqrt(2) g2 = (np.random.randn(N) + 1j * np.random.randn(N)) / np.sqrt(2) # IRS 相位对齐 Bob：目标让每个复乘积项同相 # 合法链路每一项是 conj(h1[n]) * exp(1j*theta_n) * h2[n] theta = -np.angle(h1) - np.angle(h2) # 可选：离散相位量化 if phase_quant_bits is not None: step = 2 * np.pi / (1 << phase_quant_bits) theta = np.round(theta / step) * step # IRS 在 Bob 和 Eve 处形成的等效信道 beam_b = np.sum(np.conj(h1) * np.exp(1j * theta) * h2) beam_e = np.sum(np.conj(h1) * np.exp(1j * theta) * g2) # 信噪比 gamma_b = P * np.abs(beam_b) ** 2 / sigma2 gamma_e = P * np.abs(beam_e) ** 2 / sigma2 # 保密速率 rate = np.log2(1 + gamma_b) - np.log2(1 + gamma_e) return max(0, rate)

theta 的计算是这段代码最容易出错的位置。合法链路 beam_b 的每一项是 conj(h1[n]) * exp(jθn) * h2[n]，把这个复数的相位写成 -arg(h1[n]) + θn + arg(h2[n])。为了让每一项变成正实数，θn 需要等于 arg(h1[n]) - arg(h2[n])。上面代码里的 -np.angle(h1) - np.angle(h2) 等价于 -angle(h1) - angle(h2)，因为 conj(h1[n]) 的相位是 -angle(h1[n])，所以这样写是对的。注意这里不能写成 np.angle(h1 * h2)，因为需要的是 conj(h1) 和 h2 乘积的相位，而不是 h1*h2 的相位，两者差了一个负号。

beam_e 也使用同一组 theta，因为 IRS 无法同时为 Bob 和 Eve 定制不同的相位。beam_e 的每一项是 conj(h1[n]) * exp(jθn) * g2[n]，θn 与 g2 无关，所以相位不会系统对齐，N 项相加接近随机游走。gamma_b 和 gamma_e 分别求出后，log2 差值如果为负，保密速率取 0，这是安全传输的自然结果：窃听者信道比合法信道还好时，这条链路就不应该传输秘密信息。

3.3 扫描 N 并绘制平均保密速率曲线

有了单样本函数，接下来就是对每个 N 做蒙特卡洛平均，然后画图。这一段代码很短，但它决定了“随 N 变化”的曲线形状。代码：

def run_N_sweep(N_list, mc_iters, P, sigma2): avg_rates = [] for N in N_list: # 多次信道实现取平均 rates = [secrecy_rate_single(N, P, sigma2) for _ in range(mc_iters)] avg = np.mean(rates) avg_rates.append(avg) print(f"N={N:4d}, average secrecy rate = {avg:.4f} bit/s/Hz") return np.array(avg_rates) avg_rates = run_N_sweep(N_list, mc_iters, P, sigma2) plt.figure(figsize=(7, 5)) plt.semilogx(N_list, avg_rates, 'o-', linewidth=2) plt.xlabel('Number of IRS elements (N)') plt.ylabel('Average secrecy rate (bit/s/Hz)') plt.grid(True, which='both', linestyle='--', alpha=0.6) plt.title(f'SNR = {SNR_dB} dB, {mc_iters} Monte Carlo runs') plt.show()

运行后你会看到曲线整体上升，而且 N 每增加 4 倍，保密速率大约增加 1 bit/s/Hz。这个现象对应第 2 章推导的 log2 N 增长。如果曲线不升反降，大概率是 theta 公式错了，直接跳到第 5 章排查。mc_iters=5000 在 N=256 时大概要跑十几秒到几十秒，取决于电脑性能。如果嫌慢，可以把 mc_iters 改为 1000 先跑通流程，确认趋势后再升上去。不要一上来就用 50000 次，否则你连调 bug 的耐心都会被消耗掉。

4. 参数设置与曲线校验：让“随 N 变化”曲线贴合论文趋势

代码能跑只是第一步，曲线要和论文吻合才算复现成功。这一章讲参数怎么设、怎么校验，以及什么情况下曲线看起来像论文但其实是错的。我强烈建议每跑完一轮就画一条理论近似线，不要只盯着仿真线看。很多看起来合理的平滑曲线，其实相位公式是错的，只是因为你把窃听者也一起优化了，导致曲线形状偶然接近论文。

4.1 参数表：SNR、蒙特卡洛次数、相位量化位数的合理区间

下面是我调过的参数范围和它们对曲线形态的影响，可以作为起点参考。

选 SNR_dB=20 是因为大多数 IRS 安全论文都喜欢用 10 或 20 dB 做对比图，20 dB 下曲线趋势更明显。如果你的论文背景是低信噪比，比如 0~5 dB，那 N 较小时保密速率可能全是 0，画出来是平的，容易误以为代码有 bug。N_list 不要取 [1,2,3,...]，因为对数横轴下小 N 挤在一起，看不出规律。常用做法是翻倍或乘 4：16, 32, 64, 128, 256，这样每个点之间的距离在对数坐标下相等。

4.2 用理论渐进线验证仿真结果：log2 N 和 N² 增益

在第 2 章我们推导过，无直达信道时，高信噪比下平均保密速率近似为：

R_s ≈ log2(1 + SNR·N²) - log2(1 + SNR·N)

其中 SNR 是线性值。N 很大时可以进一步简化为 log2 N。这条线不是精确解，但作为趋势校验足够用了。你可以把理论近似线画在同一张图上，代码：

SNR_linear = 10 ** (SNR_dB / 10) theo = [np.log2(1 + SNR_linear * N**2) - np.log2(1 + SNR_linear * N) for N in N_list] plt.semilogx(N_list, avg_rates, 'o-', label='Simulation') plt.semilogx(N_list, theo, 's--', label='Theoretical approx') plt.legend() plt.show()

原则上，仿真线应该略低于理论近似线。因为理论近似把随机量换成了期望值，而仿真每次都是随机实现；当 N 较小时，log 函数的凸性会让仿真均值低于近似值。如果仿真线反而明显高于理论线，或者斜率超过 log2 N，那很可能是你无意中优化了窃听链路，让窃听者也获得了不该有的增益，导致保密速率虚高。另一种快速自检是看 N=1 的结果：无直达时，单反射单元很难同时服务 Bob 并压制 Eve，平均保密速率通常接近 0。如果 N=1 时已经有 5 bit/s/Hz，那一定有问题。

4.3 曲线平滑与误差带：N 大时方差增大怎么办

蒙特卡洛次数相同的情况下，N 越大，瞬时保密速率的方差越大。因为 gamma_b 和 gamma_e 的分布更分散，取 log 后差值更容易出现极端值，所以曲线尾部经常有毛刺。解决办法有三个：增加 mc_iters、固定随机种子、用误差带展示统计特征。固定种子尤其重要，否则每次运行曲线都变，你很难判断参数调整是有效果还是随机波动。

代码示例：

np.random.seed(2024) # 固定种子，保证每次运行可复现 rates_matrix = [] for N in N_list: rates = [secrecy_rate_single(N, P, sigma2) for _ in range(mc_iters)] rates_matrix.append(rates) mean_rates = np.mean(rates_matrix, axis=1) std_rates = np.std(rates_matrix, axis=1, ddof=1) plt.semilogx(N_list, mean_rates, 'o-', label='Mean') plt.fill_between(N_list, mean_rates - std_rates/10, mean_rates + std_rates/10, alpha=0.3, label='Mean ± 1 standard error')

注意这里画的是标准误差（std/sqrt(mc_iters)），不是标准差。标准误差表示均值估计的置信度，标准差表示单次结果的离散程度。如果论文展示的是阴影带，一般用 5%~95% 分位数，而不是均值 ± 标准误差。我一般先固定随机种子，调参时曲线不会每次变样，最后出图前再关掉种子做多组平均。

5. 避坑：复现 IRS 安全仿真常见的 5 个问题

这一章是我最想写的内容。IRS 安全仿真看起来只有一个求和公式，但实际踩坑点非常多。下面五条都是我在复现过程中亲身经历过、并且给人排查过的典型问题，全部按“现象—原因—解决”展开，希望能帮你少走弯路。

5.1 现象：N 越大保密速率反而下降或剧烈波动

现象：N 从 16 增加到 128，平均保密速率没有上升，甚至下降，或者个别 N 点突然掉到 0。

原因：最常见的相位公式错误。很多初学者把 IRS 相位写成 theta = -np.angle(h1 * h2)，但没有对 h1 取共轭。正确合法链路项是 conj(h1[n]) * h2[n] * exp(jθn)，所以相位补偿需要 θn = angle(h1[n]) - angle(h2[n])。如果你写成 theta = -np.angle(h1) - np.angle(h2)，也是等价的。但如果写成 -np.angle(h1 * h2)，就会差一个负号，导致合法链路没有被对齐，反而变成随机散射。另一个常见原因是窃听链路用了独立随机相位，而不是合法链路同一组相位，导致窃听者也被“波束成形”了，保密速率被压低。

解决：先用 N=1 手算。单个单元时 beam_b = conj(h1) * exp(jθ) * h2，你应该选择 θ 让这个数变成 |h1||h2|，也就是正实数。如果 N=1 算出来还有虚部，符号肯定错了。再把 theta 公式代入 beam_e，确认 exp(jθ) 不会让窃听链路也同相叠加。

5.2 现象：曲线抖动大，每次运行形状不同

现象：同样的参数，跑两次曲线一根向上、一根向下，或者 N 越大抖动越明显。

原因：蒙特卡洛次数不够，且没有固定随机种子。N 较大时 gamma_b 和 gamma_e 的动态范围很大，log2 差对瞬时值敏感，500 次平均可能不够。此外，随机数生成顺序也会影响结果，比如先生成 h1 再生成 h2，与论文的顺序不同可能导致微小差异，但主要问题还是样本量不足。

解决：先固定随机种子：np.random.seed(42)，然后把 mc_iters 提到 5000 以上。如果你同时画多条曲线做对比，务必保证所有曲线使用同一批信道实现，否则比较没意义。做法是先生成一个大的信道矩阵保存下来，再按 N 截取前 N 列，而不是每个 N 重新生成。

5.3 现象：加了离散相位后，曲线斜率明显低于论文

现象：设定 phase_quant_bits=1 后，保密速率比起连续相位掉了不少，而且随 N 增长的斜率变小，与论文曲线相差很大。

原因：量化步长过大，或者量化方式错了。1bit 意味着相位只能取 0 和 π，相当于每个反射单元只有同相/反相两个状态，合法链路的相位误差会有随机残余，N 越大，增益损失越明显。有些论文在量化后用随机量化，有些用四舍五入，两者结果不同。如果论文没有明确说明量化方式，一般默认是取最近的有效相位。

解决：先确认论文的相位分辨率。如果是 1bit，量化公式应该写成 theta_q = np.round(theta / (2np.pi/2^bits)) * (2np.pi/2^bits)。对比连续相位曲线，确认损耗符合预期。注意：不要直接在连续相位代码上补量化，要先确认连续相位是复现正确的，再逐级加量化。

5.4 现象：相同参数但和别人/论文的数值对不上

现象：设置 SNR=20 dB、N=64，仿真结果是 3.5 bit/s/Hz，论文里是 4.5 bit/s/Hz，差 1 bit/s/Hz 左右。

原因：信道功率归一化不一致。我在 2.3 提到过 CN(0,1) 和 CN(0,1/2) 的区别。用 randn/sqrt(2) 生成的复高斯，实部和虚部各 1/2 方差，总功率为 1。但如果你用了 np.random.randn(N) + 1jnp.random.randn(N) 而没有除 sqrt(2)，总功率就是 2，相当于 SNR 抬高了 3dB。另一种可能是 SNR 定义不同：论文可能把平均信噪比定义为 PE[|h|²]/σ²，而你实际计算的是瞬时信噪比后再平均。

解决：把信道总功率打印出来：np.mean(np.abs(h1)**2)，确认是 1。再检查 SNR 换算，sigma2 = P / 10^(SNR_dB/10)。如果你用的是 dBm，要先减 30 转成 dBW，再反推线性功率。差 3dB 或 6dB 基本都是归一化问题，不要急着怀疑算法。

5.5 现象：仿真耗时太长，N=512 以上跑不完

现象：Python 里 for 循环跑 N=512、mc_iters=10000 需要几分钟，调参非常痛苦。

原因：每次循环都生成 N 个信道并做求和，Python 循环本身开销太大。其实计算量很小，瓶颈在解释器。

解决：向量化。把蒙特卡洛次数放在矩阵的第一维，一次生成 (mc_iters, N) 的复数矩阵，用 numpy 批量求和，速度能提升几十倍。代码：

def secrecy_rate_batch(N, P, sigma2, mc_iters): # 一次性生成全部信道样本 h1 = (np.random.randn(mc_iters, N) + 1j*np.random.randn(mc_iters, N)) / np.sqrt(2) h2 = (np.random.randn(mc_iters, N) + 1j*np.random.randn(mc_iters, N)) / np.sqrt(2) g2 = (np.random.randn(mc_iters, N) + 1j*np.random.randn(mc_iters, N)) / np.sqrt(2) theta = -np.angle(h1) - np.angle(h2) beam_b = np.sum(np.conj(h1) * np.exp(1j * theta) * h2, axis=1) beam_e = np.sum(np.conj(h1) * np.exp(1j * theta) * g2, axis=1) gamma_b = P * np.abs(beam_b)**2 / sigma2 gamma_e = P * np.abs(beam_e)**2 / sigma2 rates = np.log2(1 + gamma_b) - np.log2(1 + gamma_e) return np.maximum(0, rates)

注意内存。mc_iters=5000、N=1024 时，三个复数矩阵各占约 5000102416 字节，大约 78 MB，三个共 240 MB，还在承受范围内。如果 N 再大，把 mc_iters 降到 2000，或分块处理。向量化后 N=512、mc_iters=10000 大概几秒就能跑完，调参体验会好很多。

6. 进阶：让曲线更可信的验证技巧与扩展方向

到这里，你的“随 N 变化”曲线应该能稳定复现了。最后一章我讲一个每次都会用的验证技巧，以及两个值得扩展的方向。

6.1 用“小 N 枚举 + 大 N 理论线”做校准

我习惯在跑完整仿真前，先对 N=1、2、4 做一次相位穷举。连续相位下最优相位有闭式解，但小 N 可以暴力搜索：把每个单元的相位离散成 32 个测试值，遍历组合，N=4 时有 32^4 约 100 万种，可以接受。找出使保密速率最大的相位，和闭式解对比。如果闭式解不如暴力搜索，说明相位公式有 bug。小 N 校准通过后，再用第 4 章的理论线验证大 N 趋势。这个“先小后大”的流程能过滤掉 80% 的符号错误。

6.2 从平均保密速率扩展到安全中断概率

很多论文的“随 N 变化”部分不止一条平均速率曲线，还有安全中断概率曲线。定义门限 R_th，一次信道实现的瞬时保密速率小于 R_th 就算中断。扩展代码非常简单：把平均速率改成 np.mean(rates < R_th) 即可。建议同时画平均速率和安全中断概率两条线，因为它们对 N 的敏感度不同。中断概率曲线会在某个 N 出现明显拐点，这个拐点常被论文用来设计最小反射单元数。

6.3 把连续相位换成 1bit/2bit 量化并对比

最后，可以继续用 secrecy_rate_single 里的 phase_quant_bits 参数，在同一张图画出连续相位、1bit、2bit 三条 N 曲线。这个对比在实际工程里很有用，因为硬件实现的 IRS 相位分辨率通常有限。你会发现 2bit 量化与连续相位的差距远小于 1bit，这也是很多芯片设计选择 2bit 的原因。我每次复现新论文，都会先跑 N 的曲线，再跑任何一个优化算法；曲线不对就先别碰其余部分。这个习惯救过我很多次。希望帮到你。

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