高斯体会纵谈.docx

高斯体会纵谈.docx

《高斯体会纵谈.docx》由会员分享，可在线阅读，更多相关《高斯体会纵谈.docx（9页珍藏版）》请在冰豆网上搜索。

高斯体会纵谈

通过一段时刻对GAUSSIAN的利用，发觉不仅有助于加深对量子化学的了解，而且能够巩固光谱学方面的理论知识。

至少在研究光谱学理论方面，GAUSSIAN是确实一个相当不错的东西。

在利用中，也碰到很多问题，而且在GAUSSIAN的帮忙文件中没有涉及。

这些问题有的是靠自己解决的，有的是和他人探讨的，还有的是从网上或文献中查到的。

  受本人研究方向所限，文中涉及的大多数都是处置分子激发态方面的问题。

若是在分子光谱学（包括红外、拉曼、可见-紫外）方面的知识有欠缺的话，建议先看一看哈里斯和伯特卢西写的《对称性与光谱学》（胡玉才，戴寰译，高等教育出版社，1988），这是一本浅显易懂的速成书籍。

  本人水平有限，对有些问题的明白得可能有误，欢迎批评指正。

******************目录******************

（一）一些杂七杂八的小问题速度问题

（二）关于基组的设置选择（三）用GAUSSIAN处置激发态§3-3.含时方式（TD）输入方式§3-4.CIS，TDHF，TDDFT方式激发能比较势能曲线（四）求分子的光谱常数参数Te，ωe，Re[附]Gaussian98的补充1.谐振强度f2.用TDDFT做几何结构优化3.自概念赝势（模型势，有效核势）

（一）一些杂七杂八的小问题  G98W默许的内存是128兆，我以前一直用64兆，运算机是P233，因此奇慢非常。

我后来发觉，在计算之前把所有的杀毒程序都关闭（我的运算机上装了三个），能够节省时刻20%。

后来又加了一条64兆内存，时刻就减少到原先的一半。

若是把运算机换成P300，用128兆内存，只需要之前五分之一的时刻。

而且大的内存支持更大的基组。

因此运算机配置可马虎不得。

  输入的分子坐标（包括键长和键角）必需是实型，不能是整型，不然会报错。

例如，键长2.0埃，能够写成r=2.或r=2.0，不能写成r=2。

  关于含有重原子的分子来讲，光谱项的自旋、轨道角动量已经失去了意义，那个时候的光谱项应该用自旋-轨道耦合总角动量的各个分量Ω表示。

例如，一氯化铊分子的3Π电子态，在实验上观看到的是两个分得很开的Ω分量：

A0+态和B1态。

可是GAUSSIAN仅能处置含有从H原子到Cl原子的分子（因为不是重原子，因此没有必要）。

目前处置这种问题，一样是利用免费的MOLFDIR，GAMESS（USA），NWChem，DIRAC，COLUMBUS或商业的GAMESS（UK），MOLPRO，ADF等程序。

（二）关于基组的设置  基组有两种，一种是全电子基组，另一种是价电子基组。

价电子基组对内层的电子用包括了相对论效应的赝势（缩写PP，也称为模型势，有效核势，缩写ECP）进行近似。

以铟原子为例，能够输入下面的命令来观看铟原子3-21G基组的形式：

#3-21GGFPrintIn0,2In  输出结果为：

3-21G（6D,7F）S31.00N=1层，s轨道.1221454700E+05.6124760000E-01.1848913600E+04.3676754000E+00.4063683300E+03.6901359000E+00SP31.00N=2层，s和p轨道.5504422550E+03-.1127094330E+00.1523702990E+00.1197743540E+03.8344349530E-01.6096507530E+00.3866926950E+02.9696880360E+00.3970249500E+00SP31.00N=3层，s和p轨道.4702931320E+02-.2758954250E+00-.1408484750E+00.2249642350E+02.5977348150E-01.5290866940E+00.6697116970E+01.1082147540E+01.6620681110E+00D31.00N=3层，d轨道.1021735600E+03.1205559000E+00.2839463200E+02.4884976000E+00.8924804500E+01.5850190000E+00SP31.00N=4层，s和p轨道.6572360380E+01.4284830560E+00.1091304620E-01.2502157560E+01-.4633643610E+00.5036758870E+00.9420245940E+00-.8219679320E+00.5581808650E+00D31.00N=4层，d轨道.4535363700E+01.2508574000E+00.1537148100E+01.5693113000E+00.4994922600E+00.3840635000E+00  以上都是内层电子轨道，每一个轨道用三个高斯型函数拟合，也确实是3-21G中的3。

