伪代码示意

利用CST对三维超材料的能带计算程序

刚接触超材料能带计算那会儿，被三维结构搞得头皮发麻。传统平面结构的Floquet模式分析套路在立方体、金字塔这些立体结构面前直接失效。摸索了两个月，终于搞定了CST里三维能带计算的完整流程，这里把踩过的坑和关键代码整理成实战指南。

打开CST先别急着建模，先调后台参数。在Global Properties里把Solver选成"Eigenmode"，这个模式对周期性结构更友好。三维单元的边界条件最容易翻车，记得在Boundary Conditions里把X/Y/Z方向都设成"Periodic"，相位滞后参数直接用变量kx,ky,kz代替具体数值，后面跑参数扫描要用。

建模阶段有个偷懒技巧：用VBA脚本批量生成蜂窝状立体晶格。这段循环代码能自动排列六边形单元：

For i = 0 To 10 For j = 0 To 10 Component.New "UnitCell_" & i & "_" & j ' 生成六边形基体 Block.Create "Polygon", Array(0,0.5,0), Array(1,0,0), Array(0,1,0) ' 生成内部谐振结构 Cylinder.Create "SplitRing", Array(0.3,0.3,0), 0.2, 0.05 Component.Transform "Translate", Array(i*2.5, j*2.2, 0) Next Next

重点在Z轴方向的层叠处理，用Component.Stack命令实现三维堆叠时，相位匹配参数要跟着调整，否则算出的色散曲线会有断点。

能带扫描的核心是参数化k矢量。在Solver Parameters里创建三个参数变量：

StoreParameter("kx", 0) StoreParameter("ky", 0) StoreParameter("kz", 0)

然后写个三重循环遍历布里渊区路径。比如计算Γ-X-M路径时：

for kx in np.linspace(0, π/a, 20): for ky in np.linspace(0, 0, 20): for kz in np.linspace(0, 0, 20): update_parameters(kx, ky, kz) run_solver() extract_eigenfrequencies()

实际在CST里要用VBA实现，注意每次跑完仿真要Reset才能更新k值。建议用Batch Solver模式并行计算，比单线程快5倍以上。

后处理阶段最容易忽略模式追踪。用这个Python脚本匹配不同k点的模式连续性：

def mode_matching(prev_freq, current_freqs): # 寻找最接近的上一个模式频率 indices = np.argsort(np.abs(current_freqs - prev_freq)) return indices[0]

处理完的数据用Mayavi做三维等频面可视化，能清晰看到Dirac锥的位置。遇到带隙闭合的情况，重点检查结构对称性是否被破坏，或者晶格常数设置是否导致模式简并。

最后提醒几个实战经验：1）网格划分用"Hexahedral"类型更稳定 2）设置0.1%的损耗避免奇异解 3）先跑低精度测试找到关键k点区域。曾经有个案例，完整计算需要32小时，用这个技巧优化后8小时就出结果了。