计算化学4从头计算法原理ppt课件

1、袁淑萍 国家重点实验室培育基地408室 E-mail: 电话：85950690计算化学,Gaussian 程序的输入文件,chk=water.chk % Section： 行首以%开始，段后无空行 %rwf=water.rwf #p hf/6-31g scfcyc=250 scfcon=8 Route Section：行首以#开始，段后加空行 Water title：作业的简要描述，段后加空行 0 1 Molecular Specification： O 分子说明部分，段后通常加空行 H 1 R1 H 1 R1 2 a1 R1=1.04 a1=104.0,Gaussi

2、an作业的格式,Molecular Specification,分子说明部分主要用来定义分子核相对位置 分子核相对位置可以用笛卡尔坐标，内坐标(Z-matrix)，或者是二者混合表示 笛卡尔坐标是内坐标的一种特殊形式 分子坐标的格式为: 元素符号, x, y, z 元素符号 (n) 原子1 键长 原子2 键角 原子3 二面角 格式代码 (0,180) 以原子在分子中的序数表示 通过右手规则确定,几何优化用内坐标比较容易到达终点,长度单位一般默认为， 个别软件使用 nm 角度单位一般为度(,原子和分子在体系中的位置,分子的构建,描述,能量,笛卡儿坐标 内坐标,分子结构数据的获得,简单分子,大分子

3、(如蛋白质分子,通过标准键长、键角、二面角等数据直接写出 通过分子图形软件构建,复杂分子,有单晶结构数据的，直接转化为分子结构数据 采用软件构建，必须通过适当的构象搜索方法方可得到合理的结构,其结构的构建本身就是专门的计算化学研究的课题，只有通过特定的方法(如同源模建等)方可得到可能的合理结构,ENT，PDB Brookhaven Protein Data Bank MOL (Molecular Design, Ltd. MDL) HIN Hyperchem GJF Gaussian INP Gamess ZMT MOPAC Z-Matrix CAR Accelrys, Insight II

4、CSD Cambridge Structural Database Crystallographic files MOL2 Tripos MSF Accelrys, QUANTA MSI Accelrys, Cerius2 CPD Accelrys, Catalyst CHM queryAccelrys, Catalyst query CIF Crystallographic Information file XYZ genericgeneric SKC, ISIS MDL, ISIS SMI SES fileSMILES file,常见结构数据类型,笛卡儿坐标输入,C -0.25 -0.47

5、 0. 0.00 O -0.25 0.96 0. 0.00 H 0.77 -0.83 0.00 H -0.77 -0.83 0.89 H -0.77 -0.83 -0.89 H -1.16 1.28 0.00,用文本编辑器直接输入笛卡儿坐标 原子符号 x坐标 y坐标 z坐标,C 0.0 0.0 0.0 0.0 0.0 0.0 O 0.0 0.0 1.221.22 0.0 0.0 H 0.9 0.94 0.54 0.0 H -0.94 0.54 0.0,原子符号可以是元素符号，也可以是原子序数(核电荷数)， 也可以采用在元素符号后加上原子的序列编号,C10.00.00.0 O20.01.220.

6、0 H30.940.540.0 H4 -0.940.540.0,60.0 0.00.0 80.0 1.220.0 10.940.540.0 1 -0.940.540.0,内坐标输入 (Z-Matrix,当没有构建分子的图像软件可用时，与直角坐标相比使用内坐标要方便得多。 对具有较高对称性的分子，输入内坐标Z-矩阵还可减少优化变量、节省计算时间并保证优化结果与分子的对称性严格相符,C C 1 1.54 H 1 1.09 2 109.47 H 1 1.09 2 109.47 3 -120. H 1 1.09 2 109.47 3 120. H 2 1.09 1 109.47 3 180. H 2

7、1.09 1 109.47 6 -120. H 2 1.09 1 109.47 6 1120,Z-矩阵输入步骤,选择起始原子 第二个原子与第一个原子之间的距离， 确定第三个原子与前面任一原子间的距离以及三个原子所组成的 “键角” 从第4个原子开始，需要有距离，“键角”和二面角来定义,特别注意键角的取值范围： 0键角180,应当注意： 尽管在这里用了“键长”、“键角”等术语，而且在多数情况下，它们也确实是化学上的键长、键角，也是最方便的作法。 然而，从Z-矩阵的作用来讲，本质上只要求具有数学意义而不一定要有化学意义，即：它们只须代表两原子间距离及三原子间的角度，而不管这些原子是否相互键连,右旋法

8、则”确定二面角符号 1-2-3-4 右手拇指指向23的方向，14方向与其余四指方向相同， 二面角为正，否则为负,Z矩阵中的数值也可用变量表示,H O 1 R1 O 2 R2 1 A H 3 R1 2 A 1 D Variables R1 0.9 R2 1.4 A 105.0 D -120.0,内坐标和笛卡儿混合坐标,O 0 -1. 0. 0. C 0 0. yc 0. C 0 0.-yc 0. N 0 1. 0. 0. H 2 r1 3 a1 1 b1 H 2 r1 3 a1 1 b2 H 3 r1 2 a1 1 -b1 H 3 r1 2 a1 1 -b2 H 4 r2 2 a1 3 d3 V

9、ariables: yc 1. r1 1.08 r2 1.02 a1 125. d3 160. b1 90. B2 90,例1：使用HF方法，优化 H2O2分子,chk=h2o2.chk %rwf=h2o2.rwf #p hf/6-31g opt H2O2 energy calculation 0 1 H O 1 0.9 O 2 1.4 1 105.0 H 3 0.9 2 105.0 1 120.0,chk=h2o2.chk %rwf=h2o2.rwf #p hf/6-31g opt H2O2 energy calculation 0 1 H 0.000 0.0000.000 O 0.000

10、0.9000.000 O 1.350 1.2620.000 H 1.464 1.742-0.752,内坐标表示 笛卡尔直角坐标表示 上述两个作业将对H2O2分子的结构进行完全优化，包括所有的键长键角和二面角。 通过将内坐标定义成变量，可以对分子结构进行部分优化,chk=h2o2.chk %rwf=h2o2.rwf #p hf/6-31g opt H2O2 energy calculation 0 1 H O 1 r1 O 2 r2 1 a1 H 3 r1 2 a1 1 d1 r1=0.9 r2=1.4 a1=105.0 d1=120.0 通过对两个键长和键角使用同一变量定义可以控制分子的对称性

11、,chk=h2o2.chk %rwf=h2o2.rwf #p hf/6-31g opt H2O2 energy calculation 0 1 H O 1 r1 O 2 r2 1 a1 H 3 r1 2 a1 1 d1 r1=0.9 a1=105.0 变量部分 d1=120.0 r2=1.4 常量部分 通过设定常量，可以使程序 只对分子结构进行部分优化,例2： C2H4分子 平面形分子, HF/STO-3G OPT C2H4 opt 0 1 C C 1 r1 H 1 r2 2 a1 H 1 r2 2 a1 3 180.0 H 2 r2 1 a1 3 0.0 H 2 r2 1 a1 4 0.0

12、r1=1.32 r2=1.09 a1=120.0,通过同一个变量控制C-H键长 把二面角定义为180和0的常量来控制分子的平面构型,例3： CH3F分子 分子点群C3V, HF/STO-3G OPT CH3F C3v opt 0 1 C F 1 r1 H 1 r2 2 a1 H 1 r2 2 a1 3 b H 1 r2 2 a1 3 -b r1=1.38 r2=1.09 a1=110.6 b=120.0 通过常量b来控制C3V对称性,例4： NH3分子 分子点群C3V,使用虚原子,为避免使用180键角和方便定义带有对称性的 分子，可使用虚原子。虚原子用X标记,N X 1 1. H 1 nh 2

13、 hnx H 1 nh 2 hnx 3 120.0 H 1 nh 2 hnx 3 -120.0 nh 1.0 hnx 70.0,例6： 乙炔分子 直线形分子, HF/STO-3G OPT C2h2 opt 0 1 C X 1 B1 H 1 B2 2 A1 C 1 B1 2 A2 3 D1 H 4 B2 2 A3 1 D2 B1 1.20 B2 1.07 A1 120.0 A2 60.0 A3 120.0 D1 180.0 D2 -180.0,特别注意 在Gaussian的输入文件中，必须区分实数与整数形式，所有的xyz坐标以及内坐标的数值必须为实数类型(即必须带小数点)。如键长为1.0，为简单

14、可写为“1.”，但绝对不可写做“1”,作业,写出分子 H2C=C=O分子输入坐标 其中C=C: 1.35 C=O: 1.20 C-H: 1.09 写出苯分子内坐标 C=C: 1.30 C-H: 1.09, rhf/3-21g opt CH2CO opt 0 1 O1 C2 1 1.20 X3 2 1.20 1 90.0 C4 2 1.35 3 90.0 1 -180.0 H5 4 1.09 2 120.0 3 180.0 H6 4 1.09 2 120.0 3 0.0,计算化学,4.从头计算法原理及应用-2 Principle and Applications of SCF Ab Initio

15、 Method,2. 基函数和基组 Basis Functions and Basis Ste in Ab Initio Method,LCAO-MO近似,对分子体系进行从头计算时，对多原子分子,其MO必须用原子轨道的线性组合(LCAO)来近似表示,基函数,基函数的选择对从头算结果影响很大。 在模型设计正确的前提下，基函数选择得好，结果就能预测几何构型，或指出反应途径，解释实验现象反之则很难作为讨论依据。 下面介绍基函数选用的几种函数,称为基函数,1. 类氢离子波函数,早期的量化计算，人们采用类氢离子波函数为基函数。 优点是波函数与原子轨道有一一对应关系，能较好地描述电子在空间的分布,类氢离子

16、波函数的具体形式,径向部分,角度部分,关联Laguerre多项式,Laguerre函数,基函数的角度部分,基函数的角度部分,基函数的径向部分,用类氢离子波函数做分子轨道计算时，由于径向部分的关联Laguerre多项式积分难以计算，迭代收敛很慢。 Slater发现，类氢离子波函数的计算瓶颈在径向函数部分，即加以改造，提出新的径向函数，角度部分仍沿用球谐函数，后人称之为Slater函数,2. Slater Type Orbital (STO,Slater建议用于原子轨道径向部分的函数,Slater函数对径向部分的指数加以内层电子静电屏蔽校正，并改造成无径向节面的函数形式。轨道指数,Z为核电荷，s为

17、屏蔽常数，n*为有效量子数,与类氢离子轨道相比，采用Slater型轨道可以大大减少计算积分时数学处理的复杂性。 但由于STO没有径向节面，不同的STO之间一般不正交，单一的STO轨道还不能很好地逼近原子轨道,具体的Slater函数如,采用 Slater函数的线性组合，则能更好地逼近Hartree-Fock原子轨道,一般来说，每个原子轨道采用3-5个STO基函数，即能得到接近Hartree-Fock极限的精确的原子轨道。 Watson计算了第一过渡金属的原子轨道，它包含10个STO； 对于许多平衡态问题，只需用两个STO就可以得到满意的结果,半经验算法常用STO基组，从头算开始也采用STO基组,

18、1. 单x基： STO与原子轨道一一对应SZ,2. 双x基： 每个原子轨道对应一对内（I）、外（O） STO轨道DZ,3. 扩展基： 含三个以上STO轨道，在双x基基础上加极化轨道。如用p型函数描述s轨道的极化，用d型函数描述p轨道的极化。 DZP,根据一个原子轨道使用STO的个数，可将STO分成三类,双x STO较单x有较大改善。较适应由分子环境引起的变化。环境会使原子轨道变形，内层电子变化较小，外层变形较大。因此双x 把原子轨道分成内外两部分处理，更容易满足上述要求,与真实的轨道形状接近 对氢原子, 它是精确的 通常用于原子计算 有正确的渐进性质(r) 在近核部分表现很好，在核处满足尖点条

19、件(r0) 相互不正交，需要计算大量的三中心、四中心双电子积分，要用1/r12无穷极数展开，使运算变得十分复杂,STO的特点,3. Gauss Type Orbital (GTO,Boys在1950年提出用Gauss函数拟合Slater函数，计算积分。 Gauss函数用球坐标表示为： X=Rn,l(r,)Yl,m(,) 角度部分仍用球谐函数，径向部分表达为 Rn,l(r,)=2n+1(2n-1)!-1/2(2)-1/4rn-1e-r2 =D rn-1e-r2,若用笛卡儿坐标表示,X,i，j，k 的意义： i+j+k=l 角量子数,这里的i，j，k 不同于量子化学中的主、角、磁量子数的概念。 如

20、当i=j=k=0，称为1s型Gauss函数,6个Gauss函数可重新组合使之与5个d原子轨道对称性相符; 同理， i，j，k 之和为3，可构造10个4f函数，重新组合后可得到7个与原子轨道对应的f轨道,常用的GTO基组函数形式如下,GTO只有1S，没有2S，3S，而P型函数也只有2Px等，没有3P，4P等，即GTO与原子轨道没有一一对应关系。 但GTO有很明显的优点： 容易分离变量 GTO乘积定理能将双中心积分化为单中心积分，类推三中心、四中心积分也可化作双中心，最后化作单中心积分来处理。使多电子体系中许多问题得以解决,STO和GTO的比较,STO和GTO差别在于指数函数e-r 与e-r2,

21、GTO是指数2次方。因此GTO的特点是： 在近核处太平滑,即在核处没有尖点,在r比较大时, 衰减太快,这两部分都不符合原子轨道的真实行为,Gauss型轨道计算起来快速、方便。但它与电子运动的真实图象相差太远。 Slater型轨道会增加分子积分计算的难度，但它更接近电子运动的真实图象。 为了弥补这些缺陷，可用GTO的线性组合逼近STO的办法。 这种基组，以STO作为基函数，而每个STO用若干个GTO来逼近。解出的分子轨道仍然用STO的线性组合来表示，而GTO不作为原子轨道，只作为中间数学工具,从头算耗时最多的是双电子积分，若基函数个数为N，则 单电子积分数为P=N(N+1)/2， 双电子积分数为

22、: Q=P(P+1)/2=N(N+1)(N2+N+2)/8N4/8 即双电子积分数约与基函数四次方成正比。 对一个单STO，用3个GTO拟合，GTO要做20个积分，但一个STO双电子积分耗时约为GTO的102-103倍，因此总的来说计算GTO积分仍然比STO积分快得多。 目前的从头算程序大多采用GTO,基组(Basis Set,基组就是对某个元素所用基函数的集合 基组的类型 STO-NG 价层分裂基组 极化基组 弥散基组 赝势基组,STO-NG基组,STO-NG，即用N个GTO拟合一个STO，再用一个STO拟合一个原子轨道。也称为极小基组，或单zeta基组。 N=3-6。STO-NG基组的特点

23、： 从1s轨道到3d，4f轨道，每个轨道均用同样数目的GTO拟合。 应用最广泛的是STO-3G，它是极小基的首选基组，常用于几何构型优化，寻找反应过渡态的试探计算，较大体系的计算。比如C的原子轨道为1s, 2s, 2px, 2py, 2pz，先用5个STO拟合5个原子轨道，再用15个GTO拟合5个STO,极小基组包含了描述一个原子的最少的基函数 例如 第一周期原子H、He：1s 第二周期原子Li、Be、CNe： 1s, 2s, 2px, 2py, 2pz STO-3G是最小基组 一般量化程序都备有STO-NG基组，只要写出相关的关键词，程序会自动调用这些基组,价层分裂基组（N-31G和3-21

24、G,化学反应主要在价轨道进行，而STO-NG则对内层轨道与价轨道一视同仁。为了提高计算精度，需用多个GTO去拟合一个原子轨道，对于较大的体系，计算机内存与计算时间将难以承受，故有人想出两者兼顾的方法，对内层轨道用一个STO拟合原子轨道，价轨道则用两个STO拟合，则价层分为内轨（I）与外轨（O）两个，称为分裂价基。即： 只有价层的基组加倍 内层轨道用最小基组, 因为它们原子化合成分子的过程中几乎不变,双分裂价基组,也叫双基组：对于价键轨道都用两个STO函数来进行描述。常用的有：3-21G，4-31G和6-31G几种可选。另外还有三分裂价基如6-311G,6-31G即内层轨道用6个GTO拟合一个S

25、TO，一个STO拟合一个原子轨道，价层则分为内外轨，即双的I层用3个GTO拟合一个STO（I），O层1个GTO拟合一个STO（O），然后一个STO（I）和一个STO（O）拟合一个价层原子轨道。 例如：P原子用6-31G基组： 1s，2s，2p轨道各用6个GTO拟合一个STO，再用一个STO拟合一个原子轨道，3s，3p则分别用3个GTO拟合3s(I)，3p(I)，1个GTO拟合3s(O)，3p(O)，然后3s(I)与3s(O)的STO拟合3s轨道， 3p(I)与3p(O)的STO拟合3p轨道,一些元素的分裂价基函数 H、He：1s（I）, 1s（O）也可表示为： 1s, 1s Li至Ne：1s,

26、 2s, 2s, 2Px, 2Py, 2Pz, 2Px, 2Py, 2Pz Ga到Kr： 1s 2s, 2Px, 2Py, 2Pz, 3s，3Px，3Py，3Pz，3dz2，3dx2-y2， 3dxy，3dyz，3dxz 4s，4Px，4Py，4Pz，4s，4Px，4Py，4Pz 分裂价基比STO-NG仅价轨道增加一倍，但计算精度提高不少，在计算较大体系的反应生成热、解离能时，常选用分裂价基组,极化基组,分裂价基组允许轨道改变其大小(size)，但不能改变形状(shape)。极化基组则取消了这样的限制，增加了角动量 角动量更大的函数添加到基组中时, 轨道在角度方向就有了更大的可变性 比如对氢原

27、子加上p函数, 对C原子加上d函数 如果不加极化函数, 即使使用大基组, Hartree-Fock计算得到的NH3的结构是平面型的 没有极化函数时, 环丙烷的张力太大,6-31G(d) (也写作6-31G*)基组，来源于6-31G基组，并在其基础上，对于非氢原子增加了d轨道函数； 6-31G(d,p) (也写作6-31G*) ，在6-31G(d) 的基础上，在氢原子轨道中加入了p轨道,Effect of Polarization Functions,现以乙烯（C2H4）为例，说明各种基组原子轨道数目： STO-3G-C： 1s, 2s, 2Px, 2Py, 2Pz, H: 1s；共25+14=

28、14个原子轨道、42个初始的GTO； 3-21G- C： 2s, 2Px, 2Py, 2Pz为双层, H： 1s为双层，共2（1+24）+24=26个STO，42个初始GTO； 6-31G* - C：2s, 2Px, 2Py, 2Pz层外加d极化函数，共2（1+2 4+5）+2 4=36个STO，70个初始GTO； 6-31G* - 在6-31G* 的基础上在H1s外层加p极化函数，共48个STO，82个初始GTO。 由此可见，随着基组的增大，原子轨道数逐步增加，而计算机容量与时间消耗则与基函数的四次方成正比。所以大基组一般只能用于小体系,弥散基组,把指数很小的函数添加到基组中，弥散函数是s和

29、p轨道函数的大号的版本，它们允许轨道占据更大的空间 需要的情形: 阴离子, 电负性高的原子F, Cl, 计算电子亲和能和气相酸度, 弱作用(氢键, 范德华作用) 6-31+G(d)基组表示的是6-31G(d)基组在非氢原子（重原子）上加上弥散s和p函数, 6-31+G(d)基组表示对氢原子添加弥散s函数且对非氢原子添加弥散s和p函数,高角动量基组是在分裂价层基组基础上增加多个极化基函数 6-31G(2d)是在6-31G基础上增加两个d轨道的函数。 6-311+G(3df,3pd)则增加了更多的极化函数，包括三个分裂的价键基组，在重原子和氢原子上加的弥散函数，在重原子上加的三个d函数和一个f函数，在氢原子上加的三个p函数和一个d函数。 6-311+G(3df,2df,p)基组在第二周期以及以上都采用三个d函数和一个f函数的极化,而对于第一周期采用两个d函数和一个f函数的极化,高角动量基组,赝势, 有效势基组（Pseudopotential, PP，or Effective Core Potential，ECP,第三周期以上的原子的基组很难处理，由于存在非常大的核 在化学作用中, 内层轨道没有大的变化 价层轨道同时受到核和内层电子的作用 可以构造一个赝势, 来代替来自