Multiwfn支持的分析化学键的方法一览

sobereva

PS：本文可能经常更新，内容始终是对于Multiwfn官网上最新版本而言的，不要用老版本程序、看老版本手册！

Multiwfn支持的分析化学键的方法一览

An overview of chemical bond analysis methods supported by Multiwfn

文/Sobereva@北京科音

First release: 2019-Mar-19  Last update: 2024-Jul-12

文章目录
0 前言
1 AIM分析
2 键级分析
3 基于轨道的成键分析：轨道定域化、AdNDP、NBO与NAdO/BOD
4 对ELF的分析
5 其它实空间函数的分析
6 电子密度差分析
7 通过片段电荷与CDA分析定量考察成键伴随的电荷转移
8 ETS-NOCV分析
9 IGM与IBSI分析
10 态密度(DOS)曲线图分析
11 考察键的极性
12 考察键的偶极矩
13 一些其它的分析以及与成键分析有关的量
14 总结


0 前言

波函数分析程序Multiwfn在分析化学键方面有非常强大的功能，支持的分析方法极多，是考察化学键特征必不可少的利器。恰当运用，能帮助你搞清楚这些关于化学键的问题：键是什么类型？键是否形成了？键的强度如何？键的具体构成是怎样的？成键导致了何种电子结构变化？键的极性如何？等等。在本文中，笔者将把通过Multiwfn可以实现的分析化学键的方法进行全面的概述。限于篇幅和笔者的精力，每种分析方法不可能在本文说得特别详细，但本文会明确告诉读者都该具体看什么资料、怎么学相关操作，读者只要稍微花一点时间就可以掌握相关基本知识，并且能够用Multiwfn轻松、顺利地完成相应的分析。只要读者认认真真把本文读完，充分消化本文的内容，那么今后在成键分析方面就会特别得心应手，基本没有搞不清楚的问题，很容易能让文章充实起来。

另外，关于这些分析方法的原理、细节、特点以及在Multiwfn中的使用技巧，笔者会在北京科音每年都举办的“量子化学波函数分析与Multiwfn程序培训班”（http://www.keinsci.com/workshop/WFN_content.html）里非常深入、系统、全面地教授，并给出大量实际例子，十分推荐大家参加。读者还可以下载http://sobereva.com/167里提到的迄今用Multiwfn发表的超过一万篇论文的pdf合集，本文提到的大多数研究方法都已被其中的很多文章利用，可以将这里面的文章当例子库用。

如果读者对Multiwfn不了解，务必参看这些文章以了解基本知识：
《Multiwfn入门tips》（http://sobereva.com/167）
《Multiwfn FAQ》（http://sobereva.com/452）
《Multiwfn波函数分析程序的意义、功能与用途》（http://sobereva.com/184）

尚不会产生Multiwfn输入文件者仔细看《详谈Multiwfn支持的输入文件类型、产生方法以及相互转换》（http://sobereva.com/379）。如此文所示，Multiwfn可以支持几乎所有主流量子化学程序产生的含有波函数信息的文件。做本文提到的这些成键分析都不需要在太高级别下产生波函数。什么级别够用，看《Multiwfn FAQ》中Q35。有些分析方法在使用时不能有弥散函数，这点看《Multiwfn FAQ》中Q36。

Multiwfn在弱相互作用的分析上也极其强大，价值极高。但本文说的都是强相互作用，即化学键强度范畴的。弱相互作用的分析请看笔者专门写的另一篇文章《Multiwfn支持的弱相互作用的分析方法概览》（http://sobereva.com/252）。

PS：本文内容在Multiwfn手册4.A.11节也有，但是是高度简化版，可以推荐给看不懂中文的同行看。


1 AIM分析

AIM全称Atoms in Molecules，是已故的Bader发展起来的一套重要的波函数分析方法，如今已经得到极其广泛的使用，也是分析化学键最流行的方法之一。读者想了解这方面知识的话，参看《AIM学习资料和重要文献合集》（http://bbs.keinsci.com/thread-362-1-1.html）。AIM理论分析化学键主要通过拓扑分析(topology analysis)手段，在Multiwfn中的相关操作看《使用Multiwfn做拓扑分析以及计算孤对电子角度》（http://sobereva.com/108）、《使用Multiwfn+VMD快速地绘制高质量AIM拓扑分析图》（含视频演示。http://sobereva.com/445）、手册4.2.1节。一般人几分钟就能学会基本操作，就是敲几下键盘的事。用Multiwfn做AIM拓扑分析远比用其它程序方便得多、快得多，而且功能强大得多，使用上还灵活得多。

AIM理论定义了临界点(critical point, CP)的概念，相当于某个函数的梯度的模为0的点。AIM分析中通常讨论的是电子密度这个函数的临界点。临界点可分为(3,-3)、(3,-1)、(3,+1)、(3,+3)四种，差异在于临界点处这个函数的Hessian矩阵的三个本征值的符号，括号里面第二个数是本征值为正的数目减去为负的数目。本征值可以理解为这个函数在相应本征矢方向的曲率。

考察原子间的相互作用问题，自然想到要从这两个原子之间相互作用区域的电子结构特征上入手去分析。在形成化学键的原子间，肯定存在一个电子密度的(3,-1)型临界点，这被称为键临界点(Bond critical point, BCP)，在顺着键的方向上这个位置处电子密度的曲率为正，而在与之正交的方向上为负。在AIM理论中认为BCP是原子间相互作用区域中最有代表性的一个点，因此BCP的属性可以用来考察相应化学键的特征，包括强度和本质。在Multiwfn中做电子密度的拓扑分析就可以把BCP的位置，以及BCP的各种性质都计算出来。

很多文献提出了通过BCP的性质考察化学键特征的方法，有的能讨论强度，有的能判断作用类型。下面提及一些相对知名的，这些方法都可以被尝试用于分析你自己的体系。大家多看看上面提到的AIM相关资料的话还会发现更多基于BCP的分析方法。下面提及的各种函数在Multiwfn手册2.6、2.7节也都有简要介绍。

• ρ(BCP)和V(BCP)
BCP位置的电子密度ρ(r)和势能密度V(r)与化学键的强度有密切关系，以下简称ρ(BCP)和V(BCP)，后者体现在BCP位置处的电子感受到的外势，必为负值。对于同类化学键，通常ρ(BCP)越大、V(BCP)越负，则化学键强度越大。（这里说的“同类化学键”指的是特定两种元素组成的化学键，后同。比如A体系的某个C-O键和B体系的某个C-O键可以这么对比，而C-O键和C-F键就不叫同类化学键，不能靠ρ(BCP)和V(BCP)对比强度）

• ▽2ρ(BCP)
电子密度的拉普拉斯（Laplacian）函数，即▽2ρ(r)，定义是这个位置的电子密度Hessian矩阵的三个对角元的加和。负值区域对应电子凝聚区域，如果是成键区域出现了▽2ρ<0的现象，一般可以认为形成了共价键。这个函数在BCP的位置，即▽2ρ(BCP)，常被用来判断化学键的类型。很多文章将▽2ρ(BCP)>0作为化学键是非共价相互作用（也叫闭壳层相互作用）的判据，而将▽2ρ(BCP)<0作为这个键是共价键的判据。这个判据往往是对的，但反例非常多，所以并不可靠。比如对于一氧化碳里的C-O键，明显这是极性共价键，但它的▽2ρ(BCP)为正。这关键在于往往BCP的出现位置不那么合理，不能真正作为化学键作用区域最有代表性的点，关于这一点笔者之前专门写过文章讨论：《AIM键临界点处电子密度拉普拉斯值符号判断相互作用类型失败原因的图形分析》（http://sobereva.com/161）。如果用Multiwfn绘制▽2ρ的图像而不是仅着眼于BCP这一个点的话，就可以避免被▽2ρ(BCP)的符号给坑了。另外，即便是BCP的出现位置有意义的情况，把▽2ρ(BCP)<0当做是共价键的判断标准往往也太严了，特别是对于含有很重（也因此半径很大）的原子根本不适用，这是因为这些原子的成键区域涉及的空间范围很广，导致电子凝聚现象被过度“稀释”了。比如对于I2，虽然I-I是共价键，但两个I之间▽2ρ函数也完全为正值，因此▽2ρ(BCP)>0。在笔者来看，▽2ρ(BCP)<0算是共价键的充分非必要条件。

• H(BCP)
能量密度(energy density)体现的是某个点上电子的能量，是动能密度与势能密度之和，文献中常用H(r)或者E(r)表示。在Angew. Chem. Int. Ed. Engl., 23, 627 (1984)中，作者提出H(BCP)<0时这个键可以被认为是共价作用，而>0则是非共价作用。根据笔者的经验，这个标准虽然不错，但往往有点过于宽松了，比如CaO的Ca-O键的H(BCP)<0，但实际上Ca-O理应被认为是离子键。碰到一个凭直觉模棱两可的情况，如果你倾向于把它说成是共价作用为主，那就用H(BCP)当判断标准吧。在我来看H(BCP)<0适合作为共价键的必要而非充分条件。