SP21.00N=5层，s和p轨道.1659704190E+00.1189893480E+01-.9903308880E+00  这是价电子轨道，每一个轨道用二个高斯型函数拟合，也确实是3-21G中的2。

SP11.00N=6层，s和p轨道.5433974090E-01.1000000000E+01.1000000000E+01  这是空轨道，每一个轨道用一个高斯型函数拟合，也确实是3-21G中的1。

  可见3-21G是全电子基组。

再看看铟原子LanL2DZ基组的结果：

LANL2DZ（5D,7F）.4915000000E+00-.4241868100E+01.3404000000E+00.4478482600E+01.7740000000E-01.1000000000E+01.9755000000E+00-.1226473000E+00.1550000000E+00.1041757100E+01.4740000000E-01.1000000000E+01  这些都是价电子轨道。

在ECP上加上了双zeta（DZ）基组。

因此LanL2DZ是价电子基组。

  在GAUSSIAN中，CEP也是很重要的价电子基组，而且98版中对元素的适用范围扩大到了86号元素Rn。

需要注意的是，超过了Ar元素后，每一种元素的CEP-4G，CEP-31G，CEP-121G是等价的（但和加了极化*的CEP仍是有区别的）。

SDD适用于除Fr和Ra之外的所有元素，可是结果的准确性并非太好。

  需要注意的是有些方式不需要输入基组，像各类半体会方式，G1，G2等等。

  《中国科学》B30（5），2000，419页的文章给出了一种选择基组的方案，以为基函数的选择应该依照体系中不同原子的性质及实际的化学环境：

1。

描述一样体系时可依照该原子在元素周期表的位置从左到右依次增加基函数用量。

2。

对核负电原子，利用更多的基函数和适当的极化函数或弥散函数。

对核正电原子，基函数用量可适当减少。

3。

对正常饱和共价键原子不需要增加极化或弥散函数，对氢键、弱彼此作用体系、官能团、零价或低价态金属原子等灵敏体系，要增加极化或弥散函数。

  如此就能够够用适中的基组和可经受的计算量的到相当靠得住的计算结果。

方案适用于HF，MPn，DFT等计算中。

  在GAUSSIAN中有关的关键字是ExtraBasis和Gen。

二者都是在routesection指定基组的基础上，再加上附加的基函数。

不同的是Gen能够自己概念基组。

事实上可并非那么简单。

下面是一个利用ExtraBasis的例子：

#HF/3-21GExtraBasisAgClmoleculespecificationCl06-311G（d,p）****  那个地址对分子中所有的原子都利用3-21G，可是对Cl原子利用更大的基组6-311G（d,p）。

若是计算进程中报错，一样是由于下面两点缘故造成的：

1。

并非是所有的基组都能够用（哪些基组适用于哪些原子需要去参考帮忙文件），即便是别离适用于不同原子的各个基组，放在一路也不必然就行，可能有匹配或收敛的问题，也有可能是BUG。

因此需要反复多试。

2。

基组设置的顺序问题。

仍是看上面的例子，我以为routesection一行中设置的基组，应当适用于该分子中所有的原子（一样来讲是小一些的，尽管在计算中并非对所有显现的原子都利用那个基组），而在分子坐标下面的基组只要对所概念的原子适用就能够够了（一样来讲是大一些的基组）。

比如，若是把上面的输入变成：

#HF/6-311G（d,p）ExtraBasisAgClmoleculespecificationAg03-21G****  尽管都是对Ag原子用3-21G基组，对Cl原子用6-311G（d,p）基组，可是可能会报错。

固然实际情形并非老是如此，换个分子，基组哪个在前哪个在后可能并无关系，可是注意基组顺序确实能够避免犯错。

这可能是一个BUG。

（三）用GAUSSIAN处置激发态  输入文件中的多重度，指的是具有这种多重度的最低能量的电子态，不必然确实是基态，只是一样输入文件都用基态。

激发态的多重度是在输入文件的routesection行操纵的，如：

singlets，triplets。

