CREST频率计算与振动模式

想确认小分子结构是否位于局部极小点，可以在优化后计算 Hessian 和振动频率。本例用两个孤立分子 H₂O、H₂S，介绍 CREST 3.0.2 的原生 ohess 工作流；它会先优化，再计算数值 Hessian。

1. 下载与运行#

下载完整输入、参考数据和绘图包，解压后进入 crest-harmonic，用 sbatch job.slurm 提交。CPU 分区按 56 核整节点分配和计费，包内最多同时运行四个单线程小任务。在已分配的 CPU 作业内也可运行：

bash run.sh

软件为 CREST 3.0.2，内置 tblite/GFN2；电子温度 300 K、accuracy=0.0001、extreme 优化。Python 分析与绘图使用 sisl/0.16.4 模块环境。

2. 准备结构与完整设置#

水分子的起始输入 inputs/water/structure.xyz：

3
Own starting structure, Angstrom
O 0 0 0
H 0.9572 0 0
H -0.2399872084090341 0.9266272064859951 0

硫化氢输入 inputs/hydrogen-sulfide/structure.xyz：

3
Own starting structure, Angstrom
S 0 0 0
H 1.336 0 0
H -0.04895599483195902 1.3351027340882848 0

两者使用相同形式的 input.toml：

input = "structure.xyz"
runtype = "ohess"
threads = 1
[calculation]
optlev = "extreme"
[[calculation.level]]
method = "gfn2"
chrg = 0
uhf = 0
accuracy = 0.0001
rddip = true

把这两份文件放在一个新目录后，原生命令为 crest input.toml > native.out 2>&1。包内两个分子已分别准备，结果位于 cases/water 和 cases/hydrogen-sulfide。换成自己的体系时，先核对电荷与自旋，并确保优化后的力已足够小。

3. 找到输出文件#

ls cases/water
cat cases/water/vibspectrum1

crestopt.xyz 是优化结构，numhess1 是 9×9 的笛卡尔 Hessian，单位 Eh/Bohr²；vibspectrum1 给出质量加权并投影后的频率，单位 cm⁻¹；g98.out 可帮助查看振动方向。

非线性三原子分子有九个笛卡尔自由度，其中六个是整体平移和转动，只有三个内部振动。不要把接近零的六个刚体模式当作六个额外的物理振动。

图中箭头来自实际 Hessian 的内部振动特征向量，长度经过放大，只表示运动方向，不是分子的真实热振幅。

4. 区分频率和光谱强度#

3.0.2 这条数值 Hessian 路径在 vibspectrum1 中把内部模式的强度统一写为 1000；g98.out 中部分约化质量、强度等字段也使用固定占位值。不要用这些数字画真实 IR 或 Raman 强度。 本例图片只展示实际计算的频率和振动方向。

没有明显虚频是局部极小点的重要检查；存在明显虚频时，应查看其运动方向，沿该方向轻微扰动后重新优化，或进一步提高数值精度。接近零的整体模式与低频柔性运动要结合体系判断。

5. 导出并重画#

cat results/harmonic-modes.csv
python -B verify.py

CSV 同时给出原生频率，以及从原生 Hessian、独立位移梯度 Hessian 得到的频率。固定小例还在每个分子的九个坐标方向做正负 0.0025 Bohr 位移，以完整梯度核对 Hessian；它帮助理解数值二阶导数的来源，不要求每个日常任务都重复这一套步骤。

这些是 GFN2 教学数值，不代表实验频率或高精度基准。计算热化学修正继续阅读零点能与热化学。方法与版本可参阅官方3.0输入说明及3.0.2 Hessian实现。