• |V(BCP)|/G(BCP)
在J. Chem. Phys., 117, 5529 (2002)中作者提出了通过|V(BCP)|/G(BCP)判断成键类型的标准。文中认为这个量<1、>1但<2、>2这三种情况分别对应键是闭壳层作用、介中作用、共价作用。

• eta指数
J. Phys. Chem. A, 114, 552 (2010)提出了eta指数，在Angew. Chem. Int. Ed., 53, 2766 (2014)等文中也做了进一步讨论。eta指数是个函数，定义为η(r)=|λ1(r)|/λ3(r)。λ1、λ2、λ3是当前位置的电子密度Hessian矩阵的本征值，由小到大排序，后同。BCP位置的eta指数，即η(BCP)，小于1的时候被认为对应闭壳层相互作用，大于1的时候被认为是共价作用为主，而且数值越大共价性越强。不过这也不是完全靠谱的判断标准，比如Ni(CO)4中Ni-C和C-O对应的η(BCP)都小于1，但是它们都应当认为是极性共价作用。

• H(r)/ρ(r)
H(r)/ρ(r)的物理意义是BCP位置上单位电子的能量密度。BCP位置的这个函数被叫做Bond degree (BD)，是在J. Chem. Phys., 117, 5529 (2002)中提出的。文中认为对于共价作用（通常H(BCP)<0），BD越负作用越强；而对于非共价相互作用（通常H(BCP)>0），BD越正作用越弱。

• 键椭率（Bond ellipticity）
Bader在J. Am. Chem. Soc., 105, 5061 (1983)一文中提出了键椭率的概念，定义为ε(r)=[λ1(r)/λ2(r)]-1。BCP位置的键椭率，即ε(BCP)，可以一定程度上反映成键电子偏离轴对称分布的程度。假设电子密度绕着键是对称分布的，比如C-H键，那么λ1(BCP)=λ2(BCP)，正好ε(BCP)=0。而比如对于乙烯，虽然C-C键的sigma电子是轴对称分布的，但是pi电子是明显偏离轴对称分布的，这就导致了ε(BCP)显著大于0。而乙炔的ε(BCP)仍为0，因为其两个pi轨道上的电子合在一起还是轴对称分布的。所以，只要恰当理解BCP位置的键椭率的定义就可以讨论不少问题，比如化学键是否有明显pi特征等等。

• 源函数（Source function）
这个函数写为SF(r,r')，变量包含了两个三维坐标，一般研究的时候将其中r作为参考点，此时SF(r,r')就反映的是r'这个位置对r处的电子密度拉普拉斯函数的贡献（可正可负），见手册2.6节的相应条目的简介，全面的介绍和应用见Struct. & Bond., 147, 193 (2010)。实际用这个函数研究化学键时通常将参考点设在BCP，然后考察体系各个原子对化学键的影响。Multiwfn手册4.17.15给出了源函数的具体分析例子。

此外，还有很多其它文章以各种形式通过BCP的属性来考察化学键问题，比如J. Am. Chem. Soc., 120, 13429 (1998)、J. Comput. Chem., 39, 1697 (2018)等等。还有些文章通过临界点（但不限于BCP）的性质考察晶体的金属性，有兴趣可以看看，比如J. Am. Chem. Soc., 124, 14721 (2002)、Chem. Phys. Lett., 471, 174 (2009)、J. Phys.: Condens. Matter, 14, 10251 (2002)。

只要你把前面提到的Multiwfn做AIM拓扑分析的资料都稍微看了，肯定已经知道怎么计算BCP的上述这些量了。其中大部分是在显示临界点属性的时候直接输出的，有个别的没有输出，但基于上文介绍的定义以及输出信息中的其它量，可以非常容易地计算。或者，去查手册2.7节，看看那些没直接输出的量对应的是第几号用户自定义函数（比如会看到H(r)/ρ(r)对应第17号、V(r)/G(r)对应第35号），启动Multiwfn前把settings.ini里的iuserfunc设为相应的序号，然后用Multiwfn输出临界点性质的时候看显示的User-defined function的值是多少就行了。上面提到了一些基于BCP判断键是共价作用还是非共价作用的方法，但是没有哪种方法对所有体系都是普适的，应当将多种指标结合考虑，或者用后文提到的对ELF等函数来作图的方式考察。

Multiwfn还可以直接给出各个轨道对BCP属性或者三维空间中特定的点的贡献，从而将实空间函数角度的分析与轨道角度的分析有机地结合，这里说的“实空间函数”泛指所有变量是三维空间坐标的函数，后同。具体来说，计算i轨道的贡献，就是把i轨道的占据数设为0，看看属性计算结果相对于原先情况变化了多少，变化量的负值就相当于是i轨道的贡献。有的属性可以这样严格的分解为轨道的贡献，比如电子密度、电子密度拉普拉斯函数，而有的属性比如ELF则不能严格这样分解。恰当利用这种分析，特别是结合Multiwfn的轨道定域化分析，对于讨论化学键特征挺有帮助，看手册4.2.4节的例子。

键径也是AIM分析框架中的重要组成部分，它是以数学方式定义的，是从BCP处出发，分别冲着Hessian矩阵唯一的正本征值对应的本征矢的正、逆方向沿着电子密度梯度最大的方向走出来的轨迹。从实际化学意义上看，键径比较好地展现出了原子间相互作用的最主要路径。形成化学键的原子之间几乎总是会有对应的键径和BCP出现，但出现键径和BCP的两个原子间未必存在化学键作用（明显的弱相互作用都可能没有），因为它们的出现与否不体现相互作用强度的大小，比如哪怕离得很远的两个He间也照样有键径和BCP的出现。Multiwfn的拓扑分析功能不仅可以飞快地产生键径，还能将所有Multiwfn支持的函数在键径上的变化直接绘制成曲线图，看手册4.2.3节的例子。这对于考察化学键不同区域的特征差异、对比不同化学键的细节特征很有帮助，比如可以绘制键径上的键椭率的曲线图，从而了解键的不同位置上电子分布偏离轴对称分布的程度。

此外，要知道键径并不总是一条直线，也可以是弯的，对原子间连线的偏离情况也可以用来讨论键的一些特征。比如环丙烷当中C-C的键径是相对于环中心向外凸的，这体现出环张力导致成键电子主体分部区域向外膨胀（通过对后文提及的ELF、变形密度、电子密度拉普拉斯函数等函数作图也可以体现）。而对于乙硼烷B2H6，硼与体系中央的H的键径是往环中心凹的，体现此体系是个缺电子体系。


2 键级分析

键级对于化学键的分析来说是一个极其重要的概念，它是对化学键特征的一种定量化的表现形式。键级的定义很多，零零碎碎加起来得有二三十种，其中最重要的键级定义在Multiwfn里都直接支持，敲几下键盘就能算出来。Multiwfn支持的键级在手册3.11节都有十分详细的介绍，在手册4.9节有诸多分析例子，读者一定要仔细看一看！看了便知，在Multiwfn里计算键级容易至极。下面笔者对一些与键级相关的较为重要的内容简单说一下。

顺带一提，笔者有一篇重要的博文《谈谈原子间是否成键的判断问题》（http://sobereva.com/414），里面利用了键级进行了分析并做了不少讨论，必看。在《谈谈18碳环的几何结构和电子结构》（http://sobereva.com/515）中笔者也利用了键级，澄清了18碳环体系并非是一些人臆测的那样是单-三键交替的电子结构。

Multiwfn已经支持了将很多下面介绍的键级用于CP2K计算的周期性体系的波函数上面，专门的介绍和例子见《使用Multiwfn对周期性体系做键级分析和NAdO分析考察成键特征》（http://sobereva.com/719）。

• 形式键级
最朴素的键级定义是形式键级，也就是把分子画成Lewis式，两个原子间如果有n个横杠键级就是n。这种键级没有什么实际意义，因为本身画成什么Lewis式就有很大人为因素，尤其是电子高度离域的体系必须得用多个Lewis式的共振才能定性表述其电子结构，只用一个Lewis式对应的形式键级来说事则定性都是错的。NBO理论里面有个自然共振理论(NRT)，可以给出各个Lewis式的权重，把每个Lewis式的形式键级做权重平均，得到的就是NRT键级。这种键级定义虽然物理意义看起来不错，但是对稍大的体系计算极其昂贵，而且Gaussian自带的NBO3.1不支持，还得花钱买之后的NBO才行，所以也很少用，Multiwfn也没打算去支持。NRT键级有个额外优点是可以给出共价性和离子性各自对键级的贡献。

• 分子轨道理论的键级
所谓分子轨道理论定义的键级就是结构化学书里说的成键分子轨道电子数减去反键分子轨道电子数，这种定义非常狭隘。对于双原子分子这个定义没问题，但是对于多原子分子，由于分子轨道往往离域性非常强，比如可能同时主要分布在十几个原子上，缺乏与特定化学键的对应关系，因此绝大多数情况根本没法按照这种方式去算某个键的键级。后来还有人基于自然轨道的占据数去按照这种定义计算键级，叫做有效键级，见Angew. Chem. Int. Ed., 46, 1469 (2007)，当然同样也缺乏普适性。

