COMSOL光学仿真：波动光学的有限元实现与散射分析

COMSOL在光学领域的核心能力来自波动光学模块和射线光学模块两个互补的物理场接口。前者用有限元方法直接求解麦克斯韦方程组，适用于亚波长结构和近场分析；后者用几何光学近似追踪光线轨迹，适用于宏观光学系统。两者的选择依据是结构特征尺寸与光波长的比值。

波动光学模块的物理基础

COMSOL波动光学模块求解频域麦克斯韦方程组：

∇×(μ_r⁻¹·∇×E) – k₀²(ε_r – jσ/(ωε₀))E = 0

其中ε_r是相对介电常数，μ_r是相对磁导率，k₀=ω/c是真空波数。这是一个矢量亥姆霍兹方程，COMSOL用矢量基函数（Nédélec单元）离散化。

频域分析假设电磁场是时谐的（e^{jωt}），适用于单色光照明。对于宽带或脉冲光源，需要做波长扫描——COMSOL支持在单个研究步骤中扫描多个波长，每个波长独立求解。

端口与边界条件

端口边界（Port）：波动光学中最核心的边界条件。Port定义入射波的横场分布和传播方向，支持矩形波导模式、圆波导模式、平面波和自定义模式。端口同时也是出射边界——出射波被吸收，反射系数直接从端口参数提取。

完美匹配层（PML）：截断无限大计算域的人工吸收层。PML的原理是在边界区域设置复数坐标变换，使进入PML的电磁波被指数衰减吸收，无反射。PML的层数通常5-8层，每层厚度λ/10到λ/4。PML太薄会导致吸收不完全（残余反射>1%），太厚浪费计算资源。

散射边界条件（Scattering Boundary Condition）：近似吸收边界，适用于法向入射的平面波。比PML简单但精度低，只适合远场散射方向已知的情况。

周期性边界条件（Floquet Periodicity）：用于模拟无限周期结构（光栅、超表面）。源面和目标面的场用Floquet定理关联：E(r+R) = E(r)*exp(jk·R)，其中R是晶格矢量，k是布洛赫波矢。

散射场公式

对于散射问题（微粒、纳米结构的光散射），COMSOL采用散射场公式：

E_total = E_background + E_scatter

求解器计算的是散射场E_scatter，背景场E_background由解析公式给出（平面波、高斯光束等）。这种方法的优势是入射场不需要在计算域中离散化，避免了网格对入射波的数值色散误差。

散射截面和吸收截面的计算：

σ_scat = ∫∫(E_scatter × H_scatter*)·dS / I_inc σ_abs = ∫∫(E_total × H_total*)·dS / I_inc – σ_scat

消光截面σ_ext = σ_scat + σ_abs，与Mie理论解析解的对比是验证模型精度的标准方法。

实际案例：金纳米球的光散射

以直径100 nm的金纳米球在可见光波段的散射为例：

材料属性

金的介电常数用Drude-Lorentz模型或Johnson & Christy实验数据。COMSOL内置了Au的材料数据（插值函数）：

• 波长400 nm：ε_r = -1.46 + 2.41j
• 波长500 nm：ε_r = -2.54 + 3.56j
• 波长600 nm：ε_r = -10.7 + 1.2j
• 波长700 nm：ε_r = -16.2 + 1.1j

介电常数实部为负（金属特性），虚部正（吸收），且都是波长的函数。

模型设置

1. 几何：直径100 nm的球，计算域取直径2 μm的球壳
2. PML：外层5层，每层λ/8
3. 背景场：沿z方向传播的平面波，x偏振
4. 网格：纳米球内四面体网格，最大尺寸λ/(10*n)，约30 nm；PML区域扫掠网格
5. 波长扫描：400-700 nm，步长10 nm

结果分析

消光截面在530 nm附近出现峰值（局域表面等离子共振，LSPR），σ_ext≈2.5×10⁻¹⁴ m²，与Mie理论计算值偏差<2%。

近场增强|E/E_0|²在球表面达到约15倍，共振波长处增强最大。这个近场增强是表面增强拉曼散射（SERS）的物理基础。

射线光学模块

当结构尺寸>>波长时（如透镜系统、望远镜），波动光学的有限元方法计算量不可接受。射线光学用几何近似追踪光线轨迹：

• 光线在均匀介质中沿直线传播
• 在折射率突变界面遵循Snell定律
• 在连续变化介质中遵循费马原理（光程极值）

COMSOL射线光学的应用场景：

光线追踪的蒙特卡洛扩展

对于散射介质（如大气、浑浊流体），COMSOL支持在光线与散射粒子碰撞时随机改变方向（Henvey-Greenstein相函数），用蒙特卡洛方法统计光强分布。需要追踪大量光线（通常10⁴-10⁶条）才能获得统计稳定的结果。

网格要求与精度控制

波动光学的网格要求比结构力学更严格：

最小网格尺寸：电磁波在介质中的波长λ_n = λ₀/n。每个波长内至少需要5-10个单元才能准确解析波形。因此最大网格尺寸h_max ≤ λ_n/10。

金属区域：金属中电磁波快速衰减（趋肤效应），趋肤深度δ = λ/(2π*Im(n))。金在可见光波段δ≈25-40 nm。金属区域的网格需要在趋肤深度方向至少3-5层，每层厚度δ/3。

界面处网格匹配：介电常数突变的界面（金属-介质、高折射率-低折射率）两侧网格必须匹配，否则界面反射系数会引入数值误差。COMSOL默认在界面两侧生成一致网格。

常见问题诊断

结果依赖网格尺寸：典型的网格不收敛问题。逐步加密网格（h_max从λ/5到λ/20），检查关心的物理量（散射截面、共振频率）是否收敛。如果Δσ_scat/σ_scat < 1%在h_max=λ/10时满足，网格足够。

PML残余反射：如果PML吸收不完全，会在计算域内形成驻波，导致结果异常振荡。解决方法：增加PML层数、增大PML厚度、检查PML区域网格（应沿法向扫掠）。

端口模式不匹配：如果Port定义的模式与实际波导模式不符，入射波会在端口处反射，污染计算结果。用COMSOL的”Port Mode Analysis”研究步骤先求解波导模式，再将模式索引传递给Port边界。