频域鲁棒的完美匹配层边界积分方程法（PML-BIE）：攻克周期阵列散射中的瑞利-伍德（Rayleigh-Wood）异常难题

来源论文: https://arxiv.org/abs/2606.07971v1 生成时间: Jun 14, 2026 01:11

0. 执行摘要

在微纳光学、光子晶体、超材料设计以及量子化学中的电子散射理论中，周期性阵列结构的波动散射模拟一直居于核心地位。传统的边界积分方程法（BIE）通常依赖于拟周期格林函数（Quasi-Periodic Green’s Function, QGF），但在物理上十分关键的**瑞利-伍德异常（Rayleigh-Wood Anomalies, RW anomalies）**频点处，由于格林函数中晶格和（Lattice Sums）的非绝对收敛甚至直接发散，导致数值格式彻底失效。尽管基于自由空间基本解的窗口格林函数法（WGF）和完全匹配层边界积分方程法（PML-BIE）能够在一定程度上绕过复杂的QGF计算，但它们在RW异常频点或其邻域内依然面临着致命的精度崩塌。其根本原因在于：当散射波转变为掠入射的表面波（即传播模式变为渐逝模式的临界状态）时，PML截断误差的指数衰减率会趋于零，同时边界积分系统中的非物理（非传播）模式的投影方程变得极度病态。

针对这一极具挑战性的国际学术难题，本论文提出了一种一致高精度、频域鲁棒的改性PML-BIE方法。该方法创造性地将复拉伸PML截断技术与一种全新的**有限模式修正算子（Finite-Mode Correction）**相结合。通过引入互补的半无限区间积分分析，将尾部算子的贡献转化为闭合解析解，并在RW奇异点处通过洛必达法则（L’Hôpital’s rule）进行严格的解析延拓，消除了数值不稳定性。同时，利用拟周期边界条件的空间拓展技巧，有效规避了截断单元边界的“角点奇异性（Corner Singularity）”。

数值基准测试表明，本方法在极薄的PML厚度下（仅为 $T = 4\lambda$）即可实现接近机器精度的结果（能量平衡误差达到 $10^{-11}$ 量级），其精度和收敛速度相比于传统修正WGF方法提升了数个数量级，且其收敛行为完全独立于RW异常频点。这一成果不仅为纳米光子学器件的高频电磁仿真提供了高效、高精度的数学工具，也为量子化学中一维/二维周期晶格上的电子本征态散射及格林函数计算提供了全新的理论视角。

1. 核心科学问题，理论基础，技术难点，方法细节

1.1 物理与数学模型

考虑二维空间中沿 $x_1$ 方向具有周期 $\Lambda > 0$ 的无限周期阵列介质障碍物 $D_2 = \bigcup_{i\in\mathbb{Z}} \Omega_2^{(i)}$。定义外界介质区域为 $D_1 = \mathbb{R}^2 \setminus \overline{D_2}$。障碍物内外均填充均匀、各向同性的介质，分别由介电常数 $\epsilon_j$ 和磁导率 $\mu_j$（$j=1,2$ 分别对应 $D_1,D_2$）刻画。在角频率 $\omega > 0$ 下，介质 $D_j$ 中的波数为 $k_j = \omega \sqrt{\epsilon_j \mu_j}$。

假设一束平面波 $u^{\text{inc}}(x) = e^{i(\alpha x_1 - \beta x_2)}$ 自上方以入射角 $\theta^{\text{inc}} \in (-\pi/2, \pi/2)$ 入射，其中 $\alpha = k_1 \sin\theta^{\text{inc}}$ 为 Bloch 相位常数（Bloch wavenumber），$\beta = k_1 \cos\theta^{\text{inc}}$。总场 $u$ 在各个区域内满足亥姆霍兹方程：

并在障碍物边界 $\partial D_2$ 上满足如下透射边界条件（Transmission Boundary Conditions）：

其中 $\gamma^{\pm}_D$ 和 $\gamma^{\pm}_N$ 分别代表边界上的 Dirichlet 和 Neumann 迹；对于 TE 极化，比例系数 $\eta = \mu_1 / \mu_2$；对于 TM 极化，$\eta = \epsilon_1 / \epsilon_2$。