• Mayer键级
目前最常用的键级是Mayer键级，计算耗时极低，普适性极强（有机、无机、过渡金属间成键等都能用），而且结果易于理解。对于一个Lewis式就可以定性正确描述的化学键，其Mayer键级与形式键级在数值上很接近，比如乙烷、乙烯、乙炔中的C-C键的Mayer键级分别接近于1.0、2.0、3.0，而对于苯这种电子高度离域的体系，我们可以认为每个C-C键是一个sigma和半个pi键的总和，实际算出来的Mayer键级也接近1.5，比较符合化学观念。从物理本质上讲，Mayer键级反映的就是两个原子间共享的电子对数。对同类化学键，Mayer键级和键的强度有一定正相关，但不能用于不同类型化学键之间强度的比较。比如H2的键强度明显的大于Li2，但由于两个体系都可以认为两个原子间完美共享了一对电子，因此二者的Mayer键级基本都为1.0。而且Mayer键级对键的极性不敏感，比如LiF的Mayer键级也和Li2一样接近1.0，因此一般体现不出非极性共价键和极性共价键的差异。在Dalton Trans., 2001, 2095中给出了大量Mayer键级的应用实例，有兴趣可以看看。

• 模糊键级与离域化指数
Mayer键级、模糊键级（Fuzzy bond order）以及离域化指数（Delocalization index, DI）都可以被Multiwfn计算。如J. Phys. Chem. A, 109, 9904 (2005)所指出的，它们仨物理本质其实一样，都是衡量原子间共享电子对数，直接体现电子在两个原子空间之间的离域程度，只不过划分原子空间的方式不同。Mayer键级用的是Hilbert空间划分（即靠基函数与原子的对应关系来划分），模糊键级基于模糊原子空间划分（在Multiwfn里可以用Becke、Hirshfeld、Hirshfeld-I等不同方式定义），原始版本的DI用的是AIM的原子盆方式划分。这三个量对于非极性键分析结果几乎一样，但对于极性键，由于原子空间划分方式不同，结果会有差异，而且极性越大差异往往越大。一般来说，讨论化学键的时候就计算Mayer键级就行了，最省时间，例子在手册4.9.1节就有。但Mayer键级怕弥散函数，因此用了弥散函数时绝对不能用Mayer键级，而模糊键级和DI都不怕弥散函数。一般不推荐计算DI，因为需要做昂贵的AIM盆分析。

• Wiberg键级
常有人在文献中看到Wiberg键级，死活非要算Wiberg键级，到底Wiberg键级是个啥？Wiberg键级最早是Wiberg在Tetrahedron, 24, 1083 (1968)当中的第1086页的脚注中提出的，Wiberg键级提出的时候还处在半经验方法流行的年代，这个键级原本也只能用于基于正交基的半经验波函数。然而我们如今做量子化学计算用的基函数都是非正交的，那么那些文献里的Wiberg键级是怎么算的？实际上，他们是先做基函数的正交化，然后把基于非正交基函数描述的波函数转化到正交基下描述，再套用原始的Wiberg键级的计算公式。然而，基函数的正交化方法是不唯一的，比如有对称正交化、OWSO正交化、Schmit正交化，基于不同正交化方法所产生的Wiberg键级也是明显不同的。因此当我们如今说Wiberg键级的时候，必须说清楚是怎么做的正交化。Multiwfn目前版本给出的Wiberg键级是基于对称正交化的，而NBO程序给出的Wiberg键级是基于OWSO正交化的，所以这两个程序给出的Wiberg键级是不同的，但只要不用弥散函数，二者的结果一般还是相近的。另外，要知道Mayer和Wiberg键级物理本质其实也是相同的，从形式上看，Mayer键级相当于把Wiberg键级扩展到了对应于如今我们做量化计算用的非正交基的情况，因此一般情况就直接用Multiwfn算Mayer键级就完了，非要去效仿一些文献算Wiberg键级是莫名其妙。不过，也有极个别情况，Mayer键级的结果比较诡异，明显违背化学直觉，这时候可以尝试改用Multiwfn计算Wiberg键级或模糊键级。

• Mulliken键级
Mulliken重叠布居(overlap population)也被一些人叫做Mulliken键级，很多人问的“键布居数”一般也是指的这个。这种键级没有什么实际意义，因为和化学键的强度不仅缺乏较好正相关性，甚至有的时候还有负相关性，因此拿它来说事很容易误导。不过Mulliken键级有一个独特的特点是数值可以是负的，两个原子间Mulliken键级为明显负值的时候往往暗示这俩原子间整体显现反键效应、缺乏自发结合的驱动力。

• 拉普拉斯键级
笔者提出过一个键级叫做拉普拉斯键级（Laplacian bond order, LBO），在Multiwfn中的计算例子见4.9.3节。LBO的定义是在原子间模糊重叠空间内对▽2ρ为负的区域中对▽2ρ进行积分，然后乘以10。因为▽2ρ为负的区域体现电子由于共价作用而凝聚的区域，所以这部分的积分可以体现共价作用强度，而乘以10是为了让结果的数量级与形式键级接近，便于讨论。这个键级很有价值，已经被广泛使用。强烈建议大家看LBO的原文J. Phys. Chem. A, 117, 3100 (2013)（https://pubs.acs.org/doi/10.1021/jp4010345），其中对这个键级的思想、特点做了充分的阐述和展示，并且与其它键级用大量体系做了充分的对比，对了解不同键级间的特征差异很有好处。LBO有以下特点：
(1)LBO的大小和形式键级接近（对前两周期非极性键而言），因此易于考察。比如对乙烷的C-C键很接近1.0，对N2很接近3.0。
(2)LBO的大小可以体现键的强度。比如乙烷、乙烯、乙炔的C-C键的键解离能(BDE)比例为1:1.85:2.61，LBO键级的比例为1:1.96:2.56，与BDE相符非常好，而比如Mayer键级是1:1.90:2.41，可见虽然也与BDE有正相关性，但不如LBO关系那么好。再比如N2的LBO是3.05，P2是1.86，如实反映出N-N的键强显著大于P-P的事实。在《一篇最全面、系统的研究新颖独特的18碳环的理论文章》（http://sobereva.com/524）中介绍的笔者的18碳环的研究文章中就用了这种方法估计了18碳环体系中的两类C-C键的BDE，非常方便，详见文章的“Other molecular properties”那一节的相应部分，在文章的补充材料里还给出了ωB97XD/def2-TZVP级别下LBO与BDE的关系图和拟合出的关系式。在笔者的Carbon, 165, 468 (2020)中，LBO对18碳环两类的C-C键特征差异的描述和价层电子密度、Mayer键级、轨道定域化分析都是一致的，也体现出LBO对键强度衡量的可靠性。此文是LBO不错的应用例子，可以在使用LBO时连同原文一并引用。
(3)LBO的大小可以反映出键的极性。比如随着CH3CH3→CH3NH2→CH3OH→CH3F的顺序，取代基原子的电负性逐渐加大，C-X键的极性也越来越大，LBO也是逐渐降低的。因此对于极性键而言，LBO可以认为体现的是其共价部分的贡献。

不过，LBO也有局限性，就是不适合用于离子性太强的键，比如NaCl，也不适合用于含有第四周期及之后元素的成键，因为对于那些半径很大的原子，如前所述，▽2ρ这个函数本身就没法很好表现出共价作用特征了。总的来说LBO最适合考察有机体系，它虽不像Mayer键级那样具有高度普适性，但LBO独特的特征可以起到与Mayer键级互补的作用。

对于LBO适用的场合，由于LBO的大小与键能的正比关系较好，因此对于不太便于计算键能的时候可以借助LBO按比例估计。比如某个环状体系里有A-B键，由于在环里而不好按照常规方式算键能，于是可以对某个含有A-B键的模型体系照常计算键能和LBO，然后计算环状体系里的A-B键的LBO，并按照比例来折算出环里面的键能。诸如我们可以根据乙烷的C-C键键能和LBO值粗略估计环丙烷里的C-C键键能。

• 多中心键级
前面说的都是双中心键级，即用来考察两个原子间成键的，但多个原子间也可以成键，比如苯体系有显著的六中心pi键。为了研究这类问题，有人提出了多中心键指数，我习惯称之为多中心键级（multicenter bond order, MCBO），形式上类似于将Mayer键级扩展到多中心。多中心键级是个非常有用的键级定义，在讨论多中心作用上十分方便，例子见手册4.9.3节和《使用AdNDP方法以及ELF/LOL、多中心键级研究多中心键》（http://sobereva.com/138）。Multiwfn算多中心键级的时候可以用原始的定义，也就是基于基组定义的基函数来算，此时和Mayer键级一样怕弥散函数。Multiwfn也支持基于自然原子轨道(NAO)来计算多中心键级，这样虽然步骤略麻烦一点，但是不怕弥散函数，在原理上也更好点。Multiwfn的主功能15里还能基于模糊原子空间计算多中心离域化指数，见手册3.18.10节的介绍，本质和多中心键级一样，但不怕弥散函数，缺点是耗时更高。

• 键级的分解分析
Multiwfn的键级计算功能绝不仅仅是计算出键级数值那么简单，Multiwfn可以将Mayer键级和Mulliken键级分解成轨道的贡献，前者也叫作“占据数扰动的Mayer键级”分析，见手册4.9.1节和4.9.5节的实例。如果做分解前先做轨道定域化，那么可以直接得到各个成键轨道对键级的影响，这在讨论实际化学键时往往特别有用，在http://sobereva.com/380的例子中利用到了这点。