尽管在GAUSSIAN利用帮忙中没有明说，可是从给出的例子中能够体会出来。

另外，在HyperChem的帮忙文件中有类似的说明，尽管是针对HyperChem的，但我想对GAUSSIAN一样适用。

  GAUSSIAN能够用CIS，CASSCF和ZINDO方式计算激发态。

GAUSSIAN中的CIS实际是对分子激发态的零级组态相关处置。

有些从头计算程序（如GAMESS）包括计算激发态的一级组态相关（FOCI）和二级组态相关（SOCI），我不清楚是不是和CIS相对应。

从结果来看，FOCI和SOCI比CIS略有改善。

ZINDO方式适用范围窄，只能用于周期表前48种原子组成的分子，那个地址不谈。

98版还新增加了含时（TimeDepend）方式处置激发态（包括含时Hatree-Fock（TDHF）和含时密度泛函（TDDFT））。

  咱们一样计算的激发态称为Low-lyingExcitedStates，也确实是由原子价电子组成的分子电子激发态。

更高的电子态是Rydberg态。

从头算法和其他半体会方式对空轨道，专门是较高空轨道的计算结果不行，因此取得的Rydberg态的结果缺乏明确的物理意义。

§3-3.含时方式（TD）输入方式  帮忙文件对TDDFT这么重要的方式说得真是太简单了，让人去参考CIS的输入方式，连个例子都没有。

CIS的输入方式是：

#CIS（Triplets,NStates=10）/3-21GTest  那就用B3LYP泛函碰运气吧：

#TDB3LYP（Triplets,NStates=10）/3-21GTest  可是不行。

其实正确的输入方式是：

#B3LYP/3-21GTD（Triplets,NStates=10）Test  那个问题浪费了我很长时刻。

为了表示它的重要，单独分为一节。

§3-4.CIS，TDHF，TDDFT方式激发能比较  一样来讲，CIS和TDHF方式结果的准确度是比较差的，TDDFT的结果是最好的。

而各类TDDFT相较，B3LYP，B3P86，MPW1PW91是最好的（依具体的分子和利用的基组而不同）。

能够参考文献：

J.Chem.Phys.,111（7）,1999,2889Chem.Phys.Lett.,297（1-2）,1998,60J.Chem.Phys.,109（23）,1998,10180  这些文章都是计算C，H，O组成的有机物在基态平稳位置的垂直激发能。

我对含有重元素In的分子激发态进行计算，发觉若是选好基组，很多TDDFT的结果和实验值的差距乃至小于500个波数，而CIS的结果那么可能要差5000多个波数。

因此此刻利用TDDFT方式的比较多。

另外还有人用TDDFT方式计算DNA如此的大分子。

1。

用94版进行CIS计算，50-50选项不支持LanL2DZ和3-21G*基组。

98版那么可不能显现那个问题。

那个问题可能还跟计算的分子有关。

2。

同时计算单重激发态和三重激发态（也确实是利用50-50那个选项），有时候会显现致命错误：

某些电子态在某些位置的能量有不规那么起伏，而单独计算单重态或三重态（利用Singlets或Triplets选项）那么可不能显现那个问题。

我在计算中曾经遇见过两次。

下面那个例子是从网上看到的，我想这位德国的大哥也遇见了一样的问题。

  他对同一个任务执行了三种计算：

1）#Bp86/6-31G*td=（nstates=5,50-50）scf=directdensity2）#Bp86/6-31G*td=（triplets,nstates=1）scf=directdensity--link1--#Bp86/6-31G*td=（nstates=10）scf=directguess=readgeom=check--link1--#Bp86/6-31G*td=（triplets,nstates=10）scf=directguess=readgeom=check3）#Bp86/6-31G*td=（nstates=10）scf=direct#Bp86/6-31G*td=（triplets,nstates=10）scf=directdensity  第一种计算产生10个态（5S和5T），后两种产生20个。

关键字Density对结果没有阻碍。

结果是：

        1）    2）    3）  通过比较发觉，对前六个态，三种处置的结果相同，可是后面的态，利用了50-50的方式就和另外两种方式的结果不一样了。

  我个人以为，尽可能不要利用50-50那个选项。

推荐单重态、三重态分开计算，虽谈多花一倍时刻，可是结果保险。

  最笨的方式莫过于算完一个点后，更改输入文件的坐标，计算下一个点了。

