Gaussian 自旋污染与波函数稳定性

下载完整输入、实际数值、六份轨道/密度网格和绘图脚本，在新目录解压后：

bash run.sh

先用已有结果生成 analysis/、figures/，无需运行Gaussian。程序使用集群ASE 3.29 Python环境中的NumPy、Matplotlib。

要重新计算，先在自己的账号初始化Gaussian：

yskit gaussian-init

按Gaussian输入与CPU提交准备CPU作业，把执行命令设为 bash run.sh compute。已经取得CPU作业资源并加载个人环境后，也可在另一份新解压目录执行：

bash run.sh compute

它依次运行29份Gaussian输入、formchk和6次cubegen，随后从自己的 work/ 分析和画图。客户正式计算应在分配到的CPU资源内执行。已有 work/analysis/figures 时另解压一份；只有完整 work/ 而没有分析、图片时，可用 bash run.sh analyze 接着处理。

本例采用Gaussian 16 C.01 AVX厂商预编译CPU程序，个人初始化，没有公共Gaussian module；没有自编译器/MPI版本需要套用。输入每份1线程、1GB，按实际作业资源调整。本例不代表GPU或Intel Xeon 8173M节点性能。

1. 电荷、多重度与实际波函数分别检查#

输入 0 1 表示中性、请求单重态；0 3表示中性、请求三重态。多重度 2S+1对应理想纯自旋态的 S(S+1)：单重态0、双重态0.75、三重态2。UHF单行输入多重度并不能保证最终波函数就是纯自旋态。

本例以H₂从短键到断键的8点对照，再看16电子O₂三重态的UHF与ROHF，全部采用STO-3G。比较目的是学习诊断，不用于预测收敛的O₂光谱或实验能级。

2. RHF检查与UHF优化的完整输入#

4 Å H₂的RHF检查：

%nprocshared=1
%mem=1GB
%chk=result.chk
#p RHF/STO-3G SCF=Tight Stable NoSymm

Finite active-space and spin practical h2-rhf-4p0

0 1
H 0.000000000000 0.000000000000 -2.000000000000
H 0.000000000000 0.000000000000 2.000000000000

同一几何的UHF轨道稳定化：

%nprocshared=1
%mem=1GB
%chk=result.chk
#p UHF/STO-3G SCF=Tight Guess=Mix Stable=Opt NoSymm

Finite active-space and spin practical h2-uhf-4p0

0 1
H 0.000000000000 0.000000000000 -2.000000000000
H 0.000000000000 0.000000000000 2.000000000000

先在分配的CPU作业中运行，再查输出：

g16 < input.gjf > output.log
formchk result.chk result.fchk
grep -i 'instability\|stable' output.log

查看自旋行时使用固定字符串：

grep -F 'S**2' output.log

Stable只检查；Stable=Opt可以调整电子轨道。两者都不移动本例的原子。SCF收敛、轨道稳定、几何优化和纯自旋是不同问题，正式任务分别确认。

3. before与after不能混为一谈#

H₂在4 Å的UHF原生日志显示，自旋湮灭前 S**2=1.0000，打印的湮灭后为0.0000。根据未投影的α/β密度独立计算为0.99998006，与before一致。

不能因为after打印为0就把这一份UHF电子能、密度和轨道当成纯单重态结果。 断键UHF能够降低受限解的能量，同时可能明显破坏自旋对称性。若目标依赖纯自旋态、断键和近简并组态，结合活性空间方法及适合的后续理论判断。

RHF的自旋数接近0也不代表它在长键处稳定。本例RHF在较长键上恰好演示了这一点。

4. O₂三重态：UHF与ROHF的实际输入#

UHF：

%nprocshared=1
%mem=1GB
%chk=result.chk
#p UHF/STO-3G SCF=Tight Stable=Opt NoSymm

Finite active-space and spin practical o2-uhf-mult3

0 3
O 0.000000000000 0.000000000000 -0.605000000000
O 0.000000000000 0.000000000000 0.605000000000

ROHF：

%nprocshared=1
%mem=1GB
%chk=result.chk
#p ROHF/STO-3G SCF=Tight NoSymm

Finite active-space and spin practical o2-rohf-mult3

0 3
O 0.000000000000 0.000000000000 -0.605000000000
O 0.000000000000 0.000000000000 0.605000000000

O₂共16电子，0 3在本例对应9α与7β。ROHF限制成对电子共用轨道，保留未成对轨道；方法和实际支持的后续相关处理需要配套选择。

UHF的能量更低并不自动意味着某一种方法普遍更准确；本例只对同一有限基组、固定键长和具体解比较。ROHF这一输入没有请求Stable，不把UHF的稳定性结果迁移给ROHF。

5. 下载包的自旋计算适用于哪里#

basis_density.py根据实际FCHK的s/sp高斯基组独立构造AO重叠矩阵；对HF单行列式，使用α/β密度和重叠矩阵计算S²，同时核对 Tr(DS) 的电子数。没有把自旋密度的网格积分直接当作S²。

这条单行列式公式不能用于CASSCF的一阶约化密度矩阵。下载包对CASSCF只标明输入请求的自旋，没有宣称从一阶密度独立测得其总S²，也没有计算二阶约化密度矩阵。脚本仅覆盖本例真实s/sp基函数，不能直接拿去处理任意d/f/ECP任务。

head analysis/energies.csv

逐点完整数值在 analysis/RESULTS.json。

操作顺序为：确定电子数和多重度→完成SCF→检查稳定性→读取before自旋污染→结合目标选择受限/非受限/多参考方法。活性空间的具体示例见下一篇。