由于结构的周期性和平面波入射的对称性，总场 $u$（以及散射场 $u^{\text{sct}} = u - u^{\text{inc}}$）满足拟周期性（Quasi-periodicity）：

这允许我们将计算域限制在单个周期单元（Unit Cell）$U$ 内。如图 1 所示，定义两条无限长的平行虚边界：

其中 $\hat{f}(t)$ 是光滑标量函数，在 $|t| > t_0$ 时收敛于常数 $-\Lambda/2$。在单个周期单元内，散射场 $u^{\text{sct}}$ 必须满足向上和向下的瑞利展开放射条件（Rayleigh Expansion Radiation Condition）：

其中 $\alpha_n = \alpha + \frac{2\pi n}{\Lambda}$。至关重要的是，纵向传播常数 $\beta_n$ 定义为：

瑞利-伍德异常（Rayleigh-Wood Anomalies） 定义为存在某个整数 $n$ 使得 $\beta_n = 0$ 的临界物理状态，此时对应的第 $n$ 阶衍射波以掠入射（grazing）形式沿周期方向无限传播，不再向纵向辐射。

1.2 技术难点：标准 PML-BIE 方法在 RW 异常处的失效机理

为了将无界单元边界 $\Gamma_2$ 和 $\Gamma_3$ 上的边界积分进行截断，本方法在纵向（$x_2$ 轴）引入复拉伸完全匹配层（PML）。复拉伸坐标定义为：

选用经典的分段幂函数型吸收函数 $\sigma(x_2)$（在物理介质层 $|x_2| \le H$ 处 $\sigma = 0$，在 PML 层 $H < |x_2| \le H+T$ 内按 $\sigma(x_2) = S (\frac{x_2-H}{T})^P$ 渐变增加）。在复拉伸坐标系下，自由空间格林函数形式变为 $\widetilde{\Phi}(k_1, x, y) = \frac{i}{4} H_0^{(1)}(k_1 \rho(\tilde{x}, \tilde{y}))$，其中拉伸距离 $\rho$ 定义为：

然而，当直接将边界积分系统截断在 $[-H-T, H+T]$ 的有限周期边界 $\Gamma_2^b$ 和 $\Gamma_3^b$ 上时，系统的收敛行为会遭遇严重的局部退化。论文在第三节给出了严格的数学证明：

对于 $n \in \mathcal{P}$（传播模式，$\beta_n > 0$）和 $n \in \mathcal{Q}$（渐逝模式，$\beta_n \in i\mathbb{R}^+$），PML 截断产生的尾部算子误差 $\psi$ 会随 PML 厚度 $T$ 的增加呈指数衰减：

致命的技术难点在于：

1. 当频率极其接近 RW 异常频点时，最小纵向波数 $\beta_n \to 0$，这导致传播模式的衰减速率中包含项 $e^{-\beta_n \frac{ST}{P+1}}$。当 $\beta_n \to 0$ 时，该指数衰减项退化，导致标准的 PML-BIE 方法在 RW 点邻域内的指数收敛阶急剧退化，最终导致数值精度崩塌（如图 3(a)-(c) 所示）。
2. 从边界积分算子的谱理论来看，原本用于约束物理散射场的投影关系： $$\frac{1}{\Lambda} \int_{-\Lambda/2}^{\Lambda/2} \left( \partial_{x_2} u^{\text{sct}}(x_1, \pm h) \mp i \beta_n u^{\text{sct}}(x_1, \pm h) \right) e^{-i \alpha_n x_1} dx_1 = 0 \quad (n\in\mathbb{Z})$$ 由于高阶不传播模式在有限 PML 截断下对应的尾部算子变得极端微小，导致用于求解非传播模式（尤其是接近临界截止的 evanescent modes）振幅系数的线性子方程组的条件数急剧恶化。在数值离散后，这种病态性会污染整个 GMRES 迭代解的精度，产生严重的数值噪声。

1.3 核心解决方案：带有限模式修正的改性 PML-BIE 模型