如此费时费力。

在优化基态几何构型的时候，能够利用Scan关键字进行扫描，事实上Scan对激发态一样适用。

如此开上一晚上运算机，到早上所有激发态的势能曲线数据就全有了。

例如计算从R=0.6埃到1.2埃氢分子的10个单重激发态：

#RCIS（NStates=10）/6-311GScanH2ExcitedStates0,1H1H2H1r  需要注意三点：

1。

结果给的是垂直激发能，必需换算成相关于基态平稳位置的能量。

2。

扫面范围不要离基态平稳位置太远，一样从0.75Re到2Re就能够够了，太远了会报错。

3。

某些理论和某些基组的搭配可能不支持Scan，这时就只好手工操纵了。

  最后是数据处置，这一步是最枯燥的了。

数据处置需要用到Excel、Oringin、Matlab。

这些软件用不着钻研很深，明白几个最大体的命令就足够用了。

  默许的CIS计算求出最低的三个激发态。

以它为例，GAUSSIAN在初始的猜想当选择电子态数量两倍的向量（也确实是6个）进行迭代，直到有三个收敛。

关于大多数分子体系来讲，这种默许的向量数量足够了。

可是关于高对称性的分子，那么需要特殊的方案，不然可不能取得比较高的收敛态。

方式是增加NStates的数量，以增加初始向量的数量。

推荐NStates值为最大阿贝尔点群操作的数量，也确实是输出文件在Standardorientation之前的Symmetry部份（在输入文件中利用#操纵打印输出就能够够取得，#T不行）。

这是一个计算苯分子的例子：

