Gaussian ONIOM机械与电子嵌入

下载完整输入、选定原生数据与绘图脚本，解压到新目录：

unzip gaussian-oniom-practical*.zip -d oniom-example
cd oniom-example
bash run.sh

这会用随附的实际计算结果生成 analysis/energies.csv、analysis/finite-difference.csv 和 figures/ 下的5张图。脚本需要NumPy和Matplotlib，集群默认使用ASE的Python环境；其他电脑可用 PYTHON=python3 bash run.sh 指定自己的环境。

本组采用个人有权使用的 Gaussian16 C.01 AVX CPU程序，单例1核心、1GB；资料包不包含Gaussian程序。完成自己的Gaussian环境初始化、申请CPU资源后检查入口：

command -v g16
command -v formchk

需要全部重算时，在新解压目录执行 bash run.sh compute，顺序运行23个很小的计算；已有 work/ 同名目录时会拒绝覆盖。该命令不负责申请队列资源，应在已分配的CPU任务中执行。

1. 先区分整数电荷和MM部分电荷#

这里使用前述乙醇实际B3LYP/6-31G(d)计算的Mulliken电荷，固定写入原子行，展示如何输入和比较。它们是教学用固定电荷，不是已经参数化、可直接用于其他体系的UFF电荷模型。

本例原生InputFiles生成的cap部分电荷为C2的 -0.02613744e；输入链接字段写了 H-H_-0.0，也不能据此假定程序生成模型的cap电荷一定为零。应查看实际生成的模型输入和FCHK的 Atom Modified MM Charges，再核对模型部分电荷之和。

表中完整体系电荷之和约为零。模型区MM部分电荷之和不必等于模型的量子整数电荷；模型量子部分在本例仍按输入的 0 1 处理。负部分电荷在原子字段中会出现连续两个减号，如 C-C_3--0.44155841，不应删掉其中一个。

2. 机械嵌入和默认电子嵌入#

机械嵌入完整输入：

%nprocshared=1
%mem=1GB
%chk=result.chk
#p ONIOM(HF/STO-3G:UFF) Force Geom=Connectivity NoSymm SCF=Tight

Small layered-model charged-mechanical-force

0 1 0 1 0 1
C-C_3--0.44155841 0 1.171592164149 -0.403050689766 -0.000000000006 L
C-C_3--0.02613744 0 -0.006726532503 0.556916685590 -0.000000000006 L H-H_-0.0 3
O-O_3--0.61331899 0 -1.201076743568 -0.221755587510 -0.000000000001 H
H-H_-0.38932715 0 -1.953062098208 0.389286604240 0.000000000014 H
H-H_-0.11813140 0 0.039662126550 1.209262019623 0.887716250593 L
H-H_-0.11813140 0 0.039662126550 1.209262019623 -0.887716250593 L
H-H_-0.13965273 0 2.118986840378 0.147543167958 0.000000000001 L
H-H_-0.15788608 0 1.140320061155 -1.044590109792 0.886956786043 L
H-H_-0.15788608 0 1.140320061155 -1.044590109792 -0.886956786043 L

1 2 1.0 7 1.0 8 1.0 9 1.0
2 3 1.0 5 1.0 6 1.0
3 4 1.0
4
5
6
7
8
9

bash compute.sh charged-mechanical-force
bash compute.sh charged-electronic-force
grep -E 'ChrgS2:|ONIOM: gridpoint|ONIOM: extrapolated' work/charged-electronic-force/output.log

第二例仅将路线改为 ONIOM(HF/STO-3G:UFF)=EmbedCharge，其余几何、电荷和方法一致。本小模型的默认ScaleCharge=500会将所有外部嵌入电荷屏蔽掉，因而模型高层能量和机械嵌入相同；不能只看到EmbedCharge关键词就声称量子区已经受到了非零外场。

3. 显式550对照：先看实际生效了哪些电荷#

下面完整输入用 ScaleCharge=550（本身即启用电子嵌入）保留本模型的次近邻电荷；实际有6个外部位置，其中H7、H8、H9三个非零电荷保留下来，总和约0.45542489e。

%nprocshared=1
%mem=1GB
%chk=result.chk
#p ONIOM(HF/STO-3G:UFF)=(ScaleCharge=550,InputFiles) Force Geom=Connectivity NoSymm SCF=Tight

Small layered-model charged-electronic-force

0 1 0 1 0 1
C-C_3--0.44155841 0 1.171592164149 -0.403050689766 -0.000000000006 L
C-C_3--0.02613744 0 -0.006726532503 0.556916685590 -0.000000000006 L H-H_-0.0 3
O-O_3--0.61331899 0 -1.201076743568 -0.221755587510 -0.000000000001 H
H-H_-0.38932715 0 -1.953062098208 0.389286604240 0.000000000014 H
H-H_-0.11813140 0 0.039662126550 1.209262019623 0.887716250593 L
H-H_-0.11813140 0 0.039662126550 1.209262019623 -0.887716250593 L
H-H_-0.13965273 0 2.118986840378 0.147543167958 0.000000000001 L
H-H_-0.15788608 0 1.140320061155 -1.044590109792 0.886956786043 L
H-H_-0.15788608 0 1.140320061155 -1.044590109792 -0.886956786043 L

1 2 1.0 7 1.0 8 1.0 9 1.0
2 3 1.0 5 1.0 6 1.0
3 4 1.0
4
5
6
7
8
9

bash compute.sh charged-electronic-scale550-force
grep -A 16 'Point  2 -- high level on model system' work/charged-electronic-scale550-force/output.log
grep 'ONIOM: gridpoint' work/charged-electronic-scale550-force/output.log

模型高层的变化约 -0.000670384402Ha，说明这次实际存在非零嵌入外场。完整混合能量要组合三个分量；低层模型也会变化，不能只按高层一项的降低判断整体稳定性。两种边界电荷处理属于不同近似，不把它们的差直接解释为某个实物过程的反应能。

4. 读取电子嵌入能量时注意自能约定#

grep -E 'Nuclear repulsion from inactive atom pairs|ONIOM: extrapolated energy' work/charged-electronic-scale550-force/output.log
grep -E '^SCF Energy|^Total Energy' work/charged-electronic-scale550-force/result.fchk

本例外部点电荷对的核排斥自能约 0.0205520089Ha，可由保留电荷的两两库仑和独立算出。最后一个MM步骤的原始能量为 -0.008293140110Ha，减去这项后，才与原生Low(model)分量 -0.028845149054Ha 相符。完整结果仍取原生 ONIOM: extrapolated energy 或FCHK的 Total Energy。

ScaleCharge数字按距模型区的键数控制屏蔽，不代表所有真实边界都应该用550。本例强调先核对生效电荷和原子类型；真实体系应采用经过验证的力场、电荷处理与足够大的QM区。更多规则见Gaussian ONIOM手册。