ORCA同位素振动与IR实操
复用乙醇实际Hessian，比较氘、氧18和两个碳13位点，给出完整终端流程、谱图与模式对应。
阅读本篇目录
已经完成一次结构优化和频率计算后，研究同位素替换通常可以直接复用Hessian。在当前谐近似中，替换同位素改变核质量，不需要重新求一遍相同电子势能面的二阶导数。这篇用鸢算科技集群已完成的B97-3c乙醇结果，实际比较O–D、全氘、氧18和两个碳13位点。
1. 从已有Hessian开始#
原结构优化、频率与完整ORCA输入见乙醇优化和IR教程。本包直接提供那次实际计算的Hessian，乙醇有9个原子、21条内部振动。
下载同位素振动实例包，执行轻量CPU后处理：
python3 -m zipfile -e orca-isotope-vibrations-wave522-*.zip orca-isotope-vibrations
cd orca-isotope-vibrations
bash run.sh
脚本依次准备9组质量/参数文件，运行原生振动分析、两组完整HESS光谱和两个模式动画，生成数表。程序为 orca/6.1.1 的官方预编译工具，配套OpenMPI4.1.8由GCC13.3构建；本例单进程、单线程，不使用GPU，也不重新计算SCF、优化或Hessian。
完整流程需要新目录。已有自己的Hessian时，保留原件，在副本上按下面的步骤处理。
2. 只改质量那一列#
$atoms 区块每行依次是 元素、质量、x、y、z；质量为u，坐标在这个Hessian中为bohr。查看轻同位素参考：
sed -n '/^\$atoms$/,+10p' inputs/ethanol-light.hess
这组质量明确采用12C、1H、16O；原文件默认的平均原子量单独保留在 average 对照。表中同位素位移以明确的轻同位素参考为准，因此与旧教程默认质量的数值略有差别。
$atoms
9
C 12.00000000000 2.316102230731 -0.723892229615 0.000100485152
C 12.00000000000 0.108651502466 1.075897457126 0.000063870193
O 15.99491461957 -2.163553553208 -0.379385739187 -0.000161857840
H 1.00782503223 -3.571801262684 0.773215370059 0.000111264065
H 1.00782503223 0.193439551855 2.299568494978 1.670651415309
H 1.00782503223 0.193660927859 2.299744186010 -1.670386824450
H 1.00782503223 4.092959176546 0.315691618003 -0.000108278858
H 1.00782503223 2.269697026352 -1.930393766786 1.665682586621
H 1.00782503223 2.269441612156 -1.930586166521 -1.665339556002
O–H上的氢是第4个原子，即程序从0编号的原子3。做O–D替换时，把这行质量改成 2.01410177812，元素标签仍写H。只修改质量，保留原子总数、顺序和其余文件区块。
本包已准备好 inputs/ethanol-OD.hess。其他实际对照为：
| 目录 | 替换内容 |
|---|---|
work/light/ |
全部使用明确的轻同位素质量 |
work/OD/ |
仅O–H上的氢换成D |
work/allD/ |
6个氢全部换成D |
work/O18/ |
氧换成18O |
work/C13-methyl/ |
第1个碳，即甲基碳换成13C |
work/C13-alpha/ |
第2个碳，即邻氧碳换成13C |
3. 重新求振动，不重跑电子结构#
module purge
module load orca/6.1.1
cp inputs/ethanol-OD.hess OD-demo.hess
"$ORCA_ROOT/orca_vib" OD-demo.hess > OD-demo.out
grep -A 32 'VIBRATIONAL FREQUENCIES' OD-demo.out
orca_vib会更新并覆盖传入的HESS副本，写入新的频率、模式和IR数据。它也会统一平移坐标、把输出质量保留到5位小数；分子内部几何和本次Cartesian Hessian保持不变。需要再次使用精确输入质量时，仍从原始副本开始。本包在每组目录保留 before-native.hess，便于核对输入。
更新后，默认投影下仍有6个整体平移/转动零模与21条正的内部振动。显示缩放、投影选项和模式动画见振动文件处理。
4. 用更新后的HESS生成光谱#
"$ORCA_ROOT/orca_mapspc" OD-demo.hess IR -x00 -x14000 -w20 -n4001
head OD-demo.hess.ir.stk
命令产生0–4000 cm⁻¹、4001点的曲线，FWHM为20 cm⁻¹。读取的是更新后的 OD-demo.hess，其中有当前同位素的新频率和模式。独立 orca_vib 的屏幕IR表与常规ORCA主输出的列格式不同，本例绘图采用更新后的HESS文件。
原生HESS谱线文件会列出27行，其中6行整体零模的权重为零；内部振动仍是21条。强度列的含义见IR数据列说明，不要把它直接标成积分强度。
5. 看实际位移与模式对应#
本次按振动形态跟踪O–H/O–D伸缩，结果为：
| 质量组合 | 对应频率（cm⁻¹） |
|---|---|
| 轻同位素参考 | 3756.52 |
| 仅O–H换成O–D | 2735.05 |
| 全氘 | 2735.04 |
| 18O | 3744.25 |
| 甲基13C | 3756.51 |
| 邻氧13C | 3756.51 |
轻同位素参考的O–H模式号是 26，仅O–D替换后对应的是 21，两个位移方向的模式重叠约0.9979。若直接比较“替换前后第26条”，就会把新排序中的C–H伸缩当成O–D伸缩。
热图按Cartesian振动形态匹配全部模式；星号表示重叠低于0.8，需要结合动画判断混合情况。数字相同不保证是同一个物理振动，尤其是在同位素导致模式交叉时。
左图用本次计算的积分IR权重画图，每条谱各自归一化，便于比较峰位和谱形。完整数值、单位和匹配信息在CSV中。
bash analyze.sh work
head tables/frequencies-ir.csv
head tables/matched-isotope-modes.csv
/public/toolkit/miniconda3/envs/ase-3.29.0/bin/python plot_results.py
数表包含波数、模式号、原生T²、换算的积分IR强度与匹配重叠；图片写入 figures/。画图脚本默认读取包内参考数据。需要比较自己的体系时，同时调整原子质量配置和分析目标。
原理与接口见ORCA同位素振动说明。