基于FDTD的二维TM波金属方柱散射仿真

资 源 简 介

本项目旨在利用时域有限差分法（FDTD）在二维笛卡尔坐标系下，构建一个完整的电磁波散射仿真环境。主要针对横磁波（TMz模式，包含Ez、Hx、Hy分量）与理想导体（PEC）金属方柱的相互作用进行数值模拟。在该系统中，计算区域被划分为基于Yee元胞的网格，金属方柱通过将内部及边界上的电场分量强制设为零来实现PEC边界条件。为了模拟无限延伸的自由空间并消除计算窗口边缘的非物理反射，项目采用了单轴各向异性完全匹配层（UPML）作为吸收边界条件，推导并实现了UPML区域内基于拉伸坐标系的场刷新公式，确保了宽频带内的高效吸收。激励源选用高斯脉冲或调制高斯脉冲，通过硬源或软源方式加入计算区域，以模拟宽频带入射波的瞬态响应。项目不仅能够计算电磁波的传播、反射和绕射过程，还能通过总场/散射场（TF/SF）边界设置将入射波与散射波分离（可选功能），重点展示金属方柱产生的散射场分布。该代码实现了电磁场分量的时间步进迭代，并提供了动态的可视化模块，用于观测电磁波遇到障碍物后的波前变化及能量分布。

详 情 说 明

基于FDTD的二维TM波金属方柱散射仿真系统

项目简介

本项目构建了一个基于时域有限差分法（FDTD）的二维电磁波散射仿真环境。仿真主要针对横磁波（TMz模式）在包含理想电导体（PEC）金属方柱的自由空间中的传播与散射过程进行模拟。系统采用Matlab语言编写，实现了从网格剖分、激励源设置、边界条件处理到场量迭代更新及动态可视化的完整流程。

通过本项目，用户可以观察电磁波在遇到金属障碍物时的反射、绕射现象，并验证单轴各向异性完全匹配层（UPML/CPML）在截断计算区域时的吸收效果。

功能特性

• 物理模型：模拟二维TMz极化波（包含电场 $E_z$ 及磁场 $H_x, H_y$ 分量）。
• 网格算法：基于Yee元胞的空间交错网格技术。
• 边界条件：实现了基于卷积完美匹配层（CPML/UPML）的吸收边界，有效消除截断边界处的非物理反射。

* 采用多项式分布的电导率 $sigma$ 和 衰减因子 $kappa$。 * 基于复频移（CFS）参数 $alpha$ 增强低频与倏逝波吸收。

• 散射体模型：内置PEC金属方柱模型，支持参数化配置尺寸与位置。
• 激励源：高斯调制正弦脉冲（Gaussian Modulated Sinusoidal Pulse），覆盖宽频带特性。
• 实时可视化：仿真过程中动态展示 $E_z$ 场分布色图，并实时绘制监测点的时间域波形。

系统要求

• MATLAB R2016b 或更高版本（推荐）。
• 无需额外工具箱，基于MATLAB基础矩阵运算实现。

使用方法

1. 确保运行环境已配置好MATLAB。
2. 打开主脚本文件（通常为 main.m）。
3. 根据需求修改脚本头部的“仿真参数设置”区域（如网格尺寸、频率、方柱大小等）。
4. 运行脚本，程序将自动进行FDTD迭代计算。
5. 在弹出的图形窗口中观察电磁波传播动画及探针信号。

代码分析与实现细节

本项目代码逻辑清晰，主要分为以下几个核心模块：

1. 物理常数与空间网格初始化

代码首先定义了真空光速、磁导率、介电常数及自由空间波阻抗等基础物理量。计算区域被划分为 $200 times 200$ 的网格，空间步长设定为 $2text{mm}$。

2. 时间步长与稳定性控制

为了保证数值计算的稳定性，时间步长 $Delta t$ 根据Courant-Friedrichs-Lewy (CFL) 条件自动计算，并设定为理论最大值的 $0.98$ 倍，以确保在二维笛卡尔坐标系下满足数值稳定性要求。

3. 宽带激励源配置

激励源参数设定为中心频率 $2.0 text{GHz}$ 的高斯调制脉冲。代码计算了带宽、脉宽参数 $tau$ 以及时间延迟 $t_0$，确保脉冲在仿真开始时幅值平滑上升，避免高频噪声引入。源的位置被硬编码在网格的特定坐标处（$X=1/4, Y=1/2$ 处）。

4. PEC金属方柱建模

代码采用了“掩模法”来处理PEC散射体。

• 定义了方柱的物理边长（$0.08text{m}$）及中心位置。
• 将物理尺寸映射为网格索引范围。
• 创建了一个布尔类型的掩模矩阵（pec_mask），在金属方柱存在的区域标记为 true。该掩模将在后续的时间步进循环中用于强制将金属内部的电场置零。

5. UPML/CPML 吸收边界实现

这是代码中最复杂的部分，采用了标准CPML（卷积PML）策略来实现UPML边界，具体实现逻辑如下：

• 参数分布：PML层数设为10层。在边界区域内，电导率 $sigma$ 和拉伸坐标系数 $kappa$ 按照3阶多项式从计算区域内部向边界外部逐渐增加，以减少阻抗失配。
• 交错网格处理：代码分别计算了整数网格点（对应 $E$ 场位置）和半整数网格点（对应 $H$ 场位置）的PML系数，确保场分量在空间上的精确对准。
• CPML系数计算：不同于传统的裂场（Split-field）PML，代码预计算了递推卷积所需的系数 $b$ 和 $c$ 以及空间导数系数。

* 计算了 $X$ 和 $Y$ 方向上的 be, ce 系数。 * 初始化了辅助累加变量 psi 矩阵（psi_Ez_x, psi_Ez_y 等），用于存储历史场信息，通过标准的CPML更新方程实现宽频带吸收。 * 预计算了 $1/kappa$ 项，用于标准麦克斯韦方程的更新修正。

6. 可视化模块初始化

在进入时间循环前，代码初始化了绘图环境：

• 创建双子图窗口：左侧显示 $E_z$ 场的实时分布，右侧显示探针点的时域波形。
• 利用 imagesc 绘制场强，并使用 jet 色图。
• 利用 rectangle 函数在场分布图上叠加绘制出白色边框的金属方柱，用于直观展示散射体的物理位置。

7. 数据结构说明

• 场分量：Ez, Hx, Hy 矩阵存储当前时间步的电磁场值。
• 辅助变量：虽然定义了 Dz, Bx, By 等变量，但依据代码后半段的系数计算策略，实际计算主要依赖于CPML的 psi 数组和直接的 $E-H$ 更新。
• 监测探针：Ez_probe 数组用于记录源后方某一点的电场随时间变化的数据，用于分析反射波特性。