Multiwfn还可以将多中心键级分解为sigma和pi轨道各自的贡献，对于分析多中心作用本质的时候也十分有益，这里有详细介绍和示例：《在Multiwfn中单独考察pi电子结构特征》（http://sobereva.com/432）。在《通过键级曲线和ELF/LOL/RDG等值面动画研究化学反应过程》（http://sobereva.com/200）中笔者还演示了考察乙炔三聚化反应过程中sigma和pi六中心键级的变化，得到了重要信息，即反应过程中一度出现了sigma六中心离域特征，如果不做这种分析的话很难想象得到这点。在《一篇最全面、系统的研究新颖独特的18碳环的理论文章》（http://sobereva.com/524）介绍的笔者的论文中，通过键级的分解，笔者对具有双芳香性的18碳环体系中的平行于和垂直于分子平面的pi电子分别做了讨论，体现出了键级分解对电子结构特殊的体系研究上的重要价值。

Multiwfn还支持笔者独创的基于NAO对Wiberg键级进行分解的分析方法，能够清楚得知键级主要是由两个原子间的哪些原子轨道相互作用所主要贡献的，这相当于把Wiberg键级的计算从原子-原子层面细化到原子轨道-原子轨道层面，无疑可以把成键本质问题搞得明显更透彻。手册3.11.8节有原理介绍，4.9.4节有这种分析的示例，这方法很有价值，强烈建议使用。

在Multiwfn中可以将键级的计算上升到片段-片段层面，在对应键级分析的主功能8里可以看到用来定义片段的选项-1，在定义两个片段后，计算键级的时候就会把两个片段间每一对原子的键级进行加和并输出，从而衡量片段间总的成键特征。

值得一提的是，考察键级在化学反应过程中的变化，对于了解反应机理、弄清楚电子结构变化特征是极其有帮助的，仔细阅读《通过键级曲线和ELF/LOL/RDG等值面动画研究化学反应过程》后我相信读者就会了解这一点。将Multiwfn与笔者提供的批处理脚本相结合，绘制“键级 vs. 反应坐标”的变化曲线是很简单的事情。做这件事之前先得用Gaussian跑IRC，对这方面不了解的话可以参看《在Gaussian中计算IRC的方法和常见问题》（http://sobereva.com/400）。


3 基于轨道的成键分析：轨道定域化、AdNDP、NBO与NAdO/BOD

• 轨道定域化（Orbital localization）
我们平时做HF、KS-DFT计算，会得到分子轨道，严格来说叫做正则分子轨道，这种轨道往往有很强的离域性，同时分布在一大片原子上，因此一般情况完全没法通过分子轨道讨论成键问题，很多文献里拿分子轨道讨论成键往往也纯属瞎讨论。轨道定域化是一种极为重要的讨论化学键的方法，它将分子轨道通过特定的定域化算法在不失体系波函数物理意义的前提下变换成高度定域的轨道，与化学键能够对应起来，因此就可以将成键特征从轨道层面进行解读。通过这些定域化分子轨道(localized molecular orbital, LMO)，我们可以了解很多成键方面的信息，比如哪些原子间形成了共价键、原子间形成了几重键以及什么样的键、成键轨道是由哪些原子轨道以何种方式构成的、键的极性如何。轨道定域化在《Multiwfn的轨道定域化功能的使用以及与NBO、AdNDP分析的对比》（http://sobereva.com/380）做了十分详细的原理介绍和示例，手册4.19节也有例子，一定要看。用Multiwfn做轨道定域化极其简单，速度很快而且功能强大。

Multiwfn的轨道定域化功能也可以完美地用于考察化学过程。比如手册4.19.2节展示了在SN2反应过程中对应要断的键的LMO和要形成的新键的LMO的变化过程，非常好地从轨道角度描绘出了SN2反应中的电子结构变化。

• AdNDP
轨道定域化得到的LMO主要用来讨论双中心作用、三中心作用，对于更多原子间的成键，轨道定域化往往就没法讨论了，这时候需要用AdNDP方法。AdNDP相对于轨道定域化方法的优点在于能产生可以算作是半定域化的轨道，其离域程度高于定域化轨道，但是又不像分子轨道那么强，这些轨道可以用来讨论多中心作用特征，已经被用于研究金属团簇中的成键、体系的芳香性等问题上。相关介绍和示例见《使用AdNDP方法以及ELF/LOL、多中心键级研究多中心键》（http://sobereva.com/138）。AdNDP的缺点是搜索轨道过程中需要人工挑选轨道，这导致分析结果存在一定程度的主观任意性，对研究者的化学直觉也有一定要求，而且步骤比轨道定域化多很多。如果你研究的问题不牵扯多中心作用，那么用轨道定域化就够了，不需要用AdNDP。或者说，一般情况可以先用轨道定域化，如果发现所得LMO没有展现出对应自己想考察的多中心作用，那才有必要再用AdNDP来分析。

• NBO轨道
值得一提的是NBO分析框架中的NBO轨道分析和轨道定域化分析在很大程度上相似，也是找出高度定域化的轨道来讨论化学键，但算法截然不同，Multiwfn目前不支持这种分析，得用NBO程序实现。对于初学者，拿NBO轨道去讨论问题很容易会被严重误导，特别是碰到离域特征较强的键、电子结构非典型的情况。而且NBO轨道搜索算法的设计本来就有很大人为因素，物理意义并不好。如果你对NBO的理论没有透彻理解，我建议彻底放弃用NBO来讨论成键，否则容易误导自己和他人，而用Multiwfn做轨道定域化分析就没有这个问题。我在《谈谈原子间是否成键的判断问题》（http://sobereva.com/414）中也特别强调不要拿NBO来判断成键与否。不过用NBO轨道讨论也不是完全没有用，有一个优点是可以计算高占据的NBO（即Lewis型NBO）和占据数很低的NBO（即non-Lewis型NBO）之间的超共轭作用能，这在NBO分析里对应二阶微扰校正能，常用E(2)来表示。E(2)的物理意义是：假定体系的电子结构与Lewis型NBO所展现的Lewis式对应的理想电子结构完全相同时，占据的NBO轨道向相邻的非占据NBO轨道发生电子离域而令体系能量的降低，数值通过二阶微扰理论近似估计。恰当讨论E(2)可以用来讨论、解释一些化学现象。但E(2)的数值千万别去过度解释、夸大其意义，毕竟这只是基于非常粗糙的二阶微扰理论估计的。其实很多已发表的文章中的E(2)分析都是瞎分析。

• BOD函数与NAdO轨道
这种分析方法在《使用键级密度(BOD)和自然适应性轨道(NAdO)图形化研究化学键》（http://sobereva.com/535）中做了非常详细的介绍。前文已详细介绍了离域化指数(DI)，它可视为是原子间的共价键级。BOD这个实空间函数描述了三维空间各个位置对DI的贡献，因此它可以将DI通过图形化展现，从而使我们能了解哪些区域的电子对指定的两个原子间的共价作用起主要贡献。NAdO轨道则是进一步通过轨道角度对DI进行描述，所有NAdO轨道本征值的加和正好等于DI，因此通过考察本征值较大的NAdO，我们就可以清晰地搞清楚两个原子间的共价作用主要是由什么样的轨道所贡献的。对于复杂的成键情况，在讨论DI（或模糊键级）的时候，结合使用BOD/NAdO作进一步分析是相当有好处的，能了解更深层信息。在Multiwfn中还把BOD/NAdO分析扩展到了片段间相互作用的研究，极其灵活，比如可以研究二茂铁中Fe与茂环之间的作用。《18个氮原子组成的环状分子长什么样？一篇文章全面揭示18氮环的特征！》（http://sobereva.com/725）介绍的笔者的文章里通过NAdO直观揭示了18氮环中的N-N键的本质，是NAdO极好的应用实例。


4 对ELF的分析

电子定域化函数（Electron localization function, ELF）是一个研究化学体系电子结构超级重要的实空间函数。Multiwfn手册2.6节有其简要介绍，更多的介绍见“ELF综述和重要文献小合集”（http://bbs.keinsci.com/thread-2100-1-1.html），笔者写过很多与ELF相关的博文，链接都在这个合集里面，强烈建议读者观看。笔者在《电子定域化函数的含义与函数形式》（http://www.whxb.pku.edu.cn/EN/10.3866/PKU.WHXB20112786）一文中对ELF的函数形式、物理意义有深入讨论。简单来说，某个点的ELF值衡量的是这个点的电子的定域性高低，值域是[0,1]。ELF数值高的区域一般对应于形成共价作用的区域、孤对电子区域、原子内核和壳层结构区域。Multiwfn是研究ELF非常强大的工具，可以实现很多形式的分析：

