原子电荷计算方法的对比

[Review]物理化学学报(WuliHuaxueXuebao)ActaPhys.-Chim.Sin.2012,28(1),1-18JanuaryReceived:September13,2011;Revised:October25,2011;PublishedonWeb:October31,2011.∗Correspondingauthor.Email:[e-mail];Tel:[phone].TheprojectwassupportedbytheNationalNaturalScienceFoundationofChina(20773011).国家自然科学基金(20773011)资助项目ⒸEditorialofficeofActaPhysico-ChimicaSinica原子电荷计算方法的对比卢天陈飞武*(北京科技大学化学与生物工程学院,北京100083)摘要:原子电荷是对化学体系中电荷分布最简单、最直观的描述形式之一,在理论和实际应用中都有重要意义.本文介绍了12种重要的原子电荷计算方法的原理和特点,通过大量实例从不同角度比较了它们的优缺点.这些方法包括Mulliken、分子环境中的原子轨道(AOIM)、Hirshfeld、原子偶极矩校正的Hirshfeld布居(ADCH)、自然布居分析(NPA)、Merz-Kollmann(MK)、分子中的原子(AIM)、Merck分子力场94(MMFF94)、AM1-BCC、Gasteiger、电荷模型2(CM2)以及电荷均衡(QEq)方法.最后本文对如何在实际应用中选择合适的计算方法给出了建议.关键词:原子电荷;计算化学;布居分析;电负性;静电势中图分类号:O641ComparisonofComputationalMethodsforAtomicChargesLUTianCHENFei-Wu*(SchoolofChemicalandBiologicalEngineering,UniversityofScienceandTechnologyBeijing,Beijing100083,P.R.China)Abstract:Atomicchargeisoneofthesimplestandthemostintuitivedescriptionofchargedistributioninchemicalsystems.Ithasgreatsignificanceintheoryandinpracticalapplications.Inthisarticleweintroducethebasicprinciplesandspecialcharacteristicsoftwelveimportantcomputationalmethodsforthedeterminationofatomicchargesandcomparetheirprosandconsfromvariousaspectsbyconsideringalargenumberofinstances.ThesemethodsincludeMulliken,atomicorbitalsinmolecules(AOIM),Hirshfeld,atomicdipolemomentcorrectedHirshfeldpopulation(ADCH),naturalpopulationanalysis(NPA),Merz-Kollmann(MK),atominmolecules(AIM),Merckmolecularforcefield94(MMFF94),AM1-BCC,Gasteiger,chargemodel2(CM2),andchargeequilibration(QEq).Finallysomegeneralsuggestionsonhowtochooseapropermethodforpracticalapplicationsaregiven.KeyWords:Atomiccharge;Computationalchemistry;Populationanalysis;Electronegativity;Electrostaticpotential1引言原子电荷,即位于原子中心的点电荷,是对化学体系中电荷分布最简单、最直观的描述方式之一.它有很多重要意义,比如帮助化学工作者研究原子在各种化学环境中的状态、1考察分子性质、2,3预测反应位点.4另外原子电荷模型在计算化学中也有很多实用价值,如作为分子描述符用于药物虚拟筛选、5在分子对接和分子动力学/蒙特卡罗模拟中描述静电作用、6在量子力学(QM)与分子力学(MM)结合计算中表现QM区域原子对MM区域原子的静电作用势.7在其它一些理论方法中也需要借助原子电荷,如Pipek-Mezey轨道定域化方法、8doi:10.3866/PKU.WHXB20122811ActaPhys.⁃Chim.Sin.2012Vol.28电荷自洽的紧束缚密度泛函方法(SCC-DFTB)、9Truhlar课题组开发的一系列溶剂模型(SM)10-13等.原子电荷并不是可观测量,也没有客观、唯一的定义,因此获得原子电荷的方法多种多样.通过一些实验数据可以间接地考察原子电荷,14如通过分子的多极矩、红外光谱强度和频率、配体场分裂能、核磁共振(NMR)位移等,但是数据的获取相对不便,与原子电荷的对应关系存在较大经验性,也不能分析不稳定的体系和状态.计算化学的发展使原子电荷方便、快速、可靠的获得成为了可能.15-18自从1955年Mulliken电荷19-21被提出以来,迄今已有不下60种原子电荷计算方法被相继提出,近年仍有许多研究者在改进计算方法.虽然在许多计算方法的原文和部分文献22-25中都对不同类型原子电荷进行了比较和讨论,但是涉及的方法并不全面,测试分子和测试内容缺乏统一性.因此,很有必要对较为重要的、文献中常涉及的原子电荷计算方法进行综合的比较,使它们的特点、优缺点能够清晰地展现,这将有助于在实际问题研究中选择适合的计算方法.我们选取Mulliken、分子环境中的原子轨道(AOIM)、26Hirshfeld、27原子偶极矩校正的Hirshfeld布居(ADCH)、28自然布居分析(NPA)、29Merz-Kollmann(MK)、30分子中的原子(AIM)、31Merck分子力场94(MMFF94)、32,33AM1-BCC、34,35Gasteiger、36电荷模型2(CM2)37和电荷均衡(QEq)38共12种方法.在文中将首先介绍这些方法的基本原理和特点,之后通过实际体系测试和比较它们的各方面性能,最后讨论计算量的大小,并给出方法选择上的建议.2各种原子电荷计算方法的原理和特点2.1Mulliken方法Mulliken方法19-21是最古老的原子电荷计算方法.它的算法简单,计算量可忽略不计,几乎出现在所有量子化学软件中.首先考虑分子轨道波函数归一化条件(假设为实函数,后同)∫ϕi(r)2dr=1(1)r代表空间坐标.将分子轨道ϕi按以原子为中心的基函数χm展开ϕi=∑mCm,iχm(2)这里C代表系数矩阵.将式(2)代入式(1),积分后得∑mC2m,i+2∑m∑mnCm,iCn,iSm,n=1 (3)其中Sm,n=∫χm(r)χn(r)dr.上式中第一项可视为各个基函数对轨道独立贡献之和,称定域项;第二项是交叉项,体现了由于每一对基函数之间的耦合对轨道产生的联合贡献.Mulliken将分子轨道i中基函数m所占成分定义为Θm,i=C2m,i+∑n≠mCm,iCn,iSm,n(4)也就是说,定域项被完全划归到相应基函数,而交叉项被平分给对应的两个基函数.将所有轨道中属于相同原子的基函数的占据数加和就得到了原子占据数,并直接得到原子电荷qA=ZA-∑iηi∑m∈AΘm,i(5)其中Z是原子核电荷数,η是轨道占据数.Mulliken电荷存在一些严重问题:(1)交叉项平分的做法物理意义不严格,没有考虑到原子间的差异性;(2)基组依赖性非常大.尤其是使用大基组、含弥散函数基组的情况下问题更为明显.这主要是因为弥散函数延展范围大,容易破坏基组的平衡性,后有一些研究者提出了不同方法试图改进交叉项划分的不合理性,39-43其中也包括使用相对较多的Löwdin方法,39但是严重的基组依赖性却并没消除,也并未普及开来.2.2NPA方法自然布居分析(NPA)29的关键是将带有随意性的原始基组(一般为扩展基)描述的波函数转化到物理意义清晰的正交极小基下描述,使基函数与原子轨道有明确对应关系,很大程度避免了基组不平衡性问题对结果的影响,也避免了划分交叉项的困难.NPA的原理如下.首先根据原始基函数与原子的对应关系将单粒子密度矩阵中对应于各个原子的对角块依次提取并对角化,所得本征向量就是初自然原子轨道(PNAO),其本征值就是PNAO的占据数.根据占据数可将PNAO分为两类:(1)极小集轨道.它们拥有较高占据数,是对原子上的电子密度最主要、最紧凑的描述,与基态原子的内层和价层原子轨道一一对应.(2)里德堡集轨道.这是指PNAO中除极小集轨道之外的、电子占据数很少的轨道,它们对原子的电子密度的描述不起主要作用.随后对这些PNAO进行占据数权重的对称正交2卢天等:原子电荷计算方法的对比No.1化(OWSO),使它们不仅在原子内正交也在原子间正交.相对于其它正交化方法,OWSO方法可以使有意义的极小集轨道变形尽量小,而允许意义不大的里德堡集轨道自由变化,物理意义更明确.PNAO在正交化后被称为自然原子轨道(NAO).以NAO为基重新构建体系单粒子密度矩阵,将对应于相应原子的矩阵对角元加和就得到了原子占据数.值得一提的是,由于NPA是自然键轨道(NBO)方法44框架内的一个组成部分,而且NPA电荷通常由Weinhold等开发的NBO程序45计算,所以NPA电荷经常在文献中被称为NBO电荷.2.3AOIM方法和NPA的出发点相近,黎乐民等26,46提出的分子环境中的原子轨道(AOIM)方法也是将原始扩展基下的波函数转化为极小基描述后再做布居分析.AOIM有两种不同的具体实现方法:第一种方法46认为分子环境会使得原子轨道相对于原子在自由状态时收缩或扩张,分子环境的影响通过平均化的球对称有效外势来表达VʹA(r)=(1/4π)∫V(r,θ,φ)sinθdθdφ(6)其中r,θ和φ是球坐标,球坐标系以原子A的原子核为中心,V是分子势场.在这样的有效势下解原子的径向薛定谔方程,所得波函数结合角度波函数就是分子环境下的原子轨道.实际求解时使用变分法,将计算分子时用的扩展基的径向函数作为变分函数,就得到了原先扩展基与新的极小基之间的变换关系.第二种方法26利用Sanchez-Portal投影方法47,48将扩展基波函数投影到Slater极小基轨道(STO).设Löwdin正交化后的扩展基为φortho,系数矩阵为Cʹ;所期望的极小基在Löwdin正交化后为χortho,系数矩阵为Cʹ(mini),则变换关系为Cʹ(mini)m,i=∑n〈χortho􀰙φorthon〉Cʹn,i.基函数的约化会导致基组完备性下降,从而使波函数信occ息丢失,衡量丢失量的参数为Δ=∑i[1-∑m(Cʹ(mini)m,i)2].它是STO的指数与方向的函数.为了最大程度避免因约化带来的基组不完备性增加,通过有限内存大尺度限制性优化(L-BFGS-B)方法最小化Δ就得到在分子环境下膨胀、收缩和旋转后的原子轨道.通过上述过程得到极小基下的系数矩阵后,就可以照常按照Mulliken方法计算原子电荷.下文中涉及的AO-IM方法都特指第二种AOIM方法.2.4AIM方法Bader的分子中的原子(AIM)方法31,49是典型的基于实空间划分的计算原子电荷的方法.AIM方法将电子密度零通量面定义为原子间的分界面,划分出的每个原子独立的空间被称为原子盆.在分界面上没有电子密度梯度线穿过,即满