StoichiometryC6H6FrameworkgroupD6H[3C2'（HC,CH）]Deg.offreedom2FullpointgroupD2HNOp8LargestconciseAbeliansubgroupD2NOp4Standardorientation------------------------------  因此在routesection中利用CIS（NStates=8）。

（四）求分子的光谱常数  第一必需注意，原子化能和离解能不是同一个概念。

帮忙文件和ExploringChemistrywithElectronicStructureMethods一书中的说明并非确切。

计算原子化能和离解能是不同的两个进程。

以三原子分子GaCl2为例，计算离解能的进程是：

GaCl2（X2A1）---->Cl（2P）+GaCl（X1Σ+）进程中的能量转变称为离解能。

很明显，离解能还能够细分为第一步离解能、第二步离解能......，上面的是计算第一步离解能，计算第二步离解能的进程是：

GaCl（X1Σ+）---->Cl（2P）+Ga（2P）而计算原子化能的进程是：

GaCl2（X2A1）---->Cl（2P）+2Ga（2P）可见仍是有区别的。

只是关于双原子分子来讲，这两个概念是相同的。

  第二个需要注意的是，离解能一样用符号D0或De表示（前者不包括零点能，后者包括），可是在分子光谱中，还有一组D0和De。

为幸免混淆，那个地址需要作专门说明。

  第一组D0和De，也确实是离解能，数值比较大，数量级一样在0.1到1eV。

那个参数多出此刻真空-紫外、振动光谱中。

有一个公式：

      ωe^2ωeχe=-------        （^是指数符号）        4De其中的De确实是离解能。

  另一组D0和De称为离心畸变常数，数值很小，作为转动光谱中的微扰项。

如：

              h^3De=-----------------------------------    128（π^6）（μ^3）（c^3）（re^6）（ωe^2）其中的De确实是离心畸变常数。

  GAUSSIAN中有个求PH2分子的原子化能的例子，能够作为求离解能的参考。

也确实是文件e7-01。

结果是原子化能=E（P）（Hartree）+2E（H）（Hartree）-E（PH2）（Hartree）-ZPE（Hartree）      =（-341.25930）+（-0.50027）*2-（-342.50942）-（0.01322）（Hartree）  其中的零点能ZPE通过了矫正      =148.3Kcal/mol和实验值144.7Kcal/mol的差距相当小了。

  电离能缩写为I.P.。

在GAUSSIAN中有两种求I.P.的方式。

法一：

见文件e7-03，那个例子求PH2分子的电离能。

结果为：

    E      ZPEPH2+-342.144160.01347（ZPE通过矫正）PH2  -342.509420.01322（这是上面求原子化能的结果）I.P.=（E（PH2+）+ZPE（PH2+））-（E（PH2）+ZPE（PH2））    =（-342.14416+0.01347）-（-342.50942+0.01322）（Hartree）    =0.36551（Hartree）    =9.95（eV）法二：

OVGF方式。

帮忙文件里有详细的利用说明。

§4-3.求激发态的光谱学参数Te，ωe，Re  有两种方式。

1。

CIS能够用结构优化的方式。

Starl的文章里说得很详细了，那个地址再也不细讲。

仅仅需要补充的是，如何强制保留设置后的波函：

  在利用Guess=Alter设置波函以后，程序总会回到设置前的波函。

那个时候需要利用SCF关键字中的Symm选项，强制程序保留设置以后的波函。

利用方式是在routesection部份添加SCF=Symm。

2。

势能曲线扫描法。

  TDDFT方式不能计算激发态的波函，因此不能用来研究激发态的平稳构型和分子属性。

上面的结构优化方式用在TDDFT中就不适合了。

可是能够用势能曲线扫描的方式。

可是，多原子分子能量是势能曲面，操纵多个输入座标变量，扫描整个曲面是困难的，因此那个方式用于双原子分子。

如何扫描激发态的势能曲线见前面的§3-6节。

  有了势能曲线以后，就能够够找到势能最低点Re（为了精准，这一区域的点要取密一些），Re位置的能量相关于基态平稳位置能量的差值确实是电子态高度Te。

势阱的深度确实是离解能De。

  这些都简单，关键是求振动参数ωe。

求ωe也有两种方式。

一种是用Morse势函数拟合的方式。

Morse势函数形式为：

V（R）=Te+De{1-EXP[-α（R-Re）]}^2（对基态Te是0）Te、De、Re是已知的，拟合出了α后，用公式α=ωe（（2（π^2）c*μ）/（h*De））^0.5（高斯单位制）就能够够取得ωe。

  这种方式太麻烦，咱们不如此做。

实际只要把扫描出来的Morse势函数对R求二阶导数，在Re处的导数值f2和ωe的关系有：

    4（π^2）（ωe^2）c*μf2=--------------------      h*（10^18）  其中c、h、μ是国际单位制，ωe的单位是波数cm^-1。

另外在扫描的势函数中，De的单位是波数cm^-1，R和Re的单位是埃。

不然上面的f2需要乘上个系数。

  求f2用Origin软件就能够够做。

注意要取Re位置的f2。

有了f2，就能够够算出ωe，省去了计算De和拟合α的步骤，幸免引入更大的误差。

  若是用GAUSSIAN计算比较轻的分子，一样都能给出正确的激发态光谱项对称性，可是对含有重原子的分子，可能得不到对称性。

那个时候就需要依照轨道跃迁的的信息取得对称性。

  下面是用CIS方式LanL2DZ基组计算InCl分子的例子，输出文件给出InCl的电子轨道信息：

Occupied（SG）（SG）（PI）（PI）（SG）Virtual（PI）（PI）（SG）（PI）（PI）（SG）（PI）（PI）（SG）（SG）（SG）  那个地址结果只给出前4个单重激发态：

CISwavefunctionsymmetrycouldnotbedetermined.ExcitedState1:

Singlet-?

Sym4.9226eV251.87nm5->6.435875->7.54909Thisstateforoptimizationand/orsecond-ordercorrection.CopyingtheCisinglesdensityforthisstateasthe1-particleRhoCIdensity.CISwavefunctionsymmetrycouldnotbedetermined.ExcitedState2:

Singlet-?

Sym4.9226eV251.87nm5->6.54909ExcitedState3:

Singlet-SG6.7609eV183.38nm5->8.57818CISwavefunctionsymmetrycouldnotbedetermined.ExcitedState4:

Singlet-?

Sym6.8241eV181.68nm3->7.49690  由于In原子是重原子，因此大多数电子态并无给出对称性（第三个激发态除外）。

咱们发觉第一激发态是从第五个轨道跃迁到第六、七个轨道，而第六、七两个轨道是简并的（PI、DELT轨道都是二重简并的）。

也确实是说第一激发态的电子组态是：

InCl:

（内层电子）（1σ）2（2σ）2（1π）4（3σ）1（2π）1  那个组态耦合出来的单重态光谱项有一种（固然还有一个三重的）：

1Π，因此确实是1