COMSOL光子晶体六边形晶格简约第一布里渊区能带

今天在实验室折腾COMSOL的时候，突然被六边形光子晶体的能带结构卡住了脖子。每次看到布里渊区那几个对称点就头疼，特别是要把三维的布里渊区概念硬塞进二维六边形晶格的时候，简直像在玩拓扑版的俄罗斯方块。

先来点实在的，打开COMSOL新建模型，用几何序列生成六边形晶格最省事。这里有个小技巧——用参数化坐标生成六边形原胞：

a = 1e-6; % 晶格常数
r = 0.3*a; % 介质柱半径
theta = 0:60:300;
x = a*cosd(theta);
y = a*sind(theta);

这串代码画出的六边形边界线框，配合环形阵列复制就能生成完整晶格。不过要注意的是，周期性边界条件的设置必须严格对应六边形的基矢方向。有次手滑把波矢扫描方向搞反了，结果能带图直接变成抽象画，被老板当反面教材念叨了半个月。

设置材料属性时，介质柱的介电常数最好用分段函数定义：

material = (sqrt((x-x0)^2 + (y-y0)^2) < r) ? 12.0 : 1.0

这个三元运算符比几何选择更灵活，特别是处理不规则排布时。不过要小心网格划分的精细度，曾经有个案例是介质柱边缘出现锯齿导致带隙计算误差超过10%。

COMSOL光子晶体六边形晶格简约第一布里渊区能带

最刺激的还是布里渊区路径的设置。六边形的简约布里渊区边界由Γ-M-K-Γ这几个特征点构成，但COMSOL默认的直角坐标系需要换算成六边形倒格矢坐标系：

k_path = [
    0, 0;  % Γ点
    2*pi/(sqrt(3)*a), 0;  % M点
    pi/(sqrt(3)*a), pi/a;  % K点
    0, 0   % 回到Γ点
];

这个路径参数化直接决定能带计算的准确性。有次偷懒用线性插值代替实际倒格矢路径，结果在K点附近出现诡异的能带交叉，后来发现是路径取样不够导致的数值误差。

求解器参数设置有个隐藏技巧：本征频率搜索范围要跟着晶格常数动态调整。比如当a=1μm时，建议设置：

eigenfrequency = linspace(0.3e15, 0.7e15, 50)

这个量级刚好覆盖可见光波段。记得打开场监视器捕捉特定k点的电磁场分布，有次意外发现某个模态的电场集中在空气孔区域，直接推翻了之前关于带隙成因的假设。

后处理阶段建议用MATLAB直接解析COMSOL导出的.dat文件：

data = load('band_structure.dat');
scatter3(data(:,1), data(:,2), data(:,3), 20, data(:,4), 'filled');
colorbar;
view(0,90) % 俯视投影到二维布里渊区

这种三维散点图投影法能直观显示能带在动量空间的分布规律。有个冷知识：六边形晶体的狄拉克锥会出现在K点附近，当看到两个能带在该区域交叉成锥形结构时，八成是发现了类似石墨烯的光子晶体特性。

折腾完这些，终于看到漂亮的能带图跃然屏上。那些Γ-M-K路径上的起伏波纹，每个凹陷都可能藏着光子禁带的秘密。不过说实在的，每次看到计算完成的绿色进度条，都比看到实验结果还激动——这可能就是仿真工程师的奇怪乐趣吧。