为了彻底解决上述 RW 异常附近的精度崩塌，本论文提出了有限模式修正 PML-BIE 格式。其核心数学思想是将发散/慢收敛的尾部算子拆分为一个高精度的有限阶模式修正算子。整个推导逻辑如下：

步骤 1：尾部积分的精细化拆分与半无界投影

定义一组临界模式指标集 $\mathcal{C}_\delta = \{ n \in \mathbb{Z} : 0 \le |\beta_n| \le \delta \}$，其中 $\delta > 0$ 是一个用户定义的频段控制阈值。无法被标准 PML 高效吸收的慢衰减尾部场 $\psi_{p,C}$ 被近似为有限个临界 Rayleigh 模式的线性叠加：

这里，$\Psi_{n,p}^{\pm}$ 为对应第 $n$ 阶模式的半无限尾部算子：

其中 $\widetilde{\mathbb{T}}^{\pm}$ 包含在半无限区间 $[H+T, \infty)$ 和 $(-\infty, -H-T]$ 上的积分。

步骤 2：引入互补积分实现闭合解析表达

由于直接在半无限区间上对高度振荡的格林函数进行数值积分是极其困难且不稳定的，论文引入了**互补积分（Complementary Integrals）**思想。我们将无界区间 $\mathbb{R}$ 上的全积分解 $\mathbf{T}_{p,q}$（解析已知）减去有限区间 $\Gamma^b$ 上的截断数值积分 $\mathbf{T}^b_{p,q}$，从而间接、高精度地定义半无限区间积分：

通过这种代数变换，半无限尾部算子可以精确地写成有限截断算子与无界拟周期解析积分的组合（式 4.3）：

其中 $\widehat{\Psi}_{n,p}^{\pm}$ 是沿整条无界真实周期边界 $\Gamma_2$ 和 $\Gamma_3$ 的积分。根据瑞利模式在周期管道内的格林表示定理，这些无界积分可以被赋予极度简洁的代数闭合形式（式 4.6）：

由此，原本难以处理的半无限发散尾部积分，全部被转化为物理空间中的闭合解析项与计算域内的有限区间积分。这一步是本方法取得高精度突破的关键。

步骤 3：建立修正边界积分方程组（Modified PML-BIE）

将上述算子代入原积分系统，最终导出了关于改性密度 $\phi^{\text{mod}}$ 和模式系数 $C_n^{\pm}$ 的耦合矩阵方程系统（式 4.15）：

其中 $\widetilde{\mathbb{M}}$ 为一个新定义的有限秩算子（Finite-Rank Operator）：

其中 $\mathbb{L}_n^{\pm}$ 是一组作用于改性密度 $\phi^{\text{mod}}$ 上的线性泛函（由式 4.14 给出），$\widetilde{\mathbf{\Phi}}_n^{\pm}$ 则是包含 Rayleigh 迹的矢量场。由于模式系数 $C_n^{\pm}$ 被自动消去并完全隐式吸收进入算子 $\widetilde{\mathbb{M}}$ 中，这彻底避免了求解病态未知数，极大地降低了方程组的条件数。

1.4 理论延拓：洛必达法则处理 $\beta_n \to 0$ 的极限情况

当系统无限逼近物理 RW 异常点时，$\beta_n \to 0$，导致上述有限秩算子 $\widetilde{\mathbb{M}}$ 中的系数分母 $2i\beta_n \to 0$。为了保证算法在 RW 奇异点处的频域连续性与鲁棒性，论文在第四节第 4.4 部分通过洛必达法则对该极限进行了严格的数学分析。当 $\beta_n \to 0$ 时：

其中导数项 $\partial_{\beta_n} \widetilde{\mathbf{\Phi}}_n^{\pm}$ 的解析形式为：

通过这一解析极限代入，算法在 RW 奇点处的计算精度得以完美保持，完全消除了除以零的数值溢出灾难。

1.5 消除“角点奇异性”的超胞拟周期延拓方案

在改性泛函 $\mathbb{L}_n^{\pm}$ 的数值积分计算中（见式 4.14），由于积分路径在边界 $\Gamma_2^b$ 和 $\Gamma_3^b$ 与人工截断水平面 $x_2 = \pm h$ 的交界处存在几何直角，这会在数值积分中引入“角点奇异性”。

