Gaussian MK电荷与静电势cube实操

下载完整输入、实际数据和绘图脚本，解压到新目录：

bash run.sh

先从已有数据生成 analysis/、CSV和5张 figures/，不启动Gaussian。需要NumPy/Matplotlib，可使用集群ASE 3.29 Python环境，或通过PYTHON指定自己的解释器。

需要从头计算，先初始化个人Gaussian：

yskit gaussian-init

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

bash run.sh compute

顺序运行三份分子输入与formchk、九次Multiwfn电荷分析、三份cubegen静电势；首次还由Multiwfn调用个人Gaussian生成H/C/O三种ROHF原子参考。每步1线程、1GB。软件为个人授权Gaussian16 C.01 AVX厂商预编译CPU版本、内置NBO3.1，以及Multiwfn noGUI2026.10.18；没有公共Gaussian module或自编译器版本可套用。已有输出请另解压，完整work已有时可用 bash run.sh analyze。

公开的 native-fields.fchk.txt 是选定数值/基组/总密度字段，不能代替完整波函数供其他程序任意操作。需要分析自己的体系，请从自己的CHK生成完整FCHK。

1. 先区分拟合电荷与实空间静电势#

Pop=MK得到原子点电荷；cubegen Potential=SCF从自己的密度与原子核生成空间中的静电势。点电荷只是近似模型，两者不是同一种文件。下面使用与原子电荷教程相同的实际水输入：

%nprocshared=1
%mem=1GB
%chk=result.chk
#p B3LYP/6-31+G(d) Pop=(Full,NBO,Hirshfeld,MK) SCF=Tight Int=UltraFine NoSymm

Gaussian population reference water

0 1
O 0.000000000000 0.000000000000 0.000000000000
H 0.000000000000 0.758600000000 0.587100000000
H 0.000000000000 -0.758600000000 0.587100000000

2. 实际终端步骤#

在分配CPU资源并加载个人Gaussian后：

g16 < input.gjf > output.log
grep 'Normal termination' output.log
formchk result.chk result.fchk
grep 'points will be used for fitting atomic charges' output.log
grep 'Charges from ESP fit' output.log
sed -n '/ESP charges:/,/Sum of ESP charges/p' output.log

本例水原生MK使用289个拟合点，打印RMS=0.00336a.u.。水二聚体为453点、0.00390a.u.。这些是Gaussian原生拟合点上的误差；不是下面规则网格上的误差。

3. 生成实际使用的cube#

先创建网格文件：

cat > grid.txt <<'EOF'
-1 -8.125 -6.125 -6.125
-72 0.25 0.0 0.0
64 0.0 0.25 0.0
48 0.0 0.0 0.25
EOF
cubegen 0 Potential=SCF result.fchk potential.cube -1 h < grid.txt
head -n 9 potential.cube

72×64×48个网格点、每轴步长0.25Bohr。负IFlag用于文本输出；网格第一维写负数表示以Bohr给定；向量里的零写成0.0。本例起点偏移以避开核坐标恰好落在网格上。自己的大分子要扩大网格范围，本有限网格甚至没有覆盖水二聚体最上方两个H，不能当作完整分子表面。

图为最接近零平面的实际截面，坐标约−0.066Å；点为实际核坐标投影。水二聚体上方白色区域在本次网格外，投影核点仍标出。色阶只显示±0.12a.u.，超出范围的值按端点色显示，CSV/cube原始数值未截断。图不是电子密度，也不是等密度表面着色图。

4. 点电荷能复现多少静电势#

脚本用每种实际电荷重建 Vq(r)=ΣqA/|r−RA|，只比较距所有核都大于1.5Å的本例网格点；输出在 analysis/RESULTS.json 的 own_outer_grid。原生MK拟合点与此处网格不同，不能拿两个RMS直接互相替代，也不能从三个小体系宣称某种电荷普遍最好。

默认直接出图；自己重算完并具备完整work后：

bash run.sh analyze

继续看电子密度cube教程了解不同场量，别把Potential换个图标题就解释成密度。