6.2 矩量法与快速算法

6.2 矩量法与快速算法

FDTD 给整个空间都立了户口，但辐射与散射问题里，真正"活着"的只有导体与介质表面——空气是无辜的旁观者。矩量法（MoM）正是这条思路：用积分方程把未知量压到目标表面，只用给表面电流立户，空气不占编制。它天生适合开放域（天线、散射、互耦），代价是方程组变成稠密矩阵。学完本节，你将掌握"基函数展开—检验函数加权—解稠密方程组"的三步流程，明白内存账单为何随电尺寸三次方增长，以及快速多极子凭什么把这笔账砍下几个数量级。下一节的有限元法则接管另一头：复杂材料与封闭腔体。

一、积分方程：只给表面立户

思路从 3.1 的滞后位直接生长：任意导体表面的感应电流 J 会在空间产生散射场，边界条件要求"表面上总切向电场为零"（理想导体），把散射场用 J 的积分写出、代回边界条件，得到电场积分方程（EFIE）——未知函数 J 藏在积分号内。磁場积分方程（MFIE）与混合形式处理不同几何与谐振场景，原理同源。积分方程的妙处：未知量从三维体积（FDTD 的每个网格）压缩到二维表面，且辐射条件自动满足——无穷远的边界不需要人工吸收，这是它相对微分方法的先天优势。

离散化三步走：把表面电流用基函数展开（常用三角形上的 RWG 基函数，恰好保证电流连续与电荷守恒）；用检验函数对积分方程做加权（伽辽金法取基函数自身）；得到稠密矩阵方程 Z·I = V，解之即得表面电流分布，远场再由辐射积分一次算出。

# 细偶极子天线的迷你矩量法：分段正弦基 + 电场积分方程
import numpy as np

f = 300e6; lam = 1.0; k = 2*np.pi/lam
L_half, N = 0.25, 21                 # 半波偶极子剖成 N 段
z = np.linspace(-L_half, L_half, N)
dz = z[1]-z[0]

def green(r1, r2):
    r = abs(r1-r2)
    return np.exp(-1j*k*r)/(4*np.pi*r)   # 自由空间格林函数（标量示意）

Z = np.zeros((N, N), dtype=complex)
for m in range(N):
    for n in range(N):
        r = abs(z[m]-z[n])
        if r < dz/2:                     # 自奇异项：用近似距离避免发散
            r = dz/4
        Z[m, n] = green(z[m], z[n])*1j*376.73*(1 + 1j/(k*r)) * dz