为突破这一数值不连续性，论文第 4.3 节提出了一种基于 $3\Lambda$-周期超胞表示（$3\Lambda$-period supercell representation）的拟周期空间延拓方案。我们定义扩展障碍物边界族：

利用散射场满足的全局拟周期边界条件，将局部的单单元格林积分改写为在 $3\Lambda$ 超胞内三个相邻单元的联合积分（式 4.18）。通过将源点和汇点在拟周期意义下平移错开，成功地将计算点从直角交界处的奇异点转移开，极大地保证了多极离散下数值积分的精度与收敛阶。

2. 关键 Benchmark 体系，计算所得数据，性能数据

论文通过两个极具代表性的数值算例，对提出的改性 PML-BIE 方法进行了严格的精确度与收敛性验证。计算平台基于 Fortran 代码，采用 OpenMP 进行多线程并行加速。

2.1 算例一：单周期风筝形（Kite-shaped）介质障碍物阵列

物理参数设置

• 周期 $\Lambda = 2$。 障碍物几何边界 $\Gamma_1$ 的参数化方程为： $$r_1(t) = r (\cos t + 0.65 \cos 2t - 0.65, 1.5 \sin t), \quad t \in [0, 2\pi), \quad r=0.5$$
• 介质物性参数：外界 wavenumber $k_1$（扫描范围覆盖 RW 异常），内部 $k_2 = 20$，透射阻抗比 $\eta = 1.0$。
• PML 参数：半宽度 $H = 4.0$，截断高度 $h=1.0$，PML 幂次 $P=8$，吸收强度系数 $S=6$。
• 入射条件：TE 极化平面波，入射角 $\theta^{\text{inc}} = \pi / 4$。
• 临界 RW 异常波数： $$k^* = \frac{2\pi}{\Lambda (1 - \sin\theta^{\text{inc}})} = \frac{2\pi}{2 (1 - \sin(\pi/4))} \approx 10.7261$$

核心计算数据对比

论文对非改性 PML-BIE、改性 PML-BIE 以及校正窗口格林函数法（Corrected WGF）进行了详尽的三方对比。表 1 给出了在不同波数 $k_1$（偏离 RW、处于 RW 奇点、微扰 RW）下，各算法在截断宽度 $T$ 下的性能表现：

表 1：单周期风筝形障碍物下各方法的数值误差对比

数据深度解读与结论

1. 超越传统 WGF 的收敛速度：在相同的计算精度要求下，改性 PML-BIE 方法仅需 $T = 4\lambda$ 的 PML 厚度即可将自收敛误差控制在 $10^{-10}$ 级别，而校正 WGF 方法即使将窗口宽度扩大至 $T = 32\lambda$（计算自由度呈几何级数增加），其自收敛误差仍停留在 $10^{-3} \sim 10^{-5}$ 级别。这证明了本方法在极小截断尺寸下的卓越计算效率。
2. 完全的 RW 鲁棒性：图 3 和图 4 展示了能量平衡误差随 $T$ 增长的变化曲线。未改性的 PML-BIE 方法（图 3a-c）在 $k_1 = k^*$ 处误差曲线高度退化，无法收敛。而本文提出的改性 PML-BIE 方案（图 3d-f）的误差呈现出完美的线性指数下降趋势（Log 对数坐标下），且 $k_1 = 10.68, k^*, 10.76$ 三条曲线几乎重合。这无可辩驳地证明，本算法的收敛速率已经从物理 RW 异常的病态制约中完全解脱出来。

2.2 算例二：多形状分层周期障碍物阵列（复杂几何拓扑）

为了展示该算法处理复杂多介质界面、高非均匀拓扑结构的通用能力，算例二在一个周期单元内放置了五个形状各异的介质障碍物（依次为：椭圆、圆角三角形、风筝形、圆角星形、圆形），纵向半宽度设为 $H=6.0$，截断高度 $h = 5.4$。PML 参数与算例一保持一致。

非直边界截断测试