• 绘图分析
使用Multiwfn可以将ELF绘制成曲线图、各种类型的平面图以及等值面图，手册4.3、4.4、4.5节分别有相关示例。应当根据要研究的体系特征和具体情况恰当选择作图的类型，哪种形式最能充分、明确展现出你要研究的区域的ELF分布特征就用哪种作图形式，多看看上面提到的手册和相关博文里的例子就明白了。通过ELF图，我们可以立刻了解哪些地方形成了共价键，即哪些原子间形成了ELF较高的区域；也可以考察某两个原子间的键主要是共价键还是离子键，如果是共价键的话那么这两个原子的成键区域应该ELF较高，否则应当较低。ELF既可以分析双中心作用，也可以分析多中心作用。而且在Multiwfn里还可以分解为ELF-pi和ELF-sigma来分别考察pi电子和sigma电子构成的多中心作用，这点在《在Multiwfn中单独考察pi电子结构特征》（http://sobereva.com/432）专门进行了介绍。此文，将Multiwfn和笔者提供的脚本相结合，可以很容易地绘制ELF在IRC或者势能面扫描过程中的变化动画，效果又酷炫又能非常好地展现相应过程中电子结构的变化，详见《通过键级曲线和ELF/LOL/RDG等值面动画研究化学反应过程》（http://sobereva.com/200）。

• 盆分析
对于某个函数，可以通过盆分析来将整个三维空间拆分为一个一个盆(basin)，盆之间通过这个函数的零通量面分隔（垂直于这个面梯度为零），每个盆里面都有一个这个函数的极大点。Multiwfn的盆分析功能是极度普适的，可以用于任意函数。在《使用Multiwfn做电子密度、ELF、静电势、密度差等函数的盆分析》（http://sobereva.com/179）中有详细的原理讲解和实例，手册4.17节也有例子。将盆分析方法用于ELF可以得到很多有价值的定量信息。每个ELF盆都对应于一个有特殊电子结构特征的区域，比如对于每个共价键都有对应的ELF盆。分析在这个盆的空间里的特征，就可以获得与键有关的信息。比如可以对这个盆空间中电子密度进行积分了解有多少电子出现在了成键区域，可以考察这个盆里面的电偶极矩来获得键偶极矩，可以考察这个盆的定域化指数了解成键电子的定域性，等等。如手册4.17.7节示例的，Multiwfn还能给出两个原子各自对它们之间对应共价键的ELF盆中的电子数的贡献，从而考察这个键的极性。

• 拓扑分析
Multiwfn的拓扑分析模块不仅可以用于前面提到的电子密度的拓扑分析，原理上还可以用于所有Multiwfn支持的其它的函数，当然也包括ELF。一般考察的是ELF的(3,-3)和(3,-1)型临界点。前者也叫ELF吸引子(attractor)，对应ELF的极大点；后者也叫ELF二分点(bifurcation point)，是两个ELF盆之间最有代表性的点，这个位置的ELF可以在一定程度上反映在两个盆之间电子离域的难易程度。ELF的拓扑分析思想最早是在Nature, 371, 683 (1994)提出的，后来被很多文章所采用研究化学键以及化学反应过程中电子结构的内在变化，比如RSC Adv., 5, 62248 (2015)、Chem. Phys., 501, 128 (2018)、Comput. Theor. Chem., 1154, 17 (2019)。ELF的拓扑分析在手册4.2.2节有实例，在《在Multiwfn中单独考察pi电子结构特征》（http://sobereva.com/432）也结合实例介绍了将拓扑分析利用在ELF-pi上来考察pi电子共轭程度的高低的方法。值得一提的是，如果你只需要得到ELF吸引子就够了，那么使用Multiwfn做ELF盆分析也可以实现，而且不容易有遗漏。不过由于盆分析是基于均匀格点来实现的，所以给出的ELF吸引子的位置不如做ELF拓扑分析那么精确。

在Multiwfn中对于ELF还可以做一些其它分析，只要你充分理解了ELF的原理、特点，并且把Multiwfn的众多分析功能玩儿转，就可以创造新的分析手段。比如，对于同类化学键，键的强度和原子间的ELF大小往往存在正相关性，那么怎么去试图定量说明？那就可以试图对比原子间连线上或者原子间键径上ELF的最大点的数值，或者尝试用主功能15模糊空间分析功能对原子间重叠空间内的ELF值进行积分，等等。


5 对一些与ELF特征类似的实空间函数的分析

有不少实空间函数，虽然提出的思想、物理意义与ELF可能相差很大，但能说明的问题与ELF有一定程度的相似，或者函数分布特征和ELF密切相关。将上述分析ELF的方法用于分析这些函数上，往往可以获得比分析ELF更好的效果，或者结果能够与ELF分析互相补充，从不同角度更好地认识化学键的特征。这里将这些函数依次说一下，其中重要性低、不常用的就一带而过了。这些函数都可以在Multiwfn中很容易地进行分析。

• 价层电子密度
电子密度的函数分布是相当单调的，一般都是原子核位置是极大点，向四周指数型地降低。直接对电子密度去绘图是考察不出什么对化学上有意义的信息的，对其做盆分析能获得的有用信息也极为有限。笔者在《通过价层电子密度分析展现分子电子结构》（http://www.whxb.pku.edu.cn/EN/10.3866/PKU.WHXB201709252）一文中提出了价层电子密度分析的思想，也就是从总密度中扣除掉和研究化学问题关系甚微的内核电子密度，而只对价层电子密度进行分析。这篇文章强烈建议大家一读，里面不仅对价层密度的行为特征、能研究的问题进行了充分展现，还与ELF、电子密度拉普拉斯函数、变形密度（见后文）通过实例以图形方式进行了对比，对读者认识它们的共性和差异很有帮助。价层密度比ELF的定义简单得多，对大量体系，笔者发现直接绘制价层密度图形就能看出原子间是怎么成键的，还能考察化学过程中电子结构变化，这和ELF分析有不小相似性，而且对价层密度做盆分析还可以比只对ELF做盆分析获得更多的与考察成键有关信息。在Multiwfn中分析价层密度的步骤和分析总电子密度完全相同，唯一区别是在分析之前先要把内核轨道的占据数清零，这个过程可以在Multiwfn的修改波函数的界面（主功能6）里一步自动完成，实例看手册4.6.2节。

• 电子密度的拉普拉斯函数▽2ρ
▽2ρ是AIM分析中经常考察的函数，也经常以作图方式分析其分布特征。一般只需要关注▽2ρ的负值区域，即电子凝聚的区域。多数情况，共价键、孤对电子出现的区域都有对应的▽2ρ为负的区域，因此与ELF有高度的共性。▽2ρ的图像还能一定程度上还原出原子在形成化学键时出现的轨道杂化特征（虽然这只是个虚构的概念）。Bader特意在J. Phys. Chem., 100, 15398 (1996)一文中对▽2ρ和ELF和进行了细致的对比，有兴趣可以看看。▽2ρ和ELF的物理意义并不很相同，因此将它们一起分析，可以提供考察电子结构的更多视角。虽然在展现电子结构上它们往往定性相似，但在具体图形特征上还是有显著区别的，比如对于极性共价键作用区域，ELF图体现在哪里形成了共价作用，但是▽2ρ的图形则倾向于把键的极性特征也展现出来。再比如，对于孤对电子，ELF可以展现出形状显得饱满的高定域性区域，但▽2ρ对应的为负的区域只是一个薄层，没有ELF看起来那么清楚。如之前所述，▽2ρ对于含有半径很大的原子通常不合用，比如对于[Re2Cl8]2-体系中的Re-Re四重键，▽2ρ在两个Re之间数值完全为正，体现不出共价作用，但是ELF在两个Re之间则有函数值明显较高的区域出现，如实展现出了共价作用。

• LOL函数
Localized orbital locator (LOL)可以翻译为定域化轨道定位函数，在手册2.6节有简要介绍。它提出的本来目的是考察哪里是定域化轨道主要出现的地方，但从实际函数分布特征来讲，它和ELF高度相似，函数形式也很相近，而且值域都为[0,1]。基本上所有ELF能分析的问题、展现的电子结构特征用LOL也都可以展现，只不过具体特征上有所不同。虽然用LOL的人目前不如ELF多（毕竟LOL提出的晚），但是与ELF有同等实用价值。对于绘制平面图、等值面图而言，有时候ELF的图像效果更好，但也有很多时候LOL的图像效果更好、讨论问题时更清楚。用哪个更好和具体体系有关，可以都试试，然后用对当前情况效果更好的那个来分析问题。《18个氮原子组成的环状分子长什么样？一篇文章全面揭示18氮环的特征！》（http://sobereva.com/725）介绍的笔者的研究工作里用LOL考察了18氮环的成键和孤对电子分布情况，是LOL的非常有代表性的应用例子。LOL同样可以在Multiwfn中区分为LOL-pi和LOL-sigma来分别研究pi相互作用和sigma相互作用，在《在Multiwfn中单独考察pi电子结构特征》（http://sobereva.com/432）中给出了明确的分析例子。

• IRI函数
IRI是笔者在Chemistry-Methods, 1, 231 (2021)中提出的一种实空间函数，它通常以等值面图的方式来考察。IRI极具实用价值，它可以把化学体系中各种类型相互作用同时非常直观地展现出来，即一张图又可以看到所有化学键作用区域又可以看到所有弱相互作用区域。另外，基于pi电子计算的IRI的变体IRI-pi还可以清晰地展现pi作用区域、作用类型和作用强度。笔者专门写过详细的文章介绍IRI和IRI-pi，读者务必一读：《使用IRI方法图形化考察化学体系中的化学键和弱相互作用》（http://sobereva.com/598），看后你会发现IRI在实际研究中特别有用。

