Gaussian 电场梯度与全电子密度实操

下载9份完整输入、实际数值与6张图的绘图脚本，解压到新目录：

bash run.sh

先从配套数值生成 analysis/、nmr-shifts.csv 与 figures/。需要NumPy与Matplotlib，可使用集群ASE 3.29 Python环境，也可通过PYTHON指定自己的解释器。

从头重算前，先初始化个人Gaussian：

yskit gaussian-init

按Gaussian输入与CPU提交申请CPU资源，在已分配的CPU作业内加载个人环境，再运行：

bash run.sh compute

本例为个人授权Gaussian16 C.01 AVX厂商预编译CPU程序，没有公共Gaussian module或自行编译器/MPI版本可套用。整套输入最多使用4个CPU线程；TMS和显式同位素J计算4线程、2GB，其余1线程、1GB。申请资源须覆盖这些设置。默认入口只读自有数值；重算依次完成9份Gaussian、formchk、3次cubegen与分析画图。已有输出请另解压，完整work已有时可用 bash run.sh analyze。

公开的 native-fields.fchk.txt 是选定数值/基组/总密度摘录，不是完整波函数；分析自己的体系请从自己的CHK运行formchk。本组单点坐标来自已完成的优化；更换分子或优化结果后，要同步更新依赖的单点坐标。

1. 单独请求电场梯度#

本例仍用同一OH几何和UB3LYP/6-31+G(d)，单独请求Prop=EFG：

%nprocshared=1
%mem=1GB
%chk=result.chk
#p UB3LYP/6-31+G(d) Prop=EFG SCF=Tight Int=UltraFine NoSymm

OH separate EFG properties

0 2
O 0.000000000000 0.000000000000 -0.001484568972
H 0.000000000000 0.000000000000 0.981484568494

g16 < input.gjf > output.log
grep 'Normal termination' output.log
sed -n '/Center     Electric/,/Leave Link/p' output.log
formchk result.chk result.fchk

输出包含核位置的电势、电场、梯度、主值及去迹表示。先看清所读的是哪个表；不要把只请求了另一个属性的输出当成已经有EFG结果。

2. 原始矩阵和去迹矩阵要区分#

O的实际原始对角元为−1214.273711、−1218.979970、−1217.181268a.u.，它包含各向同性接触项，不是一个迹为零的矩阵。按本次原生输出，去迹表示满足 T=G−Tr(G)I/3。O的去迹主值为−2.176574、−0.369618、2.546192a.u.。

不要把表头的不同表示当成两套不同计算。脚本按实际矩阵重建主值和去迹关系，并从总密度核对核位置的接触项。这里不把EFG数值直接命名为实验四极耦合常数；那还需要正确的同位素四极矩、单位和符号约定。

3. 密度cube应明确是否包括芯电子#

本篇配套密度图来自同组OH超精细单点的检查点。进入该计算目录，从自己的CHK生成FCHK，然后创建同一网格：

formchk result.chk result.fchk
cat > grid.txt <<'EOF'
-1 -8.125 -7.125 -7.125
-64 0.25 0.0 0.0
56 0.0 0.25 0.0
64 0.0 0.0 0.25
EOF
cubegen 0 FDensity=SCF result.fchk full.cube -1 h < grid.txt
cubegen 0 Density=SCF result.fchk frozen.cube -1 h < grid.txt
cubegen 0 Spin=SCF result.fchk spin.cube -1 h < grid.txt

FDensity明确使用全电子密度，Density形式默认冻结芯处理，Spin是α−β。这里SCF计算本身没有因导出命令而重新进行；改变的是导出密度的选择。Gaussian作者cubegen手册说明了FDensity与Density的区别。

本例两个密度文件即使都带“Total SCF Density”字样，核附近场值仍不同。不能只据标题就断言已经包含所有芯电子。

4. 可视化网格不能直接当收敛电子数积分#

本组64×56×64、0.25Bohr网格用于画图。实际全电子场的有限网格积分约8.493545电子，默认冻结芯场约6.974794，自旋场约1.000101；这里没有把场值归一化成希望得到的整数。独立AO密度矩阵的Tr(DS)仍分别确认9个总电子和1个自旋电子差。

芯电子密度在核附近很尖，规则粗网格积分不能代替收敛积分。需要定量积分时应另做网格/积分方法收敛，不能把整份密度乘一个系数凑成9。图中截面也没有恰好穿过核，峰值不是核位置的完整密度。

python3 - <<'PYCODE'
import json
r=json.load(open('analysis/RESULTS.json'))
print(r['native_fields'])
PYCODE

自己的更大分子要相应扩大网格范围，并保持做差的两个场来自相同几何、相同网格和明确的密度定义。