在实际复杂仿真中，人工边界 $\Gamma_2, \Gamma_3$ 常常需要弯曲以避开复杂的障碍物。本算例特别测试了非直（Non-straight）人工截断边界的适应性（如图 8 所示）。

表 2：多形状复杂障碍物阵列下（直边界 vs 非直边界）改性 PML-BIE 方法误差 ($T = 4\lambda$)

数据深度解读与结论

该组数据表明，当人工边界由直线变为适应性弯曲线时，算法的精度没有任何实质性退化。在极高频、极复杂的介质环境以及临界物理状态下，所有自收敛、拟周期和放射边界条件误差均稳定维持在 $10^{-10}$ 到 $10^{-13}$ 的极低区间。这说明引入的“角点奇异性消除技术”与“有限模式修正”能够完美兼容任意曲线路径的网格离散，具备极强的工业实际仿真落地价值。

3. 代码实现细节，复现指南，所用的软件包及开源 repo link

3.1 数值离散核心架构

为了复现本论文提出的频域一致鲁棒 PML-BIE 方法，数值离散应采用基于切比雪夫（Chebyshev）控制点的高阶 Nyström 离散技术：

1. 边界参数化与离散：
   • 障碍物边界 $\Gamma_1$ 采用均匀网格进行 Fourier-Nyström 离散（由于其周期性边界闭合）。
   • 无界截断边界 $\Gamma_2^b$ 和 $\Gamma_3^b$ 映射到有界区间 $[-H-T, H+T]$ 上，采用**分段 Chebyshev 节点（Panels with Chebyshev nodes）**进行离散，以便在接近截断边缘和直角转折点处实现局部高密度网格加密。
2. 格林函数奇异性消解（Singularity Subtraction）： 由于二维第一类 Hankel 函数 $H_0^{(1)}(z)$ 在 $z\to 0$ 处具有 $\ln|z|$ 的对数奇异性，在计算自相互作用矩阵元（Self-interaction matrix elements）时，必须采用经典的 Alpert 积分象限法（Alpert’s quadrature rule） 或 Kress 分裂法（Kress’s splitting method） 提取并单独高阶积分解析奇异项： $$H_0^{(1)}(z) = J_0(z) + i Y_0(z) \quad \Rightarrow \quad Y_0(z) \approx \frac{2}{\pi} J_0(z) \ln\left(\frac{z}{2}\right) + \text{smooth terms}$$

3.2 线性系统求解与预条件技术

1. 线性方程组组装： 在对 Dirichlet 迹和 Neumann 迹进行高阶 Nyström 离散后，边界积分方程（式 4.15）被转化为如下致密复数线性代数方程组： $$(\mathbf{E} + \mathbf{T}^b + \mathbf{M}) \vec{\phi}^{\text{mod}} = \vec{\phi}^{\text{inc}}$$ 其中 $\mathbf{E}$ 为对角块矩阵，$\mathbf{T}^b$ 为代表局部有限积分的稠密矩阵，$\mathbf{M}$ 为由于有限模式修正产生的低秩非局部校正矩阵（Low-rank non-local correction matrix）。由于 $\mathbf{M}$ 是通过外积 $\vec{\Phi}_n \otimes \vec{\mathbb{L}}_n$ 形式定义的，在实际存储和计算中应避免直接组装整个稠密 $\mathbf{M}$，而是采用低秩快速乘法加速。
2. 快速求解器： 使用**广义极小残量法（GMRES）**进行迭代求解。设定收敛残差标准（Residual tolerance）为 $\epsilon_r = 10^{-14}$。由于有限模式修正消除了 RW 带来的病态，GMRES 迭代次数在奇点附近能保持完全稳定，不需要额外的复杂几何预条件子。

3.3 算法复现流程指南