在IRI提出之前还有个Density overlap regions indicator (DORI)函数，是在J. Chem. Theory Comput., 10, 3745 (2014)中提出的，这个函数笔者在《使用DORI函数同时考察共价和非共价相互作用》（http://sobereva.com/367）中专门进行了介绍和示例。DORI和IRI的目的一致，但DORI远不如IRI，不仅定义得特别复杂，而且图像效果巨差，因此在IRI提出之后DORI就毫无使用价值了。

• SCI函数
这个函数全称是strong covalent interaction index，是笔者的合作者刘述斌教授在J. Phys. Chem. A, 122, 3087 (2018)中提出的，之后在J. Mol. Model., 24, 213 (2018)中又做了进一步分析讨论。SCI提出当初侧重于考察金属-金属之间形成的多重键特征。SCI的定义相当于ELF函数当中的一部分，物理思想与ELF有一定差异，没有限制值域为[0,1]。只要恰当设置平面图的色彩刻度、等值面图的isovalue，SCI和ELF可以给出很相似的图像效果。读者在手册2.7节一搜SCI就知道SCI对应37号用户自定义函数，因此要分析SCI的话，把settings.ini里的iuserfunc设为37，然后用Multiwfn绘图、分析的时候被考察的函数选择user-defined function即可。

• SEDD、RoSE、PS-FID
在Chem. Phys. Lett., 582, 144 (2013)中提出的Region of Slow Electrons (RoSE)、在J. Chem. Theory Comput., 10, 3745 (2014)中提出的Single exponential decay detector (SEDD)，以及在Chem. Phys., 435, 49 (2014)中提出的Phase-space-defined Fisher information density (PS-FID)，虽然思想与ELF不同，但是分布特征与ELF比较类似，能讨论的问题也差不多。由于相对于ELF并没有什么显著额外的优点，所以实际意义不大。它们在Multiwfn里都支持，简要介绍见手册2.7节。

• EDR与D(r)函数
EDR全称Electron delocalization range function，最早在J. Chem. Phys., 141, 144104 (2014)提出，形式为EDR(r;d)，其中r是三维坐标矢量，d是个标量。此函数相当于r附近d半径内对费米穴的积分值，以此衡量r处电子在d半径内的离域程度。它相对于ELF的特点是这个半径可以自行调节，而ELF的这个值是固定的，关于ELF的这点参考笔者的《电子定域化函数的含义与函数形式》（http://www.whxb.pku.edu.cn/EN/10.3866/PKU.WHXB20112786）一文。一般来说，D(r)并不比EDR有什么特别的好处，但有时候通过恰当定义d值图像效果更好点，参看Multiwfn手册4.5.6节的例子。
在J. Chem. Theory Comput., 12, 3185 (2016)中，EDR的作者在EDR的基础上还定义了Orbital overlap distance function D(r)函数，可表示为argmaxd_EDR(r;d)，即对每个位置r，通过自动调节d得到最大的EDR值。体系中电子越软（越硬）的区域对应于越大（越小）的D(r)，将这个函数投影到范德华表面上，可以与各个位置的化学反应活性相联系，使用例子见手册4.5.7节，一个应用文章见Angew. Chem. Int. Ed., 56, 6878 (2017)。


6 电子密度差分析

形成化学键的过程中，必定伴随着电子的转移、极化现象。如果想直观地考察这点，那么最好的做法就是考察电子的密度差（electron density difference, EDD）。Multiwfn是绘制和定量分析密度差极为方便、强大的程序。关于绘制密度差笔者专门写过博文进行讨论和示例，务必仔细看看：《使用Multiwfn计算激发态之间的密度差》（http://sobereva.com/429）。密度差有很多种形式，任何形式的密度差都可以用Multiwfn计算，其中与研究化学键有关的密度差有两种：

• 整体与两个片段间的密度差
这种密度差就是用整体的密度减去两个片段各自的密度。计算整体时用的结构一般是优化后的整体结构，而计算每个片段用的结构是直接从优化后的整体结构中提取的，不再单独优化。这样的密度差专门用来考察在两个分子片段之间形成化学键造成的电子密度的变化。比如我们可以考察一个原子团簇化学吸附一个小分子由此造成的电子密度分布的变化。众所周知，共价键的形成伴随着电子往成键区域转移，因此密度差图中能发现两个原子间的成键区域中出现电子密度的显著增加，无疑是两个原子形成共价键的必要条件。

• 变形密度
变形密度（Deformation density）指的是体系的密度减去各个原子在自由状态下的密度，是密度差的一个特例。如果你想考察原本处于孤立状态的原子形成当前分子的过程中电子结构是怎么变化的，那么通过绘制变形密度就可以一目了然地考察。在《通过价层电子密度分析展现分子电子结构》（http://www.whxb.pku.edu.cn/EN/10.3866/PKU.WHXB201709252）文中笔者给出了不少变形密度图的例子，大家看此文就可以了解变形密度的基本特征。

Multiwfn对密度差的分析除了绘制成曲线图、平面图、等值面图外，还可以做盆分析，由此可以获得比如原子间电子密度增加区域的具体增加量，这点对讨论成键很有用，见《使用Multiwfn做电子密度、ELF、静电势、密度差等函数的盆分析》（http://sobereva.com/179）中的例子。

Multiwfn还能把密度差格点数据转化成局部积分曲线和电荷位移曲线（charge displacement curve），这可以非常清晰地定量展现出沿着某个方向（比如键轴方向）的电荷转移量随位置的变化，看手册4.13.6节的例子。比如在J. Comput. Chem., 35, 923 (2014)中作者就通过这个功能Multiwfn考察了AuX (X=F,Cl,Br,I,At)的成键特征，获得了极有价值的信息。


7 通过片段电荷与CDA分析定量考察成键伴随的电荷转移

• 片段电荷分析
两个不同的片段间成键必定伴随着一定程度的片段间电子转移。考察转移量最方便而且易于理解的做法就是计算在整个体系中的其中一个片段的片段电荷(fragment charge)，然后将此数值与这个片段在孤立状态下的净电荷求差。比如CH3NH2里NH2的片段电荷为-0.121，而我们假定NH2在与CH3成键之前是中性状态，即电荷为0，因此我们便知CH3与NH2成键后CH3向NH2转移了0.121个电子。这里说的片段电荷是片段中每个原子的电荷的加和。由于计算原子电荷的方法各有不同，基本相关知识见《原子电荷计算方法的对比》（http://www.whxb.pku.edu.cn/CN/abstract/abstract27818.shtml），所以不同原子电荷计算方法给出的片段间电子转移量必然也有一定差异。Multiwfn支持一大堆原子电荷计算的方法，在手册3.9节有十分详细的介绍，计算例子见手册4.7节。计算之前如果先定义片段，那么原子电荷算完了之后顺带就会输出片段电荷，都省得你去手动加和了，之前笔者录过视频演示：《使用Multiwfn计算分子的某个片段的电荷》（https://www.bilibili.com/video/av26312703/）。视频中用的是笔者提出的ADCH原子电荷（J. Theor. Comput. Chem., 11, 163 (2012), DOI: 10.1142/S0219633612500113），一般分析就用这种电荷就行。注意当体系被划分为多个片段时，基于片段电荷来分析只能得到每个片段的电荷净变化量，它是所有片段间电子转移效果的总和，并没法得知某一对片段间到底电子转移了多少。

• CDA分析
上面这种基于片段电荷的考察转移量的方法很可靠，也省事，但是没法提供片段间电荷转移方面的细节信息。如果要剖析得更深入透彻，那就需要用电荷分解分析（Charge Decomposition Analysis, CDA）了。此方法最初由Gernot Frenking在1995年提出，后来笔者提出了广义化的版本GCDA并实现进了Multiwfn。GCDA解决了CDA存在的很多局限性，详细介绍见《广义化的电荷分解分析(GCDA)方法》（http://dx.doi.org/10.12677/japc.2015.44013）。笔者在《使用Multiwfn做电荷分解分析(CDA)、绘制轨道相互作用图》（http://sobereva.com/166）中对CDA的原理和在Multiwfn中的使用做了详细介绍，给出了具体例子，很多CDA分析的要点在《广义化的电荷分解分析(GCDA)方法》中也做了细致说明。简单来说，CDA的强大之处在于：
(1)可以将片段A向B的供给量和B向A的反馈量单独给出
(2)可以告诉你电子转移是由哪些片段轨道之间的相互作用导致的
(3)当体系被划分为超过两个片段时，CDA可以告诉每一对片段间电荷转移情况
而且，在Multiwfn中，还可以告诉你整个体系的分子轨道是怎么由片段轨道构成的，并且可以直接绘制轨道相互作用图，例如笔者的一篇文章RSC Adv., 5, 78192 (2015)。恰当利用CDA分析，显然可以深入地从电荷转移的角度探究成键本质问题。


8 ETS-NOCV分析

