PySCF教程＆测评

PySCF是用Python语言写的量子化学程序。程序提供简单，轻量级的高效平台，用于量子化学计算和代码开发。当前版本可以在一个SMP节点上有效地操作中等尺寸体系（大约3000基函数，1000电子）。

（以上复制自量子化学网，这个网站上收录了很多奇怪的程序）

PySCF的安装非常简单，运行pip install pyscf即可。用PySCF进行计算需要编写Python脚本。在这里以苯乙炔的单点能为例，

将pyscf中的gto, dft等模块导入后，首先通过gto定义分子构型。gto的属性包含了基组(basis),电荷(charge)，自旋多重度(spin)，同位素(mass)等，定义方式非常自然。例如基组既可以用上述字符串定义，如果希望使用混合基组，也可以用数组定义。以下是官方文档给出的例子。

PySCF的官方文档上写分子构型中不同原子可以通过分号或换行分隔，所以如果直接读取xyz文件的话应该可以不作处理。不过由于字符串中间直接换行是不合语法的，像现在这样手动输入坐标的话就还是要在换行处补加斜杠了，因此Rita也顺便加上了分号。PySCF对笛卡尔坐标和内坐标均支持。

PySCF官方文档并没有列出其内置的基组。Rita一开始尝试在basis中指定6-311+G(d,p)报错后在报错信息中找到了内置基组的存放位置（/usr/local/lib/python2.7/dist-packages/pyscf/gto/basis/），找到了内置基组的列表。

对cc和def2系列基组支持相当完备，cc支持到5-zeta，带弥散的只支持无止尽的八月。Def2系列包含一系列RI辅助基组；完整内置pc和pcseg系列。这些基组的输入形式比较自由，如cc-pvtz, ccpvtz等均可。

Pople系列基组的文件名称就显得有些诡异。经过尝试，输入6-31G(d)可以识别，但6-311+G(d,p)就不行，要写作6-311+G**.

执行上述文件输出到test1.log里，可以看到计算得到的能量，用时80s. 默认的输出详细度下只输出了能量的数值。

为了一探Gaussian为什么比PySCF慢这么多，Rita把PySCF的输出详细程度调到最大(scf.verbose=9，一般设到4或5即可得到需要的信息)，同时也看一看PySCF的“noisy”模式会输出什么。

每一次迭代都会输出这样一大串信息。由此可以看出SCF收敛限为1E-9，比Gaussian的默认收敛限还要严格一个数量级；整个任务经过11次迭代，比Gaussian少4次，不过也不足以解释巨大的时间差。因此PySCF的速度优势应当是在算法的某处。

所以即使Gaussian的双电子积分几乎是令人绝望的优秀，但在其他地方的效率还是有很多可以优化之处的。顺带一提PySCF的DFT默认没有开RI（可以通过with_df.auxbasis属性设置），如果再开RI应该还会更快。

接下来测试UDFT。首先是苯乙炔的三重态，用PBE0/def2-SVP算单点。注意这里的mol.spin指定的是alpha电子减去beta电子，因此等于自旋多重度-1.这里输出详细度指定为4，就只输出了每一轮迭代的能量。计算在22 s完成，能量为-307.661565969255。Gaussian在相同水平下的计算用时51 s，能量-307.665874496。

在官方手册上波函数稳定性分析部分神奇地留空了。在下方给出了一个自动调节自旋态的功能，在这个例子中虽然没有指定氧气的自旋态，但如果使用UHF，通过addons的设置可以自动收敛到三重态基态。Rita很想试验一下这个功能，但发现似乎只能用于UHF，不支持将UKS的mf对象作为参数（Rita之前的代码中把scf用作了对象的名称，但现在如果import scf的话就只能用其他名称了）。由于手册里相应部分留空，Rita很好奇PySCF里是否有波函数稳定性分析的功能。

PySCF的DFT有一个好玩的特性，即可以自定义泛函。对于范围分离泛函可以通过omega属性指定omega值。不过手册里没有提到双杂化泛函，由于其使用了Libxc的泛函库，Rita去查了Libxc的手册，其中也没有双杂化相关内容。看起来是不能使用双杂化了。

PySCF的RI效率如何呢；在一开始的代码里SCF部分改成如下内容：

from pymfimport dft

mf=dft.RKS(mol)

mf.xc='pbe0'

mf.verbose=4

mf =mf.density_fit()

mf.kernel()

以默认的辅助基组进行RI.用时只降低了4 s。在mf =mf.density_fit()后面写mf.with_df.auxbasis = 'def2-tzvp-ri'指定辅助基组，耗时也只从80 s降低到了74 s。不过不加RI本身就已经这么快了，这样的成绩也算是可以接受。

然后是PySCF的CASSCF功能。仍然以苯乙炔为例，考虑到CASSCF一般比较慢，所以基组小至6-31G(d)，活性空间取(4,4)。写法如下：

相较之下Gaussian在(4,4)的活性空间下用时119 s，收敛失败。

xswl

也许我不应该用Gaussian这样的辣鸡CASSCF软件和PySCF去比？连ORCA的CASSCF都比Gaussian好得多(

可是我没有molpro谁有molpro可以给我一个吗

以及PySCF也支持NEVPT2; ORCA的垄断被打破了(

除此之外PySCF里还有一些有用的特性：

可以自定义哈密尔顿

可以生成任意原子轨道积分；此外PySCF还提供了libcint积分库的接口可以自己退（划掉）自己用

可以做轨道定域化

可以算周期性体系（手册上写了很多，太长不看，不知道这个功能怎么样）

可以提取解析梯度(手册上解析梯度只介绍了HF，而Heissian虽然有这一条目但也是空的。Rita尝试写grad.RKS但报错了)

PySCF不支持的功能

隐式溶剂化

既然梯度和Heissan那么有限应该不能用来做几何结构优化和频率计算

之前提到过的双杂化泛函

耦合簇只支持CCSD和CCSD(T)，不支持DLPNO（但支持F12）

不能并行

由此来看，Pyscf如果要与主流量化软件抗衡还远远不够。不过由于其开源以及Python易于扩展的特性，很高的计算效率，以及可以调泛函调积分等功能，作为新的泛函等的开发工具，或是为一些需求不超出其功能范围的程序提供快速计算，无疑是十分有吸引力的。

原文为PySCF官方文档。