1. 第一阶段：网格生成与参数初始化： 输入物理波数 $k_1$。检查是否存在 $n$ 使得 $|k_1 - |\alpha_n|| \le \delta$。若有，将对应的模式 $n$ 加入临界修正集合 $\mathcal{C}_\delta$。设置吸收层厚度 $T = 4\lambda$，计算 PML 各点复伸缩坐标 $\tilde{x}_2$。
2. 第二阶段：矩阵元计算： 对 $y \in \Gamma_j^b$，$x \in \Gamma_i^b$，计算拉伸距离 $\rho(\tilde{x},\tilde{y})$，进而调用 Hankel 函数计算 $\mathbf{T}^b$ 矩阵块。利用极小邻域超胞映射公式（式 4.18）计算含有拟周期相位算子 $\zeta^m$ 的高阶无角点奇异性校正积分 $\mathbb{L}_n^{\pm}$。
3. 第三阶段：奇异值洛必达处理： 在复循环中判断若 $\beta_n = 0$，激活限幅解析项 $\partial_{\beta_n} \widetilde{\mathbf{\Phi}}_n^{\pm}$，组装低秩非局部校正矩阵 $\mathbf{M}$。
4. 第四阶段：GMRES 求解与散射场后处理： 求解得到 $\phi^{\text{mod}}$，带入解析外推公式（式 4.17）计算空间任意观测点处的改性总场与物理散射场。

3.4 社区开源工具推荐与类似 Repo

虽然本论文的官方 Fortran 并行程序目前尚未整体开源，但其数学底层和离散算法与国际著名的边界积分方程库有着极高的重合度。以下开源库是进行本算法复现的极佳底层工具选择：

• FMMLIB2D (Courant Institute, NYU):
  • 链接: https://github.com/flatironinstitute/fmmlib2d
  • 用途: 提供了极其高效的二维亥姆霍兹格林函数快速多极子算法（FMM）计算、Hankel 函数的高精度复评估程序，是构建稠密边界积分矩阵的核心计算引擎。
• BIE2D (Alex Barnett, Flatiron Institute):
  • 链接: https://github.com/ahb/bie2d
  • 用途: 一个功能全面的二维边界积分方程 MATLAB 平台，内置了高阶 Nyström 切比雪夫离散、Alpert/Kress 奇异性消解积分子系统，非常适合用来快速原型化开发本论文中的三胞超胞积分延拓算法。

4. 关键引用文献，以及你对这项工作局限性的评论

4.1 关键参考文献及其承接关系

本项研究建立在国际计算波物理学界过去二十年关于无界周期边界截断技术的几项里程碑工作之上：

• [31] T. Strauszer-Caussade, L. M. Faria, A. Fernandez-Lado, and C. Pérez-Arancibia, Stud. Appl. Math., 150 (2023), pp. 277–315.
  • 承接关系：本论文的核心创新灵感——有限模式修正算子（Finite-mode correction），直接继承自该文献中针对**窗口格林函数法（WGF）**提出的模式校正思路。本论文将其成功推广至复拉伸 PML 边界积分框架，由于 PML 自身自带指数级高精度收敛，相比 WGF 的超代数收敛（Super-algebraic convergence），本方法在更小的算子尺寸下实现了更高的数值精度。
• [24] W. Lu, Y. Y. Lu, and J. Qian, SIAM J. Appl. Math., 78 (2018), pp. 246–265.
  • 承接关系：该文献奠定了层状介质中 PML-BIE 方法的基本数学构造。本论文在此基础上，进一步解决了层状边界在拟周期对称性破坏（即由于周期结构引入的非绝对收敛晶格和）以及极临界 RW 异常频点处的退化机制。
• [6] O. P. Bruno and A. G. Fernandez-Lado, Proc. A., 473 (2017), p. 20160802.
  • 承接关系：该文献代表了采用拟周期格林函数（QGF）直接求解周期散射的国际领先水平。本论文提出的自由空间 PML 替代方案，彻底规避了 QGF 方法中极其高昂且极易不收敛的 Ewald 求和技术，提供了更干净、普适的格林求解框架。

4.2 本项工作局限性的客观批判

尽管本工作在二维周期阵列散射计算中取得了开创性的成功，但在将其向更广阔的工业级电磁/声学散射应用以及量子物理仿真拓展时，仍存在以下不可忽视的局限性：

