二维光子晶体能带结构计算平台

资 源 简 介

该项目旨在利用平面波展开法（Plane Wave Expansion Method, PWE）实现对二维光子晶体光子能带结构的精确计算与可视化。主要功能模块包括几何建模、介电常数空间分布构建、傅里叶空间本征值问题求解以及能带图绘制。程序支持正方形晶格、三角形晶格等典型周期性结构，能够模拟电介质柱在背景介质中或空气孔在电介质基底中分布的物理模型。通过在倒空间中对麦克斯韦波动方程进行离散化，将波动方程转化为大型矩阵的特征值求解问题。系统支持TE（横电波）和TM（横磁波）两种极化模式的独立计算，并允许用户自定义

详 情 说 明

二维光子晶体能带结构计算平台

项目介绍

本项目是一款基于平面波展开法（Plane Wave Expansion Method, PWE）的数值模拟计算工具，专门用于分析二维光子晶体的能带结构特性。通过在倒空间中离散化麦克斯韦方程组，程序能够精确计算光子在本征模态下的色散关系，并识别光子带隙（Photonic Band Gap）。该工具为研究光子晶体的物理机制、设计光子器件以及探索新型光子拓扑态提供了基础的计算框架。

功能特性

1. 多种晶格支持：预设支持正方形晶格（Square Lattice）和三角形晶格（Triangular Lattice）两种主流的周期性结构。
2. 极化模式选择：支持横电波（TE）和横磁波（TM）两种独立极化模式的能带计算。
3. 计算参数自定义：用户可自由调节晶格常数、介质柱半径、背景及介质柱的介电常数。
4. 精度可控：通过设置倒格矢展开阶数控制平面波数量，在计算效率与精度之间取得平衡。
5. 高对称路径自动采样：自动生成第一布里渊区的不可约边界路径（如正方形的 Γ-X-M-Γ 或三角形的 Γ-M-K-Γ）。
6. 可视化与带隙分析：提供晶格结构示意图和能带图，并能自动识别带隙位置、宽度及对应能带。

实现逻辑

计算过程严格遵循以下物理与数学逻辑：

1. 空间初始化：定义实空间基矢与倒空间基矢，确定高对称点坐标。
2. 倒格矢集合构建：基于指定的展开项数，生成覆盖倒空间的双指数格点矩阵。
3. 傅里叶系数转换：计算介电常数倒数的傅里叶变换系数。对于圆柱体模型，程序利用解析公式（包含一阶贝塞尔函数）计算傅里叶分量，并构建 Toeplitz 形式的算子矩阵。
4. 波矢迭代求解：

* 沿着布里渊区路径遍历波矢 $k$。 * 针对每一个 $k$ 点，根据极化模式构建对应的特征值方程矩阵。 * TM 模式下利用电场标量方程，TE 模式下利用磁场矢量点乘特性构建矩阵。

1. 特征值提取：解算矩阵得到本征值，将其转换为归一化频率 $omega a / 2pi c$。
2. 结果输出：对能带进行排序并绘制曲线，通过搜索相邻能带间的最小/最大值差值来检测是否存在全局带隙。

关键算法与实现细节

1. 傅里叶空间变换

程序采用倒空间解析展开法。针对二维圆柱形介质柱，其介电常数的傅里叶分量通过 $2f J_1(Gr)/(Gr)$ 这一公式计算（其中 $f$ 为填充率，$J_1$ 为一阶贝塞尔函数），这种方法避免了复杂的空间离散化采样，提高了收敛速度。

2. 极化模式矩阵构建

• TM 模式：由于电场平行于柱体轴向，波动方程退化为标量形式，矩阵元素由 $eta(G-G') cdot |k+G'|^2$ 组成。
• TE 模式：磁场平行于轴向，电场在平面内。矩阵构建涉及 $(k+G)$ 与 $(k+G')$ 的点积运算，体现了矢量波的特性。

3. 本征值处理

程序默认提取前 10 条低频能带，利用 MATLAB 内置的高效特征值求解器（eig）对大型厄米矩阵（或非对称矩阵的特定形式）进行处理，并仅保留物理意义明确的正实数频率解。

4. 自动化带隙识别

在能带计算完成后，程序会自动扫描每一对相邻能带。如果第 $n$ 条能带的最大频率小于第 $n+1$ 条能带的最小频率，则判定为存在带隙，并在图表中以红色透明区域标出。

使用方法

1. 打开 MATLAB 环境。
2. 在主程序代码的“参数设置区”修改实验参数（如 lattice_type、polarization、r 等）。
3. 直接运行程序。
4. 控制台将输出带隙的具体数值范围，同时弹出包含晶格示意图和能带结构的对比图表。

系统要求

• 软件环境：MATLAB R2016b 及以上版本（需具备基础数值计算工具箱）。
• 硬件建议：建议 8GB 以上内存，由于平面波展开法涉及大型矩阵运算，较高的倒格矢展开数（N_order > 10）将显著增加计算开销。