ETS-NOCV分析是一种极为重要、有不可替代价值的考察两个或多个片段间轨道相互作用的方法。它可以把轨道相互作用能分解成为不同NOCV pair的贡献，而每个NOCV pair的密度又对应了它对电子密度差的贡献，从图形上可以判断NOCV pair对应了什么形式的轨道相互作用。因此ETS-NOCV可以清晰地解释片段之间都有什么样的相互作用、对应的相互作用强度是多少，这样的信息对于考察化学键本质，以及弱相互作用中的轨道相互作用部分，都是极具价值的。此方法目前已经很流行，并且由于Multiwfn的ETS-NOCV分析功能使用特别简单方便，日后用的人会越来越多。关于这种分析，笔者写了一篇超级详细的文章进行了全面的介绍并提供了大量例子，务必一看：《使用Multiwfn通过ETS-NOCV方法深入分析片段间的轨道相互作用》（http://sobereva.com/609）。看过之后你会深感这种分析实在太有用了。


9 IGM与IBSI分析

IGM分析是个相对比较新的分析方法，2017年才提出，笔者在《通过独立梯度模型(IGM)考察分子间弱相互作用》（http://sobereva.com/407）中对此方法的思想做了十分详细的介绍，并且给出了不少在Multiwfn中分析的例子。IGM分析方法提出的本意主要是用来分析弱相互作用，但实际上它在考察化学键方面也有用。IGM方法基于电子密度定义了一个函数δg，它可以展现出现原子间相互作用的区域，数值还能区分相互作用强度。弱相互作用的原子间δg的大小远低于化学键相互作用的原子间的情况，而且对于同类相互作用，一般相互作用越强δg也越大。因此δg可以用来展现出哪里存在化学键作用，在需要对同类键横向对比大小的情况也可以拿δg来说事。在Multiwfn中，还可以基于IGM的思想给出某一对原子作用强度的定量数值，而且还能通过着色来直观表现（不过这种分析本质上是基于几何关系估计的，并不体现电子结构，所以也不要太当回事）。

在J. Phys. Chem. A, 124, 1850 (2020)中IGM原作者在IGM框架下提出了IBSI（intrinsic bond strength index），可以翻译为内禀键强指数，介绍见3.11.9节。在JPCA文章中体现出对于共价键，IBSI与键的强度（原文通过力常数衡量）有一定正相关性，另外弱相互作用、过渡金属配位键、离子键和普通共价键的IBSI值相差很大，因此在一定程度上可以区分作用类型。其实这个特征与前述的拉普拉斯键级比较相似。通过Multiwfn计算IBSI的例子见手册4.9.6节。


10 态密度(DOS)曲线图分析

Density-of-states (DOS)图常见的类型有Total DOS (TDOS), Partial DOS (PDOS)和Overlap DOS (OPDOS)，都可以用Multiwfn基于量子化学程序输出的含有波函数信息的文件轻松地绘制。Multiwfn还能绘制相对少见的Local DOS (LDOS)，它和扫描隧道显微镜的观测结果有密切关联。各种DOS的相关原理介绍、具体的绘制和分析例子在《使用Multiwfn绘制态密度(DOS)图考察电子结构》（http://sobereva.com/482）都非常详细地给出了。

上述DOS类型中，对于考察成键有直接用处的是OPDOS，它可以考察不同能量区域的分子轨道对你指定的两个片段间起到成键作用、反键作用还是基本没有影响。实际上OPDOS就是对各个轨道计算自定义的两个片段间的Mulliken键级，然后再把数据展宽成曲线得到的。如果某个区域OPDOS曲线是明显正值，就说明这个能量区间的分子轨道对两个片段的结合是产生促进作用的，即起到成键轨道的效果；如果为明显负值就说明是对结合不利的，暗示起到反键轨道的效果；而如果数值接近0，即对结合不产生直接影响，就说明起到的是非键轨道的作用。在Multiwfn的DOS绘制界面里，片段既可以被定义为一批原子，也可以定义为比如一批原子轨道，一批壳层，或某种角动量的轨道，十分灵活，不同的定义方式可以获得不同的信息。通过OPDOS图对于考察分子特别是原子团簇的电子结构很有帮助，比如Multiwfn原文J. Comput. Chem., 33, 580 (2012)里就有用OPDOS图分析二茂铁中Fe与茂环间成键的例子，十分建议一看。另外，PDOS对讨论电子结构也是很有用的，比如在RSC Adv., 5, 78192 (2015)和J. Comput. Chem., 38, 1574 (2017)中笔者通过PDOS图考察了特殊类型的原子团簇。


11 考察键的极性

键的极性是很多人感兴趣的问题，利用Multiwfn，可以试图通过以下方式考察键的极性。

• 计算两个原子各自对它们之间的化学键对应的LMO的贡献。假设一个是贡献百分比是ΘA，另一个贡献ΘB，那么这个键的离子性百分比就是|ΘA-ΘB|，而共价性百分比就是1-|ΘA-ΘB|。显然，如果一个键是非极性的，则离子性百分比为0%。离子性百分比越大，键的极性就越高。这种考察方式通常还是比较靠谱的。比如CH3F里，对应C-F键的LMO中F贡献69%，C贡献31%，因此离子性百分比就是38%。再比如，NaCl里，对应Na-Cl的LMO中Na贡献15%，Cl贡献85%，因此离子性百分比就是70%。注意计算结果受到多方面影响，包括产生波函数的级别、做定域化轨道的算法、计算轨道成份的方法。从逻辑上，姑且可以认为离子性百分比超过一定阈值就可以将这个键指认成离子键了。但是一刀切式地设定阈值总是不科学的，总会有反例，更何况计算方法也严重影响结果，没法说哪种做法就一定最合理。在Multiwfn以这种方式分析键的极性很容易，用主功能19做轨道定域化之后，找出对应成键的轨道的序号，然后进入主功能8去计算轨道成分即可（做轨道定域化过程中也会直接输出轨道成份，但是是通过廉价、可靠性相对低一些的Mulliken方法给出的，因此建议自行用更好的方法算。轨道成份计算知识看《谈谈轨道成份的计算方法》http://sobereva.com/131）。

• 计算能量指数(energy index, EI)并进而计算键极性指数(Bond polarity index, BPI)。这个概念是J. Phys. Chem., 94, 5602 (1990)中提出的，详见手册3.200.12节的介绍和4.200.12节的例子。以这个方式考察键的极性稍微麻烦点。EI相当于原子在分子环境中每个价电子的平均能量，可以通过Multiwfn直接得到，然后通过手动将参考物质中原子的EI和实际体系中要研究的化学键涉及的原子的EI按照定义求差就可以得到BPI。BPI大小越大极性就越强。

• 做ELF盆分析，找到对应共价键对应的ELF盆，然后用Multiwfn计算两个原子各自对这个盆里的布居数的贡献，然后二者之差就可以衡量键的极性，具体例子看手册4.17.7节。这种考察方式物理意义也比较明确，但比考察LMO成份的方法更麻烦，而且键的极性比较大的话可能根本不存在对应共价键的ELF盆，届时也没法这么分析。

• 前面提到过，拉普拉斯键级体现的是原子间共价部分，属于共价键级，而Mayer键级则可视为总键级，而且对于前两周期元素形成的非极性键，二者大小几乎是一样的，因此二者之差也可以用来体现键的极性大小，差异越大极性越大。

此外，用后期收费的NBO程序计算NRT键级，用给出的离子性部分占总键级的百分比来体现极性也是一个办法。


12 考察键的偶极矩

偶极矩是分子的可观测量，键的偶极矩这个概念也挺有意义，和键的极性理应有密切关系。但键的偶极矩是不可观测的，计算方法也不唯一，不同方法给出的结果也明显有异。有以下方式可以试图考察键的偶极矩问题：

• Multiwfn的轨道定域化分析功能不仅可以给出各个LMO上的电子对体系总偶极矩的贡献，还会给出以笔者独创的方式定义的“LMO偶极矩”。对于被判断为单中心的LMO，LMO偶极矩体现这个LMO的电子云的质心对所属原子的原子核位置的偏离。而对于被判断为双中心的LMO，LMO偶极矩体现的是这个LMO的电子云质心偏离两个原子间成键区域中点的偏离，中点通过原子坐标和共价半径定义。笔者通过对一些体系的测试，发现这个量能较好地体现出键的极性。原理见手册3.22节的介绍，实际分析例子见手册4.19.4节。

• 做ELF的盆分析可以得到对应共价键的盆，对这个盆计算偶极矩即可考察键偶极矩。具体来说，如手册4.17.2节对ELF盆计算多极矩的输出所示，输出是这样的
Basin electric dipole moment:
X=   -0.104745  Y=   -0.104745  Z=    0.020728  Magnitude=    0.149575
Basin electron contribution to molecular dipole moment:
X=    0.000000  Y=    0.000000  Z=   -2.388409  Magnitude=    2.388409
其中“Basin electric dipole moment”是以这个盆的吸引子为中心算的偶极矩，这个可以作为一般意义的键偶极矩，和键的极性往往是正相关的。但实际上一个键对体系偶极矩的贡献并不等于这种方式算的键偶极矩，因为这个键自己带的净电荷也会对体系偶极矩产生贡献。这个盆里的电子对体系总偶极矩的贡献是“Basin electron contribution to molecular dipole moment”下面输出的。