1. 维度的瓶颈（难以向 3D 周期直接推广）： 本工作目前严格限定于二维（2D）空间。对于真正的三维双周期（Bi-periodic 3D）光栅或晶格结构，无界单元边界将从两条二维曲线 $\Gamma_2, \Gamma_3$ 演变为四个三维无限长曲面。在三维情况下，Maxwell 方程组具有强烈的矢量耦合特征，PML 的拉伸必须沿两个横向轴同时进行。对应的三维瑞利展展开和临界非传播模式修正矩阵的代数推导其复杂程度呈指数级增加，洛必达极限的解析推导难度极大。
2. 高频状态下的计算复杂度激增（缺乏快速求解支持）： 目前算法基于稠密 Nyström 矩阵求解。当结构的电尺寸极度增加（即进入高频波段，波数 $k_1 \Lambda \gg 1$），边界离散控制点的数量 $N$ 迅速膨胀，构造 $O(N^2)$ 的稠密矩阵和执行 $O(N^3)$ 的直接求解或 $O(N^2)$ 的 GMRES 矩阵向量乘积将面临物理内存和算力的巨大红线。虽然本方法相比 WGF 大幅减小了 PML 的截断长度，但若没有快速多极子（FMM）或 H-matrix 等快速代数加速算法的支持，仍难以应付极大规模高频周期的仿真。
3. 非均匀多层介质突变区域的处理局限： 当前的公式高度依赖于外界介质 $D_1$ 是完全均匀介质的物理假设。如果外界背景包含更加复杂的连续渐变折射率梯度（Graded-index media）或多层各向异性极化材料，自由空间格林函数将不再适用，这就瓦解了互补积分中利用闭合解析 Rayleigh 迹进行解析转换的理论大厦。

5. 其他必要的补充

5.1 与量子力学、量子化学中散射理论的跨界映射

虽然本论文直接针对的是经典电磁波与声波散射问题，但由于数学物理方程的同构性，本项工作提出的算法思想与量子化学、固体物理中周期晶格内的电子散射状态求解（Schrödinger 散射理论）有着极深的学术纽带。

1. 薛定谔方程与亥姆霍兹方程的等价性

在单电子近似下，一维/二维周期晶格中能量为 $E$ 的电子波函数 $\Psi(x)$ 满足定态薛定谔方程：

在化学物理中，由于原子势阱 $V(x)$ 的强周期性（如晶格常数为 $\Lambda$），这可以完全写为非均匀亥姆霍兹形式：

其中电子在无界真空层区域的散射满足的 Wigner 阈值定律（Wigner threshold law） 或是通道开启（channel opening）状态，其物理数学本质正对应于本论文中所研究的 Rayleigh-Wood 异常频点。当入射电子能量 $E$ 刚好跨越某个束缚态或散射通道的开启能量时，纵向动量常数 $\beta_n = \sqrt{2m(E-E_{th})/\hbar^2} \to 0$。

2. 量子输运与格林函数计算的性能提升

在紧束缚模型（Tight-binding）和基于密度泛函理论的量子输运（如非平衡格林函数法 NEGF）中，需要频繁计算无界电极（Leads）与中心散射区交界处的自能矩阵（Self-energy matrices）。自能算子的本质正是边界上的 Dirichlet-to-Neumann (DtN) 映射。传统的 DtN 映射在能带边（Band edge）附近会遭遇剧烈的奇异性发散，导致非平衡电子输运计算极难收敛。如果能将本论文提出的“基于有限模式修正的复拉伸 PML-BIE 方法”引入到量子化学的第一性原理电子输运计算中，将有望彻底打通能带边奇异点处电子透射谱的高精度仿真瓶颈。

5.2 总结与展望

本研究提出的带有限模式修正的 PML-BIE 方法，以极其精妙、严谨的数学代数处理，扫清了经典波周期散射模拟中困扰学术界数十年的瑞利-伍德异常阴霾。论文通过严格的数值理论证明与极具挑战性的多物理形状基准测试，表明该方法具有极高的精度、极致的截断效率、以及面对复杂边界几何时卓越的鲁棒性。随着研究人员进一步将其推向三维矢量 Maxwell 体系并与快速算法相结合，这一频域鲁棒的截断边界积分理论必将在光子晶体、超材料微纳光子学计算设计以及量子化学输运计算中大放异彩。