V = np.zeros(N, dtype=complex); V[N//2] = 1.0   # 中心德尔塔馈电
I = np.linalg.solve(Z, V)
print(f"输入阻抗估计 {V[N//2]/I[N//2]:.1f} 欧（参考半波偶极子约 73 加感抗）")
print("电流分布（归一化幅度）：")
for zi, ii in zip(z[::4], np.abs(I)[::4]):
    bar = "#" * int(abs(ii)/np.abs(I).max()*30)
    print(f" z={zi:+.2f} m |{bar}")
​

这段几十行的玩具内核已经能给出定性正确的电流驻波分布与阻抗量级——生产级代码会换用更严谨的基函数与奇异性积分处理，但骨架不变。MoM 的工程感就在这里：矩阵就是天线，电流就是答案。

二、稠密矩阵的账单与快速算法的谈判

MoM 的账单由稠密矩阵决定：N 个未知量，内存 O(N²)，直接求解 O(N³)（迭代法 O(N²) 每次迭代）。数量级感受：一根飞机天线几千未知量轻松，整架飞机的散射上百万未知量——平方下来内存以 TB 计，经典 MoM 在电大目标前举手投降。

快速多极子算法（FMM）重新谈判这笔账。观察：格林函数是平滑的远距函数，远距离组之间的相互作用不必两两精确计算——把表面电流分组，组内直接算、近邻直接算、远组间用"多极展开"聚合—转移—散开（科学家们戏称 ATM 机流程），内存降到 O(N·logN) 量级，百万未知量从不可能变成一晚上的事。积分方程族系由此翻新：多层快速多极子撑起了大型舰船与整机 RCS 计算，预修正快速傅里叶变换则擅长平面周期结构。代价是并行实现复杂、近远场切换的精度管理——快，是要用工程复杂度换的。

三、预条件与并行：让稠密方程组肯开口

有了快速多极子，每次矩阵向量乘从平方代价降到近线性，但迭代求解器本身还有一道坎：条件数。EFIE 在谐振频率附近条件数急剧恶化，GMRES 的残差曲线长时间趴着不动，几千次迭代都不收敛——值班室里最常见的"跑一晚上没结果"多半是它。两条常用对策：其一是换方程，用磁场积分方程与电场积分方程线性组合成混合场积分方程（CFIE），把内谐振的病态从源头压下去（代价是只适用于闭合导体）；其二是上预条件，用近相互作用矩阵（本来就是稀疏的）构造不完全 LU 或块对角预条件子，把条件数拉回迭代求解器舒服的区间，典型效果是迭代次数从上万次压到几百次。工程纪律：盯残差曲线的形状比盯绝对耗时更有信息量——平躺的残差预示病态，锯齿状残差常指向网格质量，平滑下降但收敛慢则该加预条件了。

未知量的账可以先于建模心算出来，避免"建完模才发现内存不够"的值班事故：

# 电大平台的矩量法未知量与内存账本（数量级估算）
import math

surface_area = 300.0        # 平台等效表面积，平方米
edge_len = 1.0 / 10.0       # 每波长十条边的经典剖分密度
tri_area = math.sqrt(3) / 4 * edge_len**2   # 等边三角形面积
for freq_ghz in (0.3, 1.0, 3.0):
    lam = 0.3 / freq_ghz
    n_tri = surface_area / (tri_area * lam**2)
    n_unknown = int(1.5 * n_tri)            # RWG 基函数数约为三角形数的 1.5 倍
    dense_gb = n_unknown**2 * 16 / 2**30    # 复数双精度稠密矩阵
    mlfma_gb = n_unknown * 17 * 16 / 2**30  # 多极展开近似 O(N log N)
    print(f"{freq_ghz:>4} GHz | 未知量 {n_unknown:>9,} | "
          f"稠密 {dense_gb:>8.1f} GB | 多极 {mlfma_gb:5.1f} GB")
​

数量级结论一目了然：三百兆赫兹时经典矩量法还能在工作站上喘气，到了三吉赫兹稠密矩阵以十太字节计，只有快速多极子加上并行分布内存才接得住这单出勤。这张账本也解释了为什么舰船与整机散射项目历来是计算电磁学的"重型装备竞赛"——谁先把百万未知量级的稠密方程解得又快又准，谁就能在隐身评估与平台天线布局上多出勤几轮。

四、适用边界与工程口感

MoM 的主场：开放域问题（自由空间中的天线与散射）、细线与金属结构（表面立户天然契合）、需要高精度阻抗与耦合参数的场合（天线阵互耦、馈电网络）。它的弱项同样明确：非均匀介质体（要用体积分方程或与 FEM 混合）、谐振结构的内场、以及材料极度复杂的场景——这些正是下一节 FEM 的领地。工程口感上，MoM 给出的阻抗与方向图精度高、参数化扫描快（矩阵重用方便），天线设计师的日常工具链几乎都长在它上面。

一条实用忠告：矩量法的精度陷阱常在"几何离散"而非算法本身——边缘、馈电点、缝隙处的网格剖分质量直接决定结果可信度，收敛性验证（加密一档重跑）是不可省的值班纪律。

本节要点回顾

• 表面立户：积分方程把未知量压到目标表面，辐射条件自动满足，开放域问题的天选方法。
• 三步流程：基函数展开、检验加权、稠密矩阵求解；RWG 基函数兼顾电流连续与电荷守恒。
• 账单结构：内存平方、求解立方，电大目标在经典 MoM 面前撞墙。
• 快速多极子：分组与多极展开把内存砍到 N 乘对数 N，整机级散射因此可行。
• 主场与让位：细线金属开放域是主场，复杂介质与封闭腔体让位给有限元——互补而非替代。

下一节请出变分原理：有限元法如何用能量泛函与四面体单元，把任意材料与任意几何都收进同一套离散框架。