• Multiwfn可以把体系偶极矩在Hilbert空间中分解为各个原子的偶极矩和每一对原子间的键偶极矩的加和。在Multiwfn手册里3.200.2节有详细介绍。

顺带一提，后期NBO程序里有个DIPOLE关键词可以给出占据的NLMO的偶极矩，对所有占据的NLMO的这个量加和正好等于体系偶极矩。看起来好像挺好，但实际上你将之用于实际体系就会发现很坑，比如给出的对应C-H键的NLMO的偶极矩远远大于O-H的，然而显然应该后者的极性远大于前者，可见DIPOLE分析结果和一般化学观念上的键的极性完全没有相关性。


13 一些其它的分析以及与成键分析有关的量

• Potential acting on one electron in a molecule (PAEM)
PAEM是个实空间函数，某个点的PAEM值的物理意义是一个电子在这个点感受到的总的外势，相当于静电势与交换相关势的总和。杨忠志等人提出通过分析PAEM来考察原子间是共价还是非共价作用，分析例子和详细介绍看Multiwfn手册4.3.3节。

• 键能
这通常被认为是衡量键强度最客观的量，而且可以在不同类型的键之间对比。这不是Multiwfn计算的，但可以被所有主流量化程序计算。更多相关讨论看《通过柔性力常数考察键的强度》（http://sobereva.com/364）。

• 柔性力常数
对于化学键而言，柔性力常数相当于键伸缩方向上的力常数，这被一些研究者认为是比键能更严格的考察化学键强度的方法，详见上述《通过柔性力常数考察键的强度》。这个量可以用compliance程序计算，在这篇文章里也有计算示例。

• 基于片段的能量分解

将相互作用能分解成不同的物理成分有助于理解相互作用的本质，这叫能量分解。能量分解的种类很多，在《Multiwfn支持的弱相互作用的分析方法概览》（http://sobereva.com/252）中有一节做了简要概述。我非常推荐我于2023年提出的sobEDA能量分解，又快又普适又可靠，已经被大量文章使用，见《使用sobEDA和sobEDAw方法做非常准确、快速、方便、普适的能量分解分析》（http://sobereva.com/685）。

• 基于原子的能量分解
上面这些分解是基于片段的，还有一种叫做Interacting Quantum Atoms (IQA)的能量分解，它可以给出每一对原子间相互作用的各种成分，用于考察强相互作用和弱相互作用、分子间和分子内作用都可以。但有一点很关键，也就是IQA给出的原子间总相互作用能，跟它们间的键能并没有密切关联。目前Multiwfn不支持IQA，主要也是笔者感觉这方法在实际中也没什么太大用处，得不到什么对讨论化学键特别有帮助的信息，而且耗时超级高，只能搞得动极小体系。还有一种叫Mayer能量分解的方法，和IQA在一定程度上类似，但具体实现明显不同，也是可以给出每一对原子的作用能，耗时低得多，用于大体系也可以，见《Mayer能量分解及APEX4程序使用简介》（http://sobereva.com/90）。不过笔者也不觉得此方法给出的数据有太大实际价值。

• 基于键长的分析
在《谈谈原子间是否成键的判断问题》（http://sobereva.com/414）中笔者谈到了如何基于共价半径、键长的关系十分粗略地判断是否成键，但这种做法完全依赖于阈值，人为任意性太强。在DFT-D3原文中提出了一种计算原子配位数的方法，笔者将之改造成了基于键长和原子共价半径判断原子连接性的函数I，在Multiwfn里可以用主功能100的子功能9直接计算，见手册3.100.9节的介绍。当原子间距离大致等于或小于两个原子的共价半径和的时候，这个函数接近1.0。当键长进一步增加，此函数会逐渐平滑降低，到距离达到共价半径和两倍时，此函数几乎为0。I函数不能作为键级来用，因为不体现单键、双键之类，但很适合用于考察键的存在性、完整性（或者说被破坏的百分比）。由于这种分析只依赖于几何结构，所以xyz、pdb、mol等任何Multiwfn支持的含有坐标信息的文件都可以当输入文件。


14 总结

本文全面总结了Multiwfn支持的十分丰富的化学键的分析方法，对于需要其它程序才能实现的分析方法也做了简要提及，这些方法可以说给化学键的研究提供了尽可能全面的视角。相信读者读过本文之后，足矣充分认识到Multiwfn在化学键分析上的关键地位，是实际研究几乎必不可少的程序。只要把本文介绍的分析方法充分掌握，把相关手册章节和博文里的例子都彻底搞明白，只要有个像样的文章素材，借助Multiwfn很容易就能搞出一篇内容丰富充实的文章。由于Multiwfn程序极度灵活，功能极多，显然绝对不可能各种分析、各种情况都会在手册或博文里给出例子，否则就没完没了了，所以大家务必把例子的操作思路，每一步的意图领会清楚，举一反三，这样不仅可以驾驭好各种现有的方法，甚至大家还能把Multiwfn玩出更多花样、提出新的分析方法。

本文的大部分分析方法都不仅限于基态，只要给Multiwfn的输入文件里有激发态波函数，照常分析就可以考察激发态的成键问题，对比基态和激发态的分析结果就可以明确获知电子激发（或发射）对体系成键特征产生了什么样的影响。在手册4.18.13节里明确给了一个简单的激发态波函数分析的例子。

本文的大部分分析方法也不仅限于极小点结构，也可以用于IRC上的不同的点、扫描过程中产生的各个结构，乃至从头算动力学产生的轨迹中的不同的帧。将这些分析方法用于考察这些过程上，可以获得比分析单一结构明显更多的化学上感兴趣的信息。

用于分析的波函数既可以由量子化学程序在真空中计算产生，也可以在溶剂模型下产生（一般使用隐式溶剂模型）。比如在Gaussian中，用了scrf=solvent=xxx关键词，将得到的fch/wfn/wfx文件弄到Multiwfn里分析，考察的就是在xxx溶剂环境中溶质中的成键特征。

Multiwfn的开发不会止步，笔者会不懈地令Multiwfn支持更多对于研究化学键有用的方法，届时本文也会更新，更新日期会体现在本文开头。请注意时常关注Multiwfn主页上的Update history了解最新进展。

Multiwfn主要面向孤立体系的电子结构分析，比如分子、团簇。对于通常基于平面波来计算的周期性体系，笔者打算以后另开发程序。但有很多这种体系其实可以先用第一性原理程序优化结构，然后挖出团簇模型、对边界进行恰当处理，在弄到量化程序里来算单点能得到波函数，然后照常用Multiwfn分析成键。

lijiayisjtu

受益匪浅！

第十九路猪猴

终于等来了这么详细的一篇博文了，先mark下，慢慢阅读

Tanya

太有用了，谢谢社长！

sobereva

在第三节末尾加入了BOD、NAdO轨道分析方法的说明

sobereva

今天更新的Multiwfn新支持了IBSI指数的计算，相应地给本文增加了以下内容

在J. Phys. Chem. A, 124, 1850 (2020)中IGM原作者在IGM框架下提出了IBSI（intrinsic bond strength index），可以翻译为内禀键强指数，介绍见3.11.9节。在JPCA文章中体现出对于共价键，IBSI与键的强度（原文通过力常数衡量）有一定正相关性，另外弱相互作用、过渡金属配位键、离子键和普通共价键的IBSI值相差很大，因此在一定程度上可以区分作用类型。其实这个特征与前述的拉普拉斯键级比较相似。通过Multiwfn计算IBSI的例子见手册4.9.6节。

sobereva

在此文中补充了以下内容

• EDR与D(r)函数
EDR全称Electron delocalization range function，最早在J. Chem. Phys., 141, 144104 (2014)提出，形式为EDR(r;d)，其中r是三维坐标矢量，d是个标量。此函数相当于r附近d半径内对费米穴的积分值，以此衡量r处电子在d半径内的离域程度。它相对于ELF的特点是这个半径可以自行调节，而ELF的这个值是固定的，关于ELF的这点参考笔者的《电子定域化函数的含义与函数形式》（http://www.whxb.pku.edu.cn/EN/10.3866/PKU.WHXB20112786）一文。一般来说，D(r)并不比EDR有什么特别的好处，但有时候通过恰当定义d值图像效果更好点，参看Multiwfn手册4.5.6节的例子。

在J. Chem. Theory Comput., 12, 3185 (2016)中，EDR的作者在EDR的基础上还定义了Orbital overlap distance function D(r)函数，可表示为argmaxd_EDR(r;d)，即对每个位置r，通过自动调节d得到最大的EDR值。体系中电子越软（越硬）的区域对应于越大（越小）的D(r)，将这个函数投影到范德华表面上，可以与各个位置的化学反应活性相联系，使用例子见手册4.5.7节，一个应用文章见Angew. Chem. Int. Ed., 56, 6878 (2017)。

sobereva

Multiwfn已经支持了将很多类型的键级用于CP2K计算的周期性体系的波函数上面，专门的介绍和例子见《使用Multiwfn对周期性体系做键级分析和NAdO分析考察成键特征》（http://sobereva.com/719）。