求助]如何从下面的Gaussian输出文件中找出轨道系数及轨道能!!(新手多谢),请帮忙标出来!!!
The electronic state is 1-A1.
Alpha occ. eigenvalues -- -20.58265 -11.33946 -1.39265 -0.87259 -0.69715 Alpha occ. eigenvalues -- -0.63950 -0.52294 -0.44073
Alpha virt. eigenvalues -- 0.13573 0.24842 0.33338 0.37329 0.73660 Alpha virt. eigenvalues -- 0.80783 0.84685 0.94689 1.10445 1.10700 Alpha virt. eigenvalues -- 1.13937 1.27145 1.33529 1.62050 1.78192 Alpha virt. eigenvalues -- 1.79416 1.99239 2.18347 2.23684 2.45514 Alpha virt. eigenvalues -- 2.64513 2.87165 2.97616 3.27576 4.09792 Alpha virt. eigenvalues -- 4.47637 Molecular Orbital Coefficients
1 2 3 4 5 (A1)--O (A1)--O (A1)--O (A1)--O (B2)--O EIGENVALUES -- -20.58265 -11.33946 -1.39265 -0.87259 -0.69715
1 1 C 1S 0.00000 0.99566 -0.11060 -0.16262 0.00000 2 2S 0.00047 0.02675 0.20981 0.33995 0.00000 3 2PX 0.00000 0.00000 0.00000 0.00000 0.00000 4 2PY 0.00000 0.00000 0.00000 0.00000 0.42017 5 2PZ -0.00007 0.00066 0.17259 -0.18451 0.00000 6 3S -0.00024 -0.00743 0.08051 0.31309 0.00000 7 3PX 0.00000 0.00000 0.00000 0.00000 0.00000 8 3PY 0.00000 0.00000 0.00000 0.00000 0.15760 9 3PZ -0.00048 0.00135 -0.01160 -0.07970 0.00000 10 4XX -0.00002 -0.00272 -0.01628 -0.01333 0.00000 11 4YY -0.00006 -0.00202 -0.01365 0.03019 0.00000 12 4ZZ -0.00074 -0.00123 0.03302 -0.00166 0.00000 13 4XY 0.00000 0.00000 0.00000 0.00000 0.00000 14 4XZ 0.00000 0.00000 0.00000 0.00000 0.00000 15 4YZ 0.00000 0.00000 0.00000 0.00000 -0.01394 16 2 O 1S 0.99472 -0.00038 -0.19672 0.08889 0.00000 17 2S 0.02094 0.00025 0.44184 -0.20351 0.00000 18 2PX 0.00000 0.00000 0.00000 0.00000 0.00000 19 2PY 0.00000 0.00000 0.00000 0.00000 0.32122 20 2PZ -0.00153 -0.00029 -0.13537 -0.14216 0.00000 21 3S 0.00436 -0.00058 0.37895 -0.27048 0.00000 22 3PX 0.00000 0.00000 0.00000 0.00000 0.00000 23 3PY 0.00000 0.00000 0.00000 0.00000 0.17976
24 3PZ 0.00006 0.00108 -0.04718 -0.06799 0.00000 25 4XX -0.00418 0.00015 -0.00022 -0.00041 0.00000 26 4YY -0.00383 -0.00011 -0.00073 -0.00413 0.00000 27 4ZZ -0.00356 -0.00019 0.01969 0.00906 0.00000 28 4XY 0.00000 0.00000 0.00000 0.00000 0.00000 29 4XZ 0.00000 0.00000 0.00000 0.00000 0.00000 30 4YZ 0.00000 0.00000 0.00000 0.00000 -0.02339 31 3 H 1S -0.00002 -0.00020 0.03017 0.17902 0.19082 32 2S -0.00013 0.00210 -0.00537 0.06479 0.12026 33 4 H 1S -0.00002 -0.00020 0.03017 0.17902 -0.19082 34 2S -0.00013 0.00210 -0.00537 0.06479 -0.12026 6 7 8 9 10 (A1)--O (B1)--O (B2)--O (B1)--V (A1)--V EIGENVALUES -- -0.63950 -0.52294 -0.44073 0.13573 0.24842
1 1 C 1S 0.01942 0.00000 0.00000 0.00000 -0.12212 2 2S -0.06075 0.00000 0.00000 0.00000 0.14896 3 2PX 0.00000 0.32517 0.00000 0.40259 0.00000 4 2PY 0.00000 0.00000 -0.19811 0.00000 0.00000 5 2PZ -0.37597 0.00000 0.00000 0.00000 -0.21086 6 3S 0.03971 0.00000 0.00000 0.00000 1.98096 7 3PX 0.00000 0.21231 0.00000 0.71124 0.00000 8 3PY 0.00000 0.00000 -0.04477 0.00000 0.00000 9 3PZ -0.08856 0.00000 0.00000 0.00000 -0.74977 10 4XX 0.00549 0.00000 0.00000 0.00000 -0.00273 11 4YY 0.02734 0.00000 0.00000 0.00000 -0.01265 12 4ZZ -0.01933 0.00000 0.00000 0.00000 -0.00459 13 4XY 0.00000 0.00000 0.00000 0.00000 0.00000 14 4XZ 0.00000 0.03558 0.00000 -0.03288 0.00000
15 4YZ 0.00000 0.00000 0.06035 0.00000 0.00000Sample Text
相关回复：
作者: lixiaona158 发布日期: 2008-04-03
EIGENVALUES 后面的数字就是这个轨道对应的能量，但是它的单位是HF,一般使的时候需要换成电子福特，用这个系数乘27.2116就可以了。
作者: cuihang 发布日期: 2008-04-03 那个应该是基函数对应的能量，现在量子化学中的基组已经超出了原子轨道的范围。一般来说只有HOMO和LUMO是感兴趣的
作者: cuihang 发布日期: 2008-04-03
另外，注意#后面的关键字，如果选简略输出，会少很多信息
作者: liliracial 发布日期: 2008-04-03
以 1 (A1)--O分子轨道为例，O表示已占分子轨道，其下对应的EIGENVALUES -- -20.58265，就是该轨道的本征值，即轨道能量。在该值下面对应一长串的各原子各轨道对应的数值，就是各原子轨道组成该分子轨道的系数。原子轨道系数绝对值越大，说明该原子轨道对该分子轨道的贡献越大。对1 (A1)--O分子轨道而言，很明显主要是由O 1S 原子轨道组成，其轨道系数达到0.99472。
作者: zhaohongjie 发布日期: 2008-04-03 非常感谢大家的帮助!!!!
作者: zhaohongjie 发布日期: 2008-04-03
若是我还想算一下原子轨道贡献率呢?比如说算 (A1)--O分子轨道中碳C原子的,我要把它后面对应的1S、2S、。。。。。。直到4YZ对应的系数都加在一块吗？:)谢谢！
作者: gaochao85 发布日期: 2008-04-03
那你就做一个nbo 基本上你想要的轨道信息 一般都会有
作者: zhaohongjie 发布日期: 2008-04-03 不好意思,我是新手
作者: liliracial 发布日期: 2008-04-03 QUOTE:
Originally posted by zhaohongjie at 2008-4-3 09:37:
若是我还想算一下原子轨道贡献率呢?比如说算 (A1)--O分子轨道中碳C原子的,我要把它后面对应的1S、2S、。。。。。。直到4YZ对应的系数都加在一块吗？:)谢谢！
根据原子轨道线形组合成分子轨道，各个系数的平方之和应该等于1。要算贡献也应该是计算各系数的平方，而不是简单的加和。你可以看看结构化学或者量子化学这方面的资料。
nbo，是自然键轨道，Guassian03手册里面有这个关键字，可以去查查手册，看怎么用。
NBO布局分析
思路与Mulliken布居分析一样，但在以下几个方面有所不同 (1)计算原子上的电荷时所用的原子轨道不同： Mulliken分析用的是构建分子轨道的原子轨道，因而对基组有依赖性；而NBO分析用的是NAO [natural atomic orbital]原子轨道。 NAO是在对角化密度矩阵中一中心的块矩阵得到的[原话：diagonalization of one-center blocks of the density matrix]; Since the NAOs form an orthonormal set, completely spanning the space of (generally nonorthogonal ) basis orbitals, the natural populations are inherently positive, and sum correctly to the total number of electrons. 重叠区电荷的分割是通过权重而不是均分的方式分给两个原子. (2) 给出了分子中各原子的NHO (Natural Hybrid Orbitals), 易于为普通的化学工作者掌握和解决问题. 构建NHO (Natural Hybrid Orbitals) 的步骤如下: [1]find the density matrix P in a basis set of natural atomic orbitals and diagonalize eath atomic subblock PAA to find the lone-pair hybrid on that center. NBO程序默认布居数超过1.999e的轨道是core orbital,超过 1.90e的轨道是lone pair orbital. [2]for each pair of atoms A 、L, form the two-center density matrix P (AL) . NBO程序默认只有布居数超过1.90e的轨道才被考虑,如果 电子对的数目不够, NBO程序会继续搜索三中心的block, 直到体系的总电子数被满足 [3] symmetrically orthogonalize the hybrids found in step 1 and 2 to find the final natural hybrids orbitals.
(3) 可以定量的给出轨道之间的相互作用能.
Gaussian寻找过渡态经验小结（转贴） 相关搜索: 小结, Gaussian, 经验
寻找过渡态不是一件容易的事（对于我和大多数刚涉及量化的人来说），因此我希望通过写这个经验小结能对大家有些帮助。
1.首先遇到的问题是，用哪种方法来寻找过渡态？
GAUSSIAN提供的方法是QSTN和TSN方法。两种方法各有优点和缺点。QSTN方法特别QST3方法要求输入反应物，过渡态的猜测结构，产物这三者的结构。特别麻烦。但很管用，一般不会出现不收敛的情况。对于TSN（对应关键词为OPT=TS）方法，只要求输入过渡态的初始结构，但这个初始结构非常的关键，如果结构不好，则很容易出现不收敛的情况。所以我建议，如果是刚开始做过度态的话，用QSTN方法是好的选择，等有了“感觉”之后，再用TSN方法。
2. 怎么解决经常出现的错误？
在找过度态的时候，经常碰到的一些问题就是（1）不收敛，（2）有一个错误的本征值（错误信息为：there is a wrong sign eigenvalue in hessian matrix.....)，（3）和LINK9999错误导致退出。对于不收敛的情况，可以分为两类，比如提示信息里的CONVERGENCE FAILER 提醒收敛到了10（-5），而此时你设定的SCF循环次数也仅仅是64步，那么完全有希望通过加大SCF循环次数来达到收敛的目的。倘若只收敛到10（-3）或10（-2），此时加大循环次数可能就没用了。结果还是CONVERGE FAILER。 此时可采用SCF=QC，来达到强制收敛的目的 。因为SCF=QC（LINK508）的计算量比默认的L502要大，所以不到万不得以就不用它了。出现第二个错误可以直接用 关
键词OPT=NOEIGEN 来实现。LINK9999出错是因为已经走完了默认的步数，但还未完成。系统会自动跳出。出现这种情况大多数就是因为优化步数和SCF步数超过了默认值。可用OPT（MAXCYCLE=100）和SCF（MAXCYCLE=300）来改错。
3.怎么样控制过渡态的优化，使得过渡态不至于收敛到其他的分子结构中去？
我用GAUSS VIEW 可以解决这个问题 ，当刚开始运行GAUSSIAN时，你用GVIEW去打开输出文件时，你可以看到你的过渡态的初始输入结构，当一个循环过后（从上一个LINK502到下一个LINK502），你再打开输出文件，你就可以清晰地看到优化一步后分子的构型，这样就可以随时监控过度态分子的结构，倘若已经有收敛到其他分子构型的趋势时，你就可以把它给KILL了，而不至于需要等全部工作结束后，打开输出文件才知道已经不是想要的过渡态了。 如果收敛到其它的构型上去，可以考虑缩小OPT的步长.iop(1/8=2或3）即可。
4. 还需要加其它的关键词吗？
建议在OPT中加入CALCFC。这样可以加大找到过渡态的几率。本人深有体会！
5. 如何控制寻找过渡态的步长？
这个可以用IOP来实现。具体相应的IOP是iop(1/8)=2，此语句说明在找过度态的时候，以2A（A为基本单位长度）为单位来寻找过渡态。（转自量子化学网，作者：elizerbeth）
用Gaussian寻找过渡态（Transition State）
如何寻找transition state
Answer: A sample route section #gfinput iop(6/7=3)
#B3LYP/6-31G(d) Opt(TS,Noeigen)
In order to increase the efficiency of the saddle point search,we could cal culate
the force constants by adding \"CalcFC\" keyword. #gfinput iop(6/7=3)
#B3LYP/6-31G(d) Opt(TS,Noeigen,CalcFC)
We can also ask Gassian to automatically generate a gues structure for
the reaction
by using keyword \"QST3\" or \"QST2\" #gfinput iop(6/7=3)
#B3LYP/6-31G(d) Opt(QST3,Noeigen,CalcFC) A+B-->C Reactant //title section 0 1
structure of A+B A+B-->C Product 0 1
structure of C A+B-->C TS 0 1
guess structure for the TS
Note:the corresponding atoms need to appear in the same order within all th
e molecule specifications.
发信人: ghb (Never is a long time), 信区: Gaussian 标 题: Re: 如何寻找transition state?
发信站: BBS 大话西游站 (Sat Jan 5 11:10:02 2002)
找TS好像也不是那么简单
我试了一下用QST2优化一个光环化反应的TS 用PM3方法，竟然out文件有240M！ 而且link died at L9999 ft死了
仔细想想，其实也就是对反应物分子和产物分子的Redandunt coordinate
按照设定的path＝N做N等分，从这条路径找到一个近似的TS 如果用Z-matrix一般都不会收敛
然后再逐步调整，根据Hessian判断是否到达了真正的TS
我试过这样做，不知道是否有用，大家讨论讨论吧：
用WinMopac2.0也可以做半经验的IRC计算，但不一定要从TS开始
可以选定Reaction Coordinate，让他从反应分子变到产物分子 例如对某个键断裂反应，选键长为反应坐标
可以写成这样的形式：1.5,1.6,1.7,1.725,1.75,„,3.4,3.6,4.0,„ 在可能的TS附近可以写多一点
然后开始IRC，这样能够得到每个反应坐标处的一个初始构型 当然也有能量，键级等等
Winmopac做这个非常方便，只要吧要做IRC的设为-1，不优化的设为0，优化的设为1
利用这些初始构型，固定每个反应坐标，用比较高的基组进行优化 比如B3LYP/6-31G*
这样得到一系列构性，其基态能量最高者是否就是TS？ 当然许多情况下还应该计算激发态，比如光反应。
寻找transition state是量子化学中比较难的一个问题吧， 因为体系的波函数不如基态来的那么比较容易描述 ， 关于你具体的说的那个反应，我觉得做法上可能有些技巧可以利用， 如在做SQT2是，一般都不会直接输入反应物和产物的坐标， 而是会相应作些调整，比如将active bond伸长15%~20%左右， 这样做有利于加快计算速度，
还有，transition state对基组和理论的选取非常敏感 ， 所以要谨慎 。
Gaussian中关于IRC的输入和输出分析 %Chk=irc
#p hf/sto-3g irc(forward,calcfc) guess=read geom=allcheck
****************************************** %nproc=2 %Chk=irc2
#p hf/sto-3g irc(reverse,calcfc)
guess=read geom=allcheck scf(maxcyc=200) (results
SUMMARY OF REACTION PATH FOLLOWING: 与scan看结果一样是看OPti上的那个orient。
如果是单独作irc(forward)和irc(reverse)，取出最后opt结果上面的那个orient来看irc到底是个什么样的过程。第一个是过渡态的orient了。 （作者：elizerbeth ）
如何解决Gaussian算过渡态时Error termination in NtrErr Dear Everyone,
I would like to thank who get me the following usefull answer to my question.
I started an eigenvalue-following transition state optimization in g98 but the run failed giving the following error message:
------------------------------------------------------------------------
dumping /fiocom/, unit = 3 NFiles = 1 SizExt = 524288 WInBlk = 512
defal = T LstWrd = 66048 FType=2 FMxFil=10000 Number 0
Base 20480
End 66048
End1 66048
Wr Pntr 20480
Rd Pntr 20480
Error termination in NtrErr:
NtrErr Called from FileIO.
------------------------------------------------------------------------
I used the following route for TS search:
# PM3 test opt=(EF,TS,calcfc,ModRed) scf=direct freq
Can anybody tell me how to solve this problem?Thanks a lot! ---------------------
Hey,
you ask for 2 jobs. the first one is the TS optimization, the 2nd the FREQ calculation. Are you sure that the first job doesn't terminated correctly? It could be only the 2nd job that can not start. Hope this helps.
Ast. Pr. Xavier Assfeld Xavier.Assfeld@lctn.uhp-nancy.fr
Laboratoire de Chimie theorique (T) 33 3 83 68 43 74
Universite Henri Poincare (F) 33 3 83 68 43 71
F-54506 Nancy BP 239 http://www.lctn.uhp-nancy.fr ----------------------------------------------------------------------------- Hello Daniele,
maybe (I'm not sure) you have an old chk-File in your directory and your g98
jobs tries to use this, because its given in the input: %chk=... As you don't need this file to restart (your command line has no such option), you should remove it.
(If you are not sure, please show me the upper part of your input, containing the %chk line) Hope that helps. Bye
Elmar
Elmar Gerwalin , University of Kaiserslautern,Germany
Dept. of Theoretical Chemistry
and IT Service Team
elg@chemie.uni-kl.de 0631-205-2749
----------------------------------------------------------------------------------
Dr. Bellocchi,
The EF options are not consistent with redundant internal coordinates. If you want to use EF you will need to define a Z-matrix and use Z-matrix constraints. Alternatively you can leave out EF and use Modredundant.
Douglas J. Fox
Technical Support
Gaussian, Inc.
help@gaussian.com
---------------------------------------------------------------------------------- Hi, can you forward the answers to me too, cause I have problems with eigenvalue following.Thank you in advance. Best regards. Nurcan Senyurt [senyurtn@boun.edu.tr]
----------------------------------------------------------------------------------
Best Regards
Daniele Bellocchi, PhD Student Molecular Modeling Section
Istituto di Chimica e Tecnologia del Farmaco
University of Perugia
Italy （作者：elizerbeth ）
用Gaussian进行溶剂计算 （转贴）
相关搜索: 溶剂, Gaussian 用SCRF关键词，可以指定它给出的系列溶剂（默认为H2O） 可供选择的有19种溶剂，可以在SCRF关键词里面查到，也可以写成
SCRF=N (N=1~19,分别对应N所代表的溶剂）
也可以用自己设定的溶剂
但是要在输入它的介电常数和溶剂半径A0
介电常数可以查化学手册，A0可以用Volume关键词进行计算
计算方法可以选择PCM和Dipole两种模型
e.g.
# RHF/6-311++G(d,p) OPT=Vtight Freq SCRF(PCM,Solvent=4) Test
Gaussian03常见出错分析 最初级错误
1. 自旋多重度错误 2. 变量赋值为整数 3. 变量没有赋值
4. 键角小于等于0度，大于等于180度 5. 分子描述后面没有空行
6. 二面角判断错误，造成两个原子距离过近
7. 分子描述一行内两次参考同一原子,或参考原子共线 运行出错
1. 自洽场不收敛 SCF a. 修改坐标，使之合理 b. 改变初始猜 Guess c. 增加叠代次数 SCFCYC=N d. iop（5/13=1） 2. 分子对称性改变
a. 修改坐标，强制高对称性或放松对称性 b. 给出精确的、对称性确定的角度和二面角 c. 放松对称性判据 Symm=loose
d. 不做对称性检查 iop(2/16=1) 3. 无法写大的Scratch文件RWF a. 劈裂RWF文件 %rwf=loc1,size1,loc2,size2,……..,locN,-1
b. 改变计算方法 MP2=Direct可以少占硬盘空间 c. 限制最大硬盘 maxdisk=N GB
4. FOPT出错 原因是变量数与分子自由度数不相等。 可用POPT 或直接用OPT
5. 优化过渡态只能做一个STEP 原因是负本征数目不对 添加 iop（1/11）=1
6. 组态相互作用计算中相关能叠代次数不够，增加叠代次数 QCISD（Maxcyc=N） Default.Rou设置
• 在Scratch文件夹中的Default.Rou文件中设置G03程序运行的省缺参数： • -M- 200MW • -P- 4
• -#- MaxDisk=10GB
• -#- SCF=Conventional or Direct • -#- MP2=NoDirect or Direct • -#- OPTCYC=200 • -#- SCFCYC=200
• -#- IOPs 设置 如iop（2/16=1） Default.Rou设置中的冲突
• Default route: MaxDisk=2GB SCF=Direct MP2=Direct OPTCYC=200 SCFcyc=100 iop(2/16=1) iop(5/13=1) • ------------------ • # ccsd/6-31G** opt • ------------------
• L903/L905 and L906 can only do MP2.
问题在于，ＭP2=Direct ！ 去掉这个设置，CCSD的作业就能进行了。
因此，建议在Default设置中只设置，内存，最大硬盘，等项 振动分析的关键词 振动分析的选项
拉曼光谱 Raman
非谐性 Anharmonic 势能函数不是简谐的
内转动，受阻转子 Hinderedrotor 力常数过小的转动
振-转耦合 VibRot
振动圆二色性 VCD 用于光活性分子
同位素效应，温度，压力 ReadIsotopes 不影响频率数值
算法选项：
解析频率 Analytic 数值频率 Numeric 二阶导数数值 Enonly
FREQ的缺省设置是：有解析方式就计算解析频率，不计算Raman性质
MP3，MP4，CISD，CCSD，QCISD等方法没有解析频率，只能算数值频率。
一般各种理论方法下，计算的频率和实验测得频率之间有系统误差，需要用校正因子校正（Scale Factor） 各种方法计算频率的Scale因子
方法： 约化因子 约化因子 （频率） （ZPE）
HF/3-21G 0.9085 0.9409 HF/6-31G(d) 0.8929 0.9135 MP2(Full)/6-31G(d) 0.9427 0.9646 MP2(FC)/6-31G(d) 0.9434 0.9676 BLYP/6-31G(d) 0.9940 1.0119 B3LYP/6-31G(d) 0.9613 0.9804 SVWN/6-31G(d) 0.9833 1.0079 Freq=ReadIsotopes输入 # 。。。Freq=ReadIsotopes 。。。
Title 0 2 H O 1 0.98
298.15 1.5 0.8929 温度，压力 约化因子
H 2 用质量数是2（程序读取精确值）的同位素 O 18 用质量数是18（程序读取精确值）的同位素 空行
Freq=（ReadFC，ReadIsotopes）
同一方法基组下计算不同同位素的力常数矩阵（Hessian）相同，因此从前面频率计算的chechpoint文件里读取力常数，计算同位素效应的频率和热力学性质，可以节省许多时间。
注：同种元素的不同同位素，几何结构是相同的，电子结构也是相同的。
也可以用Freqcheck.exe工具做不同温度，压力，同位素的热力学分析。 数值频率问题
• 数值频率计算不需要很大内存和硬盘空间，但是精度较差，时间也长。
Freq=(Numeric,(Step=N))
Step用于指定微分步长，缺省为0.001*5 Angstrom Freq=Enonly 计算解析一阶导数。数值二阶导数
振动分析中的问题
• 振动分析原则上只对稳定点有意义
• 振动分析一定要在结构优化相同的方法和基组上进行。 • 理论计算的振动频率和试验结果有系统误差，需要乘以Scale因子校正
• 考虑相关能方法的解析频率计算需要大的硬盘和内存 Warning -- explicit consideration of 3 degrees of freedom as vibrations may cause significant error
• 如果出现这个警告， 一般来说是出现两个方面的问题： 1. 优化的结构不够精确，需要用更严格收敛判据优化。 2. 有3个振动模式实际上应该按内转动处理。 包括：自由转子模型，受阻转子模型。
内转动模式对热力学性质计算有影响，特别是低频模式 稳定点优化中出现小的虚频问题
1。优化的结构还不是局域极小值点，需要用更严格的收敛判据优化。OPT（tight or verytight）
2。对应转动模式的虚频，可以用稍微修改转动模式对应的二面角，使之合理。
3。用OPT=CalcAll来计算 结构优化和频率计算中的常见错误
1. 自洽场不收敛 SCF Convergence Failure a. 修改坐标，使之合理
b. 改变初始猜 Guess c. 增加叠代次数 SCFCYC=N
d.先用较松的SCF判据优化结束后再用正常SCF判据优化 e.不管不收敛，继续计算，iop（5/13=1） f. 尝试用SCF=qc算法
2. 分子对称性改变a. 修改坐标，强制高对称性或放松对称性b. 给出精确的、对称性确定的角度和二面角如对Td对称性，二面角值给为，109.47122 度c. 放松对称性判据 Symm=loose d. 不做对称性检查 NoSymm或iop(2/16=1)
3. FOPT出错 原因是变量数与分子自由度数不相等。 可用POPT 或直接用OPT
4. 超出最大优化圈数Step Exceed
• 检查分子结构是否合理，如合理可增加循环圈数，OPTCYC=N
• 读取Checkpoint文件最后几何，继续优化。
• 选取上次优化过程中最好的几何结构继续优化。 GEOM=(step=n)
5. 无法写大的Scratch文件RWF
a. 劈裂RWF文件 %rwf=loc1,size1,loc2,size2,„„..,locN,-1 b. 改变计算方法
c. 限制最大硬盘 maxdisk=N GB 6. 组态相互作用计算中相关能叠代次数不够，增加叠代次数 QCISD(Maxcyc=N) CCSD
(Maxcyc=N)等。
本章讨论应用电子结构理论研究化学反应.我们将从电子密度开始,然后
回顾第四章中有关反应势垒的讨论,再讨论反应研究中的更复杂的技术,
最后,通过对相应反应的计算,来研究未知体系的反应热. 本章将引入两种新的计算方法 * 势能面 * 反应路径分析 8.1 预测电子密度
将电子密度或静电势可视化是研究一个分子体系的反应性的重要的第一步.
例8.1 文件 e8_01a, e8_01b 取代苯的电子密度
在有机化学中,亲电芳香取代反应的定位效应是已经被深入研究的课题.
在这里,我们采用电子密度对这一现象进行研究.
已经知道氯苯和硝基苯的硝化是基于同样的反应机理:苯环首先受NO2+的
攻击,产生各种异构体的阳离子异构体.当硝化完成后,产物分布如下. 邻位 间位 对位
氯硝基苯 29% 1% 70% 二硝基苯 7% 88% 1%
我们在这里检验间位和对位异构体的中间体.
分子采用B3LYP/6-31G(d)进行优化,电子密度在HF/6-31G(d)等级计算.将
电子密度按照平行苯环平面的方向切片,得到不同厚度位置的电子密度图.
间位的氯硝基苯和对位的二硝基苯的电子密度分布显示,其保留了有较大
共振范围的电子结构,相反,另两个构型的电子密度分布显示其电子分布相
对局域化,并且向苯环外的方向集中.
通过电子密度的图形,可以定性的理解电子密度和反应性的关系,在得到结论
之前,检查这个体积的电子密度是必要的.关于这方面的进一步资料可以参见
Gaussian出版的白皮书Visualizing Results from Gaussian. 8.2 计算反应焓变
例8.2 文件e8_02 水解反应
现在分析水解反应 H+ + H2O --> H3O+
目的是计算标准反应焓变dH298.其计算方法可以表示为 dH298 = dE298 + d(PV)
dE298 = dEe0 + d(dEe)298 + dEv0 + d(dEv)298 + dEr298 + dEt298 其中
dEe0: 0K时产物与反应物的能量差;
d(dEe)298: 0K到298K电子能量的变化.对于这个反应,这一项可以忽略;
dEv0: 0K时反应物和产物的零点能之差; d(dEv)298: 0K到298K振动能量的变化; dEr298: 产物和反应物的旋转能之差; dEt298: 产物和反应物的平动能之差; d(PV): 由于有一摩尔分子消失,PV=-RT.
dEe0由单点能得到,本例采用的计算方法是B3LYP/6-311+G(2df,2p).其他的各项
都要考虑内能校正,通过频率分析得到.这样,所要做的工作就是进行优化然后进
行频率分析得到所需数值.采用B3LYP/6-31G(d)就能够得到足够精确的结果.
这里注意我们不用计算H+,由于没有电子,它的电子能量显然是0;由于只有一个
原子,其振动,转动能显然也是零,这样,其只有平动能,其值为 1.5RT = 0.889kcal.mol.(详见统计热力学).
最终计算得到dH298=-163.3kcal.mol.实验值为-165.3+-1.8kcal/mol. 两者符合的相当好.
8.3 研究势能面
对于势能面的研究对反应路径分析来讲,可能产生出人意料的好的结果.本节中,
我们通过实例研究势能面的应用方法. 考虑丙烯基正离子的旋转异构体的变化, H1 H1 | |
H2a C H2b H2a C H3b \\ / \\ / <---> \\ / \\ / C C C C | | | |
H3a H3b H3a H2b (I) (I')
曾经认为两个异构体之间的变化是通过一个具有Cs对称性的过渡态完成的,在该
构型中,H2b-C-H3b组成的平面垂直与碳原子平面.采用HF/6-311++G(d,p)能够找
到这样的过渡态,但是进一步的采用MP2和QCISD以及同样基组的研究却没有得到
过渡态,而得到了极小值!这个新的具有Cs对称性的结构中,H1迁移到了端位的碳
原子上.这个新结构的能量比势能面中平衡结构的能量高10kcal/mol.
这样就有了另一条反应路线: 中间碳上的氢迁移到端位的碳原子上; 新形成的甲基旋转;
旋转后的甲基上的一个氢原子迁移回中间碳原子.
在这个例子中,应用了IRC计算来确定过渡态的确是连接产物与反应物的.
本章后面将对这一方法进行讨论.
HF方法的研究得到了假的过渡态,原因是,由于HF方法本身的限制,其计算的亚甲
基旋转的势垒要低于氢原子迁移的势垒. 8.4 势能面扫描
势能面扫描可以研究一个区域内的势能面.一般的扫描都是由一系列的在不同结
构上的单点能计算组成的.当进行势能面扫描时,要设置分子结构的变量,设置需
要变化的结构的范围和步长.
在Gaussian中,势能面扫描是自动进行的,下面是一个进行势能面扫描的算例.
#T UMP4/6-311+G(d,p) Scan Test CH PES Scan 0 2 C
H 1 R R 0.5 40 0.04
该算例要求一个对于CH的势能面扫描,所用的关键词是scan,变量的设置格式是: 名称 初始值 [点数 步长]
当只有一个参数时,变量在整个扫描中是不变的,当三个参数都设定时,变量将在 一定范围内变化.
当有多个变量时,所有的可能构型都要计算.
所有等级的计算结果都在输出文件中列出,比如进行的MP2的势能面扫描也将列出 HF方法的结果.
根据得到的扫描结果,可以得到所要的势能面,通过它,可能得到极小值的可能位
置.势能面扫描过程中不进行几何优化.
8.5 反应路径分析
在第四章中我们提到,得到一个过渡态机构不能说明它就是连接产物和反应物的
结构.分析其是不是所需过渡态的一个方法是分析虚频的简正振动状态.有时,对
振动的分析也不能够确定.本节讨论更为精确的方法.
IRC方法检验过渡态分子的趋势.计算从过渡态开始,根据能量降低的方向来寻找
极小值,就是说,寻找过渡态所连接的两个极小值.
反应路径是连接反应物与产物的,但是连接反应物和产物可以有不止一条路径,
通过不同的过渡态连接,通过IRC计算,我们可以寻找真正的反应路径,也就是能量 最低的反应路径.
反应路径计算可以确认得到的过渡态就是连接反应物和产物的过渡态,一旦确认,
还可以计算活化能(注意零点能校正). 运行IRC
在Gaussian中,运行IRC的关键词是IRC.需要注意的是,IRC计算是从过渡态开始的,
在两个反应方向上各进行固定步骤的计算(默认是6步). IRC计算的方法是这样的: * 优化过渡态
* 进行频率分析,确认所得到的是过渡态,计算零点能,生成进行IRC
计算的力矩阵
* 运行IRC,在鞍点的能量下降方向,寻找极小值.一般的,需要增大寻找的次数,
从而尽可能的接近极小值.方法是设置MaxPoints. 确定反应势垒,一般还要进行更多的工作, * 对过渡态的高等级的能量计算
* 对反应物和产物进行优化和频率分析,得到零点能,进行高等级能量计算
8.6 势能面研究实例
我们现在用Gaussian的反应路径分析来研究甲醛的势能面.这个势能面上有很多
极小值,包括甲醛,羟基卡宾,以及H2和CO.每一组之间都可以组合成不同的反应物
产物对.这里研究两个反应 H2CO <--> CO + H2 H2CO <--> HCOH 甲醛的解离
我们要确定反应过渡态的结构,预测反应的活化能.为此,我们需要以下信息:
* 甲醛,氢分子,一氧化碳分子的考虑零点能的能量. * 过渡态的几何构型和零点能校正的能量. 计算在HF/6-31G(d)水平进行,结果如下
SCF能量 零点能 总能量
H2 -1.12683 0.00968 -1.11716 CO -112.73788 0.00508 -112.7280 H2 + CO -113.84996 H2CO -113.86633 0.02668 -113.83966 计算过渡态的能量,方法是
* 过渡态几何构型优化,计算SCF能量 * 频率分析,计算零点能 * IRC计算,确认过渡态 优化过渡态
例8.3 文件 e8_03 CH2O --> H2 + CO IRC
首先考虑氧原子垂直于CHH平面的构型,同时增大OCH夹角.计算中设置
Opt=(TS,CalcFC). CalcFC一般对于过渡态的优化是有帮助的. 得到的该点几何构型与猜测的结构接近,SCF能量-113.69352 频率分析
频率分析表明其有一个虚频,零点能0.01774(校正后),总能量-113.68578 IRC计算
IRC计算需要优化好的过渡态和相应的力矩阵,得到的方法是
* 从临时文件中获得(IRC=RCFC),或 * 在IRC计算的初始进行计算(IRC=CalcFC)
IRC计算在输出文件末尾对计算进行总结,列出能量和优化的变量的值.第一个
值和最后一个值是整条路径的起点和终点.
在起点上,我们得到了一个类似甲醛分子的结构,可以认定该反应路线是通向
甲醛的,在终点上,得到了一个C-H键伸长的结构,C-O键略微缩短,也表明这条
反应路线是通向解离分子的. 计算活化能
IRC计算确认了所得到的就是我们所要的过渡态,下面就可以计算活化能了. 能量 活化能(kcal/mol)
过渡态 -113.67578
反应物 -113.83966 102.8(正向) 产物 -113.84996 109.3(反向) 计算表明两个反应方向的势垒相似.
注意IRC得到的产物的能量不一定等于两个单独的体系的和,因为当IRC计算得到
分子配合物的极小值,与两个分离体系的能量和有些差别
1,2氢迁移反应
现在用同样的步骤研究第二个反应.
反式氢基卡宾的包含零点能的总能量是-113.75709,计算方法是RHF/6-31G(d). 寻找过渡态
猜测过渡态在碳原子上的一个氢原子象氧原子方向迁移,处于同时与碳原子和
氧原子作用的位置,对其进行的频率分析表明其位一阶鞍点,包含零点能的总能 量为-113.67941 反应路径分析
IRC分析得到的两个结构,一个类似于HCOH,一个类似于H2CO,说明该结构为该反 应的过渡态. 活化能预测
计算得到的活化能为100.6(正向)和48.7(反向)kcal/mol IRC的注意事项
虽然实际的反应结构如极小点,极大点,鞍点等在势能面上存在几何的和数学的
意义,但不能简单推广到物理的和化学的意义.实际的分子是有动能的,这样它就
可以不遵循反应路径.当然,计算得到的结果提供了最经济的反应途
径.
8.7 等构反应(Isodesmic Reactions)
等构反应是指反应前后各种键的数量不变的反应,比如乙醛与乙烷生成丙酮和甲
烷的反应,反应前后,各种价键的数量都没有变化.
由于这一特点,对这样体系的研究可以得到相当精确的结果. 例8.5 文件e8_05 等构反应的反应焓变化 本例计算上面提到的反应的反应焓变,步骤如下: * 在HF/6-31G(d)水平优化结构 * 进行频率分析计算零点能
* 在B3LYP/6-311+G(3df,2p)水平计算单点能
结果得到反应焓变-9.95kcal/mol,实验值-9.9+-0.3kcal/mol. 例8.6 文件e8_06 通过等构反应确定二氧化碳分子生成焓
本例中利用等构反应CO2 + CH4 --> 2H2CO 来确定二氧化碳分子的生成焓. 原理如下:
dHcalc = 2 E0(H2CO) - (E0(CO2) + E0(CH4)) dHf(CO2) = -(dHcalc - dHf(CH4) + 2dHf(H2CO))
甲烷和甲醛的生成焓实验值分别为-16.0和-25.0kcal/mol(0K) 计算方法同上例,所得到的反应焓变为60.64kcal/mol,这样计算得到的 二氧化碳的生成焓为-94.64kcal/mol.实验值为-93.96kcal/mol.
等构反应的局限
等构反应对于反应体系和生成焓的研究非常重要,但其也有明显的缺点,如下
* 必须依赖良好的实验数据,实验数据不准确,得到的生成焓自然也不可信
* 该技术不能推广到活化势垒方面
* 该技术不能用于实际上发生不了的等构反应 * 不同的等构反应可以得到不同的生成焓数值. 例8.7 文件e8_07 等构反应的局限
本例通过两个不同的等构反应计算乙烷的生成焓,研究理论 MP2/6-31G(d)//HF/6-31G(d) 反应一
丙烷 + 氢 --> 乙烷 + 甲烷 反应二
乙烷 + 氢 --> 2 甲烷
另外通过6个反应计算SiF4的生成焓,研究方法 MP2/6-31G(d,p)//HF/6-31G(d) SiH4 + F2 --> 2H2 + SiF4 (i) SiH4 + 4HF --> 4H2 + SiF4 (ii) SiF3H + HF --> SiF4 (iii) SiF2H2 + 2HF --> 2H2 + SiF4 (iv) SiFH3 + 3HF --> 3H2 + SiF4 (v)
SiH4 + 4F2 --> SiF4 + 4HF (vi) 结果如下
乙烷生成焓实验值 -20.0+-0.1 kcal/mol 反应一 -17.361 反应二 -22.258
SiF4 实验值 -386.0+-0.3 kcal/mol (i) -375.983 (ii) -366.588 (iii) -380.506 (iv) N/A (v) -376.112 (vi) -385.379
虽然得到的乙烷生成焓和实验值比较吻合,但两个结果之间竟然有高达5kcal/mol
的差距,对于这样的简单的烃的体系仍然有这样大的差异,在使用这一方法是就 不得不小心了.
对于SiF4,不同反应得到的数据差别也很明显,虽然也有很符合的结果,但最大的
差距竟然达到20kcal/mol.同时注意这些硅化合物的生成焓本身的可靠性不高.
8.8 练习
练习8.1 文件 8_01a~c 水和反应 计算两个反应的反应焓
Li+ + H2O --> H2OLi+ at 298.15K H2O + H2O --> (H2O)2 at 373K
两个反应的实验值分别是-34.0+-0.2kcal/mol和-3.6+-.5kcal/mol 计算方法是B3LYP/6-311+G(2df,2p)//B3LYP/6-31G(d)
由于本例中要计算298.15K和373K两个温度的水分子的热力学数据,所以,可以在
一个计算中完成,输入文件如下 %Chk=water
#T RHF/6-31G(d) Freq=ReadIso ... --Link1-- %Chk=water %Nosave
#T RHF/6-31G(d) Geom=Check Freq=(ReacFC, ReadIso) Guess=Read Test ...
373.0 1.0 0.9135 16
1 1
计算得到的第一个反应反应焓变为-34.0kcal/mol.注意对Li+不能进行频率分析,
其焓的校正为1.5RT(只有平动能)
对于第二个反应,注意由于DFT方法在处理弱作用体系上有些困难,要确认所得到
的分子是真正的极小值,而不是过渡态.本例中需要采用Opt=CalcAll来进行优化.
计算结果为反应焓-2.9kcal/mol 练习8.2 文件8_02a, 8_02b 键的解离
本例中通过势能面扫描研究键的断裂过程,研究的体系是CH键 CH UMP4/6-311+G(d,p) 键长范围 0.2-2.5A CH4 RQCISD(T)/6-311++G(d,p) 和 键长范围 .75-3.15A
UQCISD(T,E4T)/6-311++G(d,p)
UQCISD的E4T选项进行的是在MP4(SDTQ)水平进行MP4计算,而不是默认的MP4(SDQ).
计算中还需要设置的关键词有
IOP(2/16=1)表示在扫描中忽略对称性变化
Guess=(Always,Mix)表示将HOMO和LUMO轨道混合来消除自旋对
称性的影响,并且
在每一个点都重新计算新的猜测波函数. 下面是检查势能面扫描数据的方法
* 将HF, MP2, MP4, QCISD(T)不同水平计算的键能对键长作图 * 利用这些图来分析键解离: * 限制性和非限制性方法的比较 * 电子相关
* 三重态对QCISD水平的贡献
对于CH,UHF曲线在其他曲线之上,相对而言,HF的结果是最差的.MP2的曲线比
MP3和MP4(SDTQ)曲线要高些,但三者相差不多.
对于CH4,UHF得到的曲线远高于其他结果,MP2的结果也比其他结果要高些,而
其他各个计算方法得到的曲线基本相似.每一次提高微扰的次数,能量曲线就
向低的方向移动.
QCISD曲线与MP4曲线很接近,在对于整个曲线的描述上,QCISD(T)比MP4(SDTQ)的 结果要好很多.
限制性方法和非限制性方法得到的结果只在键分裂之后才显现出来,RHF方法
得到的曲线趋势就是错误的,因为其得到的结果是,当原子进一步远离
时,能量
继续上升,实际上当两个原子距离达到一定程度后,能量将基本保持不变.而对
于QCISD(T)方法,两者的差别不大,限制性方法得到的曲线在键解离后要略高于
非限制性方法的结果.这说明,在描述键解离的过程时,非限制性方法的结果要
比限制性方法好,特别是在采取低等级计算模型时. 练习8.3 文件 8_03 H2CO的势能面
我们以前研究过H2CO的两个反应的势能面,得到了五个稳定点:其中三个极小值,
甲醛,反式羟基卡宾以及一氧化碳和氢分子;两个过渡态分别连接甲醛和另外两个
产物.现在要做的,就是研究这两个产物之间的变化途径.反应分两步: 反式羟基卡宾 <--> 顺式羟基卡宾 <--> 一氧化碳 + 氢 在HF/6-31G(d)水平研究这一问题,研究步骤是: * 寻找过渡态
* 确定得到的稳定点是过渡态,计算零点能 * 确定这一过渡态连接的两个极小值
一个可能的连接顺式和反式结构的过渡态就是HCO平面与COH平面成90度的二面
角.我们采用Opt=QST3,分别给出两个异构体和这一猜测的构型.然后
进行IRC计
算,得到的是接近顺式和反式异构体的构型,说明这个过渡态结构就是所要找的 过渡态.
用同样的方法,可以得到有顺式异构体到一氧化碳和氢的反应的过渡态.
反应的活化能为 正向 反向
反式 <--> 顺式 25.8 20.5 顺式 <--> 解离 68.0 131.6
这样就得到了关于H2CO的新的势能面.
点还是鞍点，或者要得到相关的热力学性质，经常需要对优化后的几何结构进行振动分析。首先，原则上说，振动频率分析只对稳定结构有意义。这里所说的稳定结构包括是势能面上的局域极小点和鞍点。实际上，一个分子可以有很多的自由度，如果在所有自由度上分子都处在稳定平衡，就是稳定的分子。频率分析得结果是所有频率都是正的，表明这是一个局域的极小点。如果分子只在一个自由度上处于不稳定平衡位置，其他自由度上都处在稳定平衡位置，说明该结构是一阶鞍点。分子在稳定自由度方向上的振动才是真实的振动，在不稳定
自由度方向上的实际上是不会有振动的。不过我们可以对不稳定方向上的运动也按振动来做数学处理，会的到负的振动频率，我们称它为虚频。虚频的出现表明该结构为鞍点。
第二，Gaussian计算中，频率的计算一定要在和分子结构优化相同的方法，基组下进行，否则计算的结果是没有意义的。我们知道，任何理论水平下的计算，都是在一定的近似下进行的，不同的理论水平的近似程度是不同的。在一种理论水平A下优化的稳定结构Geom_A会和另一种理论水平B下优化的稳定结构Geom_B有差别，也就是说Geom_A不会是理论水平B下的稳定结构。根据前面我们所讨论的，在理论水平B下对一个不稳定的结构进行频率分析是没有意义的。
第三，频率计算中可以考虑同位素效应（Freq=ReadIsotopes）。在波恩-奥本海默近似下，对于同一种元素采用不同的同位素对几何优化和电子结构计算没有影响，频率计算所需的力常数矩阵（Hessian矩阵）也不会变化，变化的只是约化质量。容易理解，重的同位素会导致低的振动频率。实际上，原子序数大的元素的同位素效应非常不明显，一般只需考虑H原子的同位素效应。
第四，各种方法计算的频率和实验结果之间存在系统误差，需要乘以一个约化因子来进行校正（Scale=f）。一般来说，理论计算的频率值
会比实验结果大。下面是一些理论水平下的约化因子。注意频率和零点能的约化因子是可以不同的。更多水平下的约化因子需要查文献获得。
方法： 约化因子 约化因子 （频率） （ZPE） HF/3-21G 0.9085 0.9409 HF/6-31G(d) 0.8929 0.9135 MP2(Full)/6-31G(d) 0.9427 0.9646 MP2(FC)/6-31G(d) 0.9434 0.9676 BLYP/6-31G(d) 0.9940 1.0119 B3LYP/6-31G(d) 0.9613 0.9804 SVWN/6-31G(d) 0.9833 1.0079
第五，Gaussian的频率计算有时会遇到下面的警告：
Warning -- explicit consideration of 3 degrees of freedom as vibrations may cause significant error
这时一般有两种可能：一种可能是，优化的几何结构不够精确，还没有达到稳定点。对于这种情况，需要考虑用OPT=tight或OPT=VeryTight，结合Int=fine或者Int=VeryFine进行更加精确的优化。第二种可能是，在计算的振动模式中，有部分的低频模式对应于内转动模式，在热力学分析中应该按自由转子或受阻转子模型处理，如果按谐振子模型处理会造成较大的误差。采用受阻内转子模型，可以用
Freq=Hindered-Rotor计算或者自己手动计算这些振动模式对热力学性质的贡献。
第六，对于计算出的力常数较小的振动模式，其势能面形状可能不满足谐振子模型的要求，力常数函数需要考虑非谐性的成分，这时可以用Freq=Anharmonic来进行非谐性处理。
第七，在频率输出前，会有这样的输出:
Low frequencies --- -0.0008 0.0003 0.0013 40.6275 59.3808 66.4408 Low frequencies --- 1799.1892 3809.4604 3943.3536
上面行的六个模式实际对应于平动（前三个）和转动（后三个），理论上，这六个数都应该是零。如果是对过渡态进行的频率计算，第一行会先输出虚频，然后才给出平动转动的数值。一般来说，平动的三个数都是接近于零的。如果转动的三个数不接近于零（一般解析频率10个波数以内，数值频率50个波数以内），说明需要进行OPT=Tight或OPT=VeryTight计算。第二行的几个实频率应该和随后输出中的相应频率进行对比，如果基本一样，说明频率计算结果很好；如果差别较大，说明这些频率受到的平动和转动的污染较大。
第八，频率计算后程序会进行一步优化计算，正常会得到四个判据都是YES的结果，这说明优化的结构很好。但是有时也会遇到四
个判据不全是YES的情况，这时需要仔细分析，一般有三种情况： 1、Force和RMS Force都是收敛的，Displacement和RMS Displacement虽然是NO，但是都比较接近收敛判据；
2、如果Force和RMS Force都是收敛的，Displacement和RMS Displacement不收敛且数值远大于收敛判据； 3、Force和RMS Force不收敛。
对于第一种情况，我们可以不管它，仍然可以认为计算结果是可靠的。对于第二和第三种情况，一般说明优化过程中估算的Hessian是不准确的，优化的结构可能还没有达到稳定点，需要重新优化。如果反复优化仍然无法解决这个问题，建议在优化过程中使用OPT=CalcAll。
第九，优化一个局域极小点，收敛后从力常数本征值看应该是局域极小点（没有负的本征值），但是频率计算中会出现虚频和负的本征值。这种情况下虚频一般是由转动模式造成的，说明分子中两个基团之间的相互位置不是很合适，需要绕转动的键相对转动到一个合适的位置，重新优化。绕某个键转动两个基团，有时可以很方便地用修改二面角的方法实现：OPT=Modredundant结合分子描述后输入：* m n * [+=]value。其中m，n为连接两个基团的键的顶端原子。有时，用OPT=CalcAll也可以解决这个问题。
第十，Gaussian程序中，频率计算是和温度、压力无关的。因为Gaussian所考虑的是原子核在一定核构型和电子运动状态构成的势场内振动，一般来说，温度和压力的变化不会影响分子的结构和电子运动状态，所以振动的势能函数是与温度、压力无关的。Freq=Isotopes关键词要求输入的温度和压力影响的只是平动，转动和振动的平衡分布和它们对热力学性质的贡献
优化第一步：
确定分子构型，可以根据对分子的了解通过GVIEW和CHEM3D等软件来构建，但更多是通过实验数据来构建（如根据晶体软件获得高斯直角坐标输入文件，软件可在大话西游上下载，用GVIEW可生成Z-矩阵高斯输入文件），需要注意的是分子的原子的序号是由输入原子的顺序或构建原子的顺序决定来实现的，所以为实现对称性输入，一定要保证第一个输入的原子是对称中心，这样可以提高运算速度。我算的分子比
较大，一直未曾尝试过，希望作过这方面工作的朋友能补全它。以下是从本论坛，大话西游及宏剑公司上下载的帖子。 将键长相近的，如
B12 1.08589 B13 1.08581 B14 1.08544 键角相近的，如
A6 119.66589 A7 120.46585 A8 119.36016 二面角相近的 如
D10 -179.82816 D11 -179.71092 都改为一致，听说这样可以减少变量，提高计算效率，是吗？ 在第一步和在以后取某些键长键角相等，感觉是一样的。
只是在第一步就设为相等，除非有实验上的证据，不然就是纯粹的凭经验了。 在前面计算的基础上，如果你比较信赖前面的计算，那么设为相等，倒还有些依据。 但是，设为相等，总是冒些风险的。对于没有对称性的体系，应该是没有绝对的相等的。
或许可以这么试试： 先PM3，再B3LYP/6-31G.（其中的某些键长键角设为相等），再B3LYP/6-31G（放开人为设定的那些键长键角相等的约束）。
比如键长，键角，还有是否成键的问题，Gview看起来就是不精确，不过基本上没问题 ，要是限制它们也许就有很大的问题，能量上一般会有差异，有时还比较大
如果要减少优化参数，不是仅仅将相似的参数改为一致，而是要根据对称性，采用相同
的参数。例如对苯分子分子指定部分如下： C
C 1 B1
C 2 B2 1 A1
C 3 B3 2 A2 1 D1
C 4 B4 3 A3 2 D2
C 1 B5 2 A4 3 D3
H 1 B6 2 A5 3 D4
H 2 B7 1 A6 6 D5
H 3 B8 2 A7 1 D6
H 4 B9 3 A8 2 D7
H 5 B10 4 A9 3 D8
H 6 B11 1 A10 2 D9
B1 1.395160
B2 1.394712
B3 1.395427
B4 1.394825
B5 1.394829
B6 1.099610
B7 1.099655
B8 1.099680
B9 1.099680
B10 1.099761
B11 1.099604
A1 120.008632
A2 119.994165
A3 119.993992
A4 119.998457
A5 119.997223
A6 119.980770
A7 120.012795
A8 119.981142
A9 120.011343
A10 120.007997
D1 -0.056843
D2 0.034114
D3 0.032348
D4 -179.972926
D5 179.953248
D6 179.961852
D7 -179.996436
D8 -179.999514
D9 179.989175
参数很多，但是通过对称性原则，并且采用亚原子可以将参数减少为： X
X 1 B0
C 1 B1 2 A1
C 1 B1 2 A1 3 D1
C 1 B1 2 A1 4 D1
C 1 B1 2 A1 5 D1
C 1 B1 2 A1 6 D1
C 1 B1 2 A1 7 D1
H 1 B2 2 A1 8 D1
H 1 B2 2 A1 3 D1
H 1 B2 2 A1 4 D1
H 1 B2 2 A1 5 D1
H 1 B2 2 A1 6 D1
H 1 B2 2 A1 7 D1
B0 1.0
B1 1.2
B2 2.2
A1 90.0
D1 60.0
对于这两个工作，所用的时间为57s和36s，对称性为C01和D6H，明显后者要远远优于前者。
好，我真的很想知道这句话“你手工把坐标写好输入，输入为C2对称性”是什么意思？老是见到版面上有人说对输入施加对称性限制，就是不明白是怎么回事。您能详细解释一下吗？有例子就更好了。麻烦了谢谢！
比如，你要输入CO2， 那么你这么输入 C 0.0 0.0 0.0
O 0.0 0.0 1.203212313123 O 0.0 0.0 1.203212313123
写得很精确，那就是对称性输入了。这样g98可自动判定为D∞h 如果你输入的差不多，比如： C 0.0 0.0 0.0
O 0.0 0.0 1.203212313123 O 0.0 0.0 1.204212313123
精度不够，那么g98就只能判定为 C∞h了。
再给你举个例子，比如甲烷的Td点群的输入坐标，考虑对称性，把C放在直角坐标的原点，四个氢放在四个角的顶点。 C 0. 0. 0. H 1. 1. 1. H -1. -1. 1. H 1. -1. -1. H -1. -1 1.
这样出来的结果点群判断就会是正确的，当然这里的1我是随便取的一个值，只是来说明这个问题而矣。 谢谢二位指教，也就是说第一个输入的原子应该是对称中心，其它输入的原子应保持数值上的对应性。那对于用Gview或3D生成的分子构型，也可以根据这个原则将处于对称位置的原子的坐标改为一致了？如果分子很大（如我这次算的C14H14N4S20）有三个环，整体性的对称输入就有些麻烦，是不是可以先生成分子构型，然后利用对称性对它的笛卡儿坐标或矩阵坐标进行修改，将相近的坐标的数据改为一致？
\"第一个输入的原子应该是对称中心，其它输入的原子应保持数值上的对应性\"没错,做到这点就很方便了,\"对于用Gview或3D生成的分子构型，也可以根据这个原则将处于对称位置的原子的坐标改为一致\"那还不如用手动输入,大体系很难办了,如果体系太大,不用保持严格的对称性吧. 有机分子，原子多了以后构型就没有什么对称性了。
但对于有些系统，比如 C60, 你当然得按照它的实际的对称性来输入。 如果你希望你的计算按照某种对称性来做，那就一定要手工输入严格的 满足对称性的坐标.用软件生成的坐标，不是严格满足对称性，即使 用symm=loose,有时也不一定能判断出来。
各位大虾是用什么软件建模的呀，我是用SGI里sybyl软件进行构像搜索后取最低能量构像导入高斯程序的；有的是用chem3D直接建模，这样子有没有什么问题啊，应该怎么做，谢谢
通常也就是这些方法,建议还是认真的核实分子结构的是否合理,然后再试着增大迭代圈数,降低收敛判据.
用什么软件建模并不是问题，对于我们量化计算，如vliant所说，核实分子的合理构型才是最重要的，包 括对个原子间的成键分析。就算不同的软件，也不可能保证给出的构型合理，因此我们一般都提倡手动输入坐标，就是为了加强对构型的理解，只有在这个基础上使用软件建模才能尽可能保证不出错。一个合理的计算结果不能寄希望于软件。 我们都知道常用的分子坐标输入有： 内坐标 直角坐标 混合坐标
分子构型输入方法的不同会对计算过程产生影响么？或者说，对于哪些体系，用哪种方法更好优化？ “混合坐标”怎么讲？
我通常用内坐标，它能很好的控制对称性，且时间要短。
“混合坐标”是直角坐标和内坐标的混合使用，高斯手册上有讲这部分内容。有看过讲哪种体系使用哪种方法更好优化，好像这都是凭经验的结果。
在“免费资源”论坛中有本《计算化学》，其中第八章（67－72页）好像提到你所关心的问题，你可以下来看一下，或许有所帮助。
如果你使用缺省的优化方案且初始几何构型完全相同，那么，不论你使用何种坐标输入方式，优化结果应该是一样的。简言之，优化结果只取决于优化方案和初始几何构型，与坐标录入方式无甚关联。
于对称性限制的优化来说，显然以内坐标输入，且选用opt=z-matrix可以减少变量，使得优化效率提高。免费资源中《gauss入门》有介绍。
您冗余坐标是指有过剩的内坐标。冗余内坐标算法可以用于结构优化。
优化第二步：
对于大分子，一般要先预优化（如用PM3或AM1），有时是用部分优化，然后才用能达到实验精度方法和基组进行优化，就我所了解大多数情况下，6-31G（D）基组就差不多。
不好意思，我是新手，刚刚接触guassian，最近要算两个较大有机分子中原子的静电荷急相关能量和几何性质，分子中已经超过了80个原子（包括氢），用chemdraw画出图形，在chem3d中生成的input文件，可是，软件报错——2070~~
好象说是连接（last link）错误，怎么回事，哪里出错了，各位指点一下~ 怎么，先用pm3预优化就可以了？？
我的意思是，初试构象对计算量有这么大的影响么？？
我以前也遇到过相似情况，同一个分子，相同的计算条件，用gview画出的模型进行计算报2070的错，用chemdraw画的模型，chem3d生成的input文件进行计算就OK，完成可以计算完成。 我还以为你画分子式的时候出什么问题了呢~~
只想优化H原子，其他重原子不做优化。怎样进行分子坐标输入，下面的坐标输入方法对吗？就是仅H原子的参数设置成Vairable. %chk=sample #p HF/6-31G* opt
Title line 0 1
cu 7.577854 5.926989 4.935050 cu 9.091644 8.489386 4.644861 cu 10.605857 5.766733 4.988313 p 10.432217 3.600341 5.598077 p 12.551342 6.888289 4.775263 p 7.164457 9.584280 4.209578 p 11.125158 9.457819 4.661391 p 7.503791 3.715690 5.278502 p 5.670551 7.057133 4.652208 s 9.100712 6.979012 6.564150 s 9.088921 6.492596 3.113104 n 8.884097 3.048614 5.950712 n 5.827205 8.602518 4.003874 n 12.406788 8.511755 5.225239
H 1 B14 2 A13 7 D12 H 1 B15 2 A14 7 D13 H 1 B16 2 A15 7 D14 H 1 B17 2 A16 7 D15 H 2 B18 1 A17 3 D16 H 2 B19 1 A18 3 D17 H 2 B20 1 A19 3 D18 H 2 B21 1 A20 3 D19 H 3 B22 1 A21 2 D20 H 3 B23 1 A22 2 D21 H 3 B24 1 A23 2 D22 H 3 B25 1 A24 2 D23 H 1 B26 2 A25 7 D24 H 2 B27 1 A26 9 D25 B14 2.990841 B15 2.990841 B16 2.985833 B17 2.985833 B18 3.009872
B19 3.009872 B20 2.998051 B21 2.998051 B22 3.015712 B23 3.015712 B24 3.002042 B25 3.002042 B26 3.489189 B27 3.477113 A13 161.290094 A14 146.871331 A15 108.517371 A16 96.114237 A17 106.549047 A18 100.057030 A19 155.415250 A20 145.871488 A21 163.425916 A22 138.991013 A23 105.076579 A24 99.387878 A25 66.120736 A26 65.475681 D12 64.540783 D13 -53.098055 D14 -163.347541 D15 152.573174 D16 152.164746 D17 -162.671553 D18 -56.766147 D19 42.781046 D20 46.122424 D21 -46.534404 D22 150.210575 D23 -164.670929 D24 -78.866142 D25 109.837455
你如果要这么做的话，必须加opt＝Z-matrix,这样的话它就只优化后面的H，前面的直角坐标就不会优化了。你只用了opt的话，那就是全优化了 你也可以通过molredunt来实现部分优化 固定原子的方法有许多,最简单的是这样的： #rhf/sto-3g opt nosymm test
the two C atoms are frozen 0,1
8 0 -1.000000 0.000000 0.000000 6 -1 0.000000 1.000000 0.000000 6 -1 0.000000 -1.000000 0.000000 8 0 1.000000 0.000000 0.000000 我不知道你的那种方法行不行，但我这个肯定可以。 前面我都贴出来了：） 这可是我的心血噢。 ＝＝＝＝＝＝＝＝＝ opt=z-matrix的部分优化。 %Chk=reac-hf3
#p RHF/STO-3G opt(z-matrix,maxcyc=500) freq=noraman scf(maxcyc=200) reactant fopt 1 1
N -1 -1.75389 -0.48105 0.31424 H 0 x2 y2 z2 H 0 x3 y3 z3 H 0 x4 y4 z4
O -1 2.16296 0.18089 -0.15503 H 0 x6 y6 z6 H 0 x7 y7 z7 H 0 x8 y8 z8 x2=-2.2088 y2=0.39705 z2=0.46255 x3=-2.14589 y3=-0.9279 z3=-0.4899 x4=-1.88451 y4=-1.06386 z4=1.11629 x6=1.36256 y6=0.28882 z6=-0.67401 x7=2.32973 y7=0.98241 z7=0.34634
x8=2.90378 y8=0.01009 z8=-0.74115
*********************************************** opt=modredundant的部分优化。
分子坐标部分可以为具体的数值，而不用为参数形式。 %Chk=reac-hf3
#p RHF/STO-3G opt(modredundant,maxcyc=500) freq=noraman scf(maxcyc=200) reactant fopt 1 1
N -1 -1.75389 -0.48105 0.31424 H 0 H 0 H 0
O -1 2.16296 0.18089 -0.15503 H 0 H 0 H 0 1 5 F
感谢各位大侠的帮助！
看了之后有个不明白的地方，还请大侠不吝赐教！ 在opt=modredundant的部分优化时最后一行的 1 5 F 是什么意思？
就是将第一和第五号原子的距离固定
问一个菜鸟的问题：N 后面的-1和H后面的0是什么意思？是像手册上所说的把坐标增加所需要的值吗？ 见笑了！！
－1指的是固定，0指的是优化，前面helpme贴的输入方法也是这样的
部分优化输入格式之二 opt=z-matrix的部分优化。 %Chk=reac-hf3
#p RHF/STO-3G opt(z-matrix,maxcyc=500)
freq=noraman scf(maxcyc=200) reactant fopt 1 1
N -1 -1.75389 -0.48105 0.31424 H 0 x2 y2 z2 H 0 x3 y3 z3 H 0 x4 y4 z4
O -1 2.16296 0.18089 -0.15503 H 0 x6 y6 z6 H 0 x7 y7 z7 H 0 x8 y8 z8 x2=-2.2088 y2=0.39705 z2=0.46255 x3=-2.14589 y3=-0.9279 z3=-0.4899 x4=-1.88451 y4=-1.06386 z4=1.11629 x6=1.36256 y6=0.28882 z6=-0.67401 x7=2.32973 y7=0.98241 z7=0.34634 x8=2.90378 y8=0.01009 z8=-0.74115
*********************************************** opt=modredundant的部分优化。
分子坐标部分可以为具体的数值，而不用为参数形式。 %Chk=reac-hf3
#p RHF/STO-3G opt(modredundant,maxcyc=500) freq=noraman scf(maxcyc=200) reactant fopt 1 1
N -1 -1.75389 -0.48105 0.31424
H 0 H 0 H 0
O -1 2.16296 0.18089 -0.15503 H 0 H 0 H 0 1 5 F
冻结坐标的目的是什么了?我虽然会做但不太清楚Gaussian为什么会有这个功能.
还有最近做了一个体系.有一个键如果不冻结的话.就优化不出那个反应中间体.那个键会断.我这样冻结键优化出来的结构计算频率也无虚频.它能看作为一个minima吗?这个结构可以作为我反应机理的中间体吗?我用溶剂化优化结构可以在不冻结坐标的情况得到一个和这个结构很接近的结构. 冻结坐标有时是有用的，举几个例子：
1，如优化溶质分子和nH2O的supermolecule，周围的水分子，可能结构变化不大，这时就可以固定水的键长或键角，减少工作量。
2，如用cluster模拟大块固体，这时也许就要按照晶格参数固定原子键长和键角。 3，又如不固定NH3－H－H2O中N－O键长是得不到质子传递势垒。
4，还有在优化稳定构型或过渡态得时候也有时用到部分优化得到接近得结果，再放开优化，这样逐步得到你想要得结果。
冻结坐标对一些大的体系,超分子体系好象是很有用的.好象做生物分子的反应时氨基酸残基必须先固定.但我这个体系如果不固定,只要一放开,结构就跑了会优化到其它构型上去.就算stepsize设得很小也没有用.
优化第三步：
优化结果分析，一般要进行频率计算，无虚频，表示是极小值；有虚频，表示是过渡态。如果出现虚频，可根据虚频的简正模式进行调整以获得极小值，我还没有碰到这种情况，希望各位朋友能给些详细的资料。
如前信所说，键长键角有差别是很正常的，因为气相分子和晶体内分子的结构不可能完
全吻合，只可能近似相同。
另外，对不同元素的限制主要看基组是否能计算该类元素，例如你所采用的方法中用了
6－31G(D)基组，该基组支持H-Kr的元素，即可以计算Br原子但是不能够计算I原子。如果你所
用的基组不支持你所计算的某个原子可以采用扩展基组来进行计算，或者降低基组的精度。
我用Gaussian优化一个分子结构，由于目前实验上没有办法得到纯净的样品来实验验证，那么我怎么判断优化得到的结构是正确的？
一些主要的键长和键角的实验值或者别人计算的值是否和与你计算的分子类似？ 这个判断不好做。
所谓正确不正确，如果计算出来的频率没有虚频，就表示你找到了一个局域稳定点。要确定是不是全局稳定点很难，只有通过别的方式比如HOMO-LUMO能级差等信息看你算出来态是不是稳定。
也觉得高级的方法和大基组的结果更可信,谢谢各位高手的指点!!再问一个问题,如果出现虚频如何修改得到分子稳定的结构
用GV看一下，根椐产生虚频的振动模式调整结构
消去虚频有什么技巧？我有一个计算结果，虚频有三个，都只有－30~－20左右，如果结构改动过大，一不小心，虚频就会更大。头痛ing～～～
那你就把结构朝振动的方向改小一点试一下
确实据我了解，消除虚频的主要方法是改变构型。
其次在计算上还可以尝试：nosymm；加大循环次数；提高收敛度；iop(1/8=1)等。 其实我自己也是有这个问题，对于势能面很平，较小的虚频很难消除。 首先，
1在优化时采用Scf(tight)的选项，增加收敛的标准。再去计算频率。如果 还有虚频，参见下一步。
2. 对称性的影响，很多情况下的虚频是由于分子本身的对称性造成的。这样， 在优化时，如果必要，要将对称性降低，还有，输入文件有时是用内坐标。 建议如果有虚频的话，将内坐标改成直角坐标优化。
3.如果上述方法还有虚频，看一下虚频，找到强度较大的，将在频率中 产生的原子的振动坐标加到相应的输入文件中。这样，重新计算。直到 虚频没有。
4 实际上，如果分子柔性较大，很难找到最低点，这是电子结构计算的问题， 这种情况下，需要动力学的 东西，用构象搜寻的办法解决。 如：模拟退火，最陡下降法，淬火法等。将得到的能量最低的构象做 一般的电子结构计算，这样，应当没有问题。
不要讲你还没有得到最稳定的结构，那么，是你的分子有问题， 要么计算错了，要么就是游离出现代计算的范畴 zixia上大虾的回复:
将内坐标改成直角坐标优化”这一方法没有多大作用，我试过，还是一样的结果 另外如果是根据振动频率的模式来调整分子结构的话，还要注意调整完后的分子结构
对应的能量是不是比未调整之前有所降低，如果没有，说明还是不行。
alwens兄所讲的第三点，经过我的猜测和与alwens的探讨，终于知道应该怎么做： 分子结构输入用直角坐标，将这些产生虚频的振动坐标直接与原始坐标相加。 alwens兄进一步提到： 这个方法有时候不太好。
不一定将坐标全部加上。有时候可以选择加1/2等。 主要的目的是去掉由分子对称所形成的虚频。
象上面的那个计算，完全可以再利用scf(tight)的选项消掉1－2个虚频。 强度很小。
又多学了一招，自我感觉很幸福中„„^_^
怎么把内坐标改成直角坐标啊 看你的输出文件
Full point group CS NOp 2 Largest Abelian subgroup CS NOp 2 Largest concise Abelian subgroup CS NOp 2 Standard orientation:
--------------------------------------------------------------------- Center Atomic Atomic Coordinates (Angstroms) Number Number Type X Y Z --------------------------------------------------------------------- 1 46 0 -0.242506 1.028743 2.204311 2 46 0 -0.476137 2.619863 0.000000 3 46 0 -0.242506 1.028743 -2.204311 4 46 0 -0.242506 -1.618028 -1.382208 5 46 0 -0.242506 -1.618028 1.382208 6 46 0 -1.725240 0.349678 0.000000
这一部分就是，如果没有，看 Input orientation:也行。
在opt的选项中加入“tight”或“saddle＝0”，若是计算过度态，应是“saddle＝1”。但是，也不一定能消除虚频！！最好，将结构重新做一下，注意分子的构型，如CH3的三个H与相邻基团是交叉，还是重叠！！ 请教虚频
请问如何知道虚频的结构是什么，是不是频率计算负值都可以认为是虚频，那我用gaussianview 看有三个负值，接下来我想做过渡态，用IRC，那我应该怎么弄呢？谢谢指导 你的输出文件中这部分就是看是否有虚频的
1 2 3 A\" A' A\"
Frequencies -- -65.9781 56.7053 59.7893 Red. masses -- 105.9032 105.9032 105.9032 Frc consts -- 0.2716 0.2006 0.2231 IR Inten -- 4.2436 0.0151 1.4107 Raman Activ -- 0.0000 0.0000 0.0000 Depolar -- 0.0000 0.0000 0.0000 Atom AN X Y Z X Y Z X Y Z 1 46 -0.49 0.13 -0.15 -0.41 -0.16 0.08 -0.10 0.15 0.22 2 46 0.00 0.00 0.16 0.64 -0.01 0.00 0.00 0.00 -0.18 3 46 0.49 -0.13 -0.15 -0.41 -0.16 -0.08 0.10 -0.15 0.22 4 46 -0.22 0.05 0.08 0.05 -0.06 0.09 -0.26 -0.19 0.00 5 46 0.22 -0.05 0.08 0.05 -0.06 -0.09 0.26 0.19 0.00 6 46 0.00 0.00 0.44 -0.06 0.23 0.00 0.00 0.00 -0.31 7 46 0.00 0.00 -0.23 -0.11 0.01 0.00 0.00 0.00 0.52 8 46 0.00 0.00 -0.24 0.23 0.20 0.00 0.00 0.00 -0.48
Frequencies 那一行数值为负，即有虚频。有三个负值，说明找的不是过渡态。但是如果另外两个虚频很小，而另外一个虚频很大，而且它的振动模式对应于你所计算的反应模式（比如反就前后某一键的拉长，缩短等），在计算所用的方法充许校正范围内，也还可以认为是找到过渡态。但是一般要求能够把另外两个虚频去掉。方法是要据振动模式，调整结构，再优化。
般认为，过渡态仅有一个虚频，中间体没有虚频！！建议仔细阅读《Exploring Chemistry with Electronic Structure Methos》和《G98W程序使用入门》，后者稍微简单一点，厦大的免费资源里有！！
有一个虚频是过渡态的必要条件；有几个虚频肯定不行了，有时候有一个虚频也不一定就找准了，做IRC成功就可靠了。 另：复合物也会有小的虚频，只要不大于100，也可认为是一极小值。
没错,有些体系可能存在several or many local minima,那么找到的具有一个虚频的过渡态可能连接两个loccal minima,而非连接反应物和产物,因此必须做IRC来判断.
关于虚频的物理意义
自恨理论功底浅薄，所以特地在此向诸位高手求教：
在寻找过渡态计算中opt=ts, 需要计算频率，如果有一个虚频就说明优化所得的构型是个 过渡态。有时在做部分限制优化opt=minimal时，也能得到一个虚频（一般发生在frozen 部分的化学键上）。我想问一下，虚频的数值大小本身具有什么物理意义。后一种情况是否也可以认定为过渡态。 我的理解：
虚频就是在这个振动方向上力常数是负的
过渡态之所以有一个虚频，是因为它是一个鞍点，鞍点上有一个负的梯度的方向 分子沿这个方向振动时将转化为反应物或产物。就象从山上掉下来，受到的力可以 认为是负的。
但并不是有一个虚频就是过渡态，还要在这个虚频的振动方向上分别指向反应物和 产物
有一个虚频并不一定是你想要的过渡态，你要看虚频的振动方向是不是对应着反应的方向。 从概念上说，虚频就是化学反应方向上的简正坐标的二次偏导小于零。
另外，根据我的经验，这个虚频的负值应该大一些，如果这个值太小，那也不足信。不知对不对，呵呵！
优化输出信息：
这说明优化失败了，没有收敛，你可以尝试：
1。降低收敛标准（算一般的有机分子不要紧，如果是过渡金属，那就不要了） 2。增大循环次数 3。换基组和算法 4。改变初始猜测
是有机分子阿，呵呵，用的是ＡＭ１，还是太大吗ｈｓｉ
这种不收敛的情况可能是你最初用软件建模的时候构型不合理，仔细考虑一下你的分子的构型是否合理，是否还存在其它的构型，单纯地降低收敛标准意义不大，不过你可以先尝试一下增大循环次数，最多scf(maxcyc=300),再不收敛的话就不是循环的问题，是构型不对。
S能不能说说，π轨道和σ轨道的组成上有什么差别？这方面我比较欠缺。 ure. You can directly analyze the components of the molecule orbitals from the population in the *.log file. 我在版面上说过了，当初的讨论还真是热烈：）你可以翻看一下
我这里简述几句：
π轨道主要是由npx，npy，nd－1，nd+1组成
σ轨道主要是由ns，npz（键轴方向），nd0
还有，nd＋2，nd－2主要是组成δ轨道。
opt 是可以restart，前面的中间文件是少不了的。 1, %chk=... # .. opt
Stopped 2, %chk=...
#.. opt(restart) guess=read geom=allcheck -----------------
对弈freq的restart，是要采用数值计算方法才可以restart。 （以下我省略chk和方法的输入，自己加上） 1, freq(numerical)
2, freq(restart,numerical) guess=read
OPT=restart不能进行，why
因为死机OPT被中断，我采用OPT=restart继续进行，但是出现以下情况： Restoring state from the checkpoint file \"Mn7O14V11O2-PCOPT1-2.chk\". FileIO operation on non-existent file.
FileIO: IOper= 2 IFilNo(1)= -997 Len= 120004 IPos= 0 Q= 8710320
dumping /fiocom/, unit = 1 NFiles = 16 SizExt = 524288 WInBlk = 512 defal = T LstWrd = 667136 FType=2 FMxFil=10000
Number 0 0 0 0 0 508 522 598 Base 144896 149504 142848 22016 444416 143360 143872 172032 End 146944 172032 143360 142336 667136 143375 144572 172034 End1 146944 172032 143360 142336 667136 143872 144896 172544 Wr Pntr 143360 149504 142848 22016 199168 143360 143872 172032 Rd Pntr 143360 149504 142948 142020 199168 143360 143872 172032 Number 634 991 992 994 995 996 998 999 Base 199168 147456 146944 20480 142336 21504 20992 172544 End 444168 149097 146949 20510 142346 21554 21042 198796 End1 444416 149504 147456 20992 142848 22016 21504 199168 Wr Pntr 199168 147456 146944 20480 142336 21504 20992 172544 Rd Pntr 199168 149097 146949 20510 142346 21554 21042 172545
dumping /fiocom/, unit = 2 NFiles = 5 SizExt = 0 WInBlk = 512
defal = F LstWrd = 267264 FType=2 FMxFil=10000 Number 0 508 522 634 998 Base 266245 20480 20545 21245 20495 End 267264 20495 21245 266245 20545 End1 267264 20495 21245 266245 20545 Wr Pntr 266245 20480 20545 21245 20495 Rd Pntr 266245 20480 20545 21245 20495
dumping /fiocom/, unit = 3 NFiles = 1 SizExt = 524288 WInBlk = 512 defal = T LstWrd = 66048 FType=2 FMxFil=10000 Number 0 Base 20480 End 66048 End1 66048 Wr Pntr 20480 Rd Pntr 20480
Error termination in NtrErr: NtrErr Called from FileIO.
为什么我的file IO会不存在？是什么参数设漏了吗？请各位大虾支持，谢谢！
在你前一次的计算中，第一单点的计算未完成，所以，checkpoint文档中尚无优化过程的信息，无法restart.
过渡态优化过程中遇到一个这样的问题
我在过渡态优化过程中用这样的方法 # hf/6-31g(d) opt=(ts,calcfc) freq当运行到L999时Link died请问可能是什么原因？是不是我的过渡态输入有问题？另外，寻找过渡态有哪些技巧，请指教！
你的计算在设定的步数里面没有找到过渡态。 这有两个可能：一是初始输入的构型不好。 二是循环次数不够 (可以maxcycle=大一点的数)。 具体还是等高手来指点吧。
对于决定收敛结构的标准，做了一些很小但很重要的变动。当力比截断值小两个数量级时（即，极限值的1/100），结构即被认为是收敛的，即便是位移比截断值还大。这个变动有利于大分子的优化，因为它们的
最小值周围可能有非常平缓的势能曲面。 (G98manual p106)所以你的opt有2个yes
频率计算完成后还有一次优化，可能这次就优化到3个了
！做频率分析的时候只需要用Freq就可以，不需要加Raman选项，对吗？ freq=noraman不作拉曼强度计算,可以提高10%-20%的计算速度.
算一个分子的电离能，里面有I，CL 在优化起正离子的时候老是没有没有结果 体系的能量趋势不固定，有时大有时小，
MAXCYCLE=80后出现2个YES，然后DIE LINK 999 但是COPY 。OUT最后的构型再去优化后 结构4个全是NO
现在还在算，但是还是NO
怎么办？ 要OPT=CALCALL吗？有什么用？ 另外，还有一个题外的问题 这样写对吗？
OPT（MAXCYCLE=80，STEPSIZE=10）可以吗？
增大优化次数，optcyc=N; N给大一些。 opt=(maxcycle=N)也一样
前面的老帖子已经讲到这些，好象在免费资源里的Gaussian98 tips. 一部分引用如下：
First, whenever you encounter SCF convergence failures you need to look at the progress of the SCF by using #P in place of # which will print out the status at each iteration. It is useful to know if the SCF is converging slowly, oscillating, taking large jumps etc. Slow convergence can be a sign that a higher order converger like SCF=QC would be helpful. Oscillations can be a sign that a small HOMO/LUMO gap is causing a flip in the occupied orbitals and SCF=VShift=n is useful. Large jumps can indicate that a poor initial guess or a change in the electronic state from the initial guess has occured and normal convergence will eventually achieved but after a larger number of cycles and a better initial guess would help.
Second, whenever possible it is helpful to start with a smaller basis set, i.e. STO-3G, LANL1MB, etc for single points. This is cheaper and so you can often generate a better initial guess which can be used with GUESS=READ and the larger basis set and without the more expensive options. In your case the SCF converges initially quite smoothly and then bounces around without making progress. This makes SCF=QC a likely approach to clean up this problem, at least at the initial point.
You might want to look at adding IOP(5/13=1) which allows the SCF to print out the MO's even in the event that convergence was not reached. This should only be used with single points, i.e. remove OPT.
As to optimization constraints. To freeze atomic coordinates you should turn on NoSymm. This method of freezing structures causes problems for redundant internals which are being looked into but NoSymm fixes the most serious ones.
下面是血红素的化学合成的模型，由于分子有点大，又加上有个Fe，所以优化非常困难， 请求斑竹支 高招：
（1）这是一个高度对称的分子，如何利用加虚原子，提高计算速度？
（2）优化Fe，我采用了LANL2DZ，请问有没有更好的基组？怎么手工加基组？ （3）由于变量太多，您认为我通过怎样的办法，减少变量来提高计算速度？ 非常感谢！！！做出文章了，我一定挂你的大名，呵呵。 S
O 1 B1
Fe 2 B2 1 A1
N 3 B3 2 A2 1 D1 N 3 B4 2 A3 1 D2 N 3 B5 2 A4 1 D3 N 3 B6 2 A5 1 D4 C 4 B7 3 A6 2 D5 C 8 B8 4 A7 3 D6 C 9 B9 8 A8 4 D7 C 4 B10 3 A9 2 D8 C 5 B11 3 A10 2 D9 C 12 B12 5 A11 3 D10 C 13 B13 12 A12 5 D11 C 5 B14 3 A13 2 D12 C 6 B15 3 A14 2 D13 C 16 B16 6 A15 3 D14 C 17 B17 16 A16 6 D15 C 6 B18 3 A17 2 D16 C 7 B19 3 A18 2 D17
C 20 B20 7 A19 3 D18 C 21 B21 20 A20 7 D19 C 7 B22 3 A21 2 D20 C 23 B23 7 A22 3 D21 H 24 B24 23 A23 7 D22 C 11 B25 4 A24 3 D23 H 26 B26 11 A25 4 D24 C 16 B27 6 A26 3 D25 H 28 B28 16 A27 6 D26 C 20 B29 7 A28 3 D27 H 30 B30 20 A29 7 D28 H 1 B31 3 A30 2 D29 H 9 B32 8 A31 4 D30 H 10 B33 9 A32 8 D31 H 13 B34 12 A33 5 D32 H 14 B35 13 A34 12 D33 H 17 B36 16 A35 6 D34 H 18 B37 17 A36 16 D35 H 21 B38 20 A37 7 D36 H 22 B39 21 A38 20 D37 B1 4.21561119 B2 1.64053835 B3 2.01170682 B4 2.00223534 B5 2.01140630 B6 2.01948030 B7 1.38502525 B8 1.45027869 B9 1.36206816 B10 1.38792662 B11 1.38546948 B12 1.45027336 B13 1.36191526 B14 1.38377173 B15 1.38899517 B16 1.45046091 B17 1.36116892 B18 1.38527197 B19 1.39574494 B20 1.44813170 B21 1.36193008 B22 1.39696571 B23 1.37779525
B24 1.08374598 B25 1.38063627 B26 1.08377153 B27 1.37911143 B28 1.08372124 B29 1.37753366 B30 1.08361344 B31 1.37864549 B32 1.07802759 B33 1.07801256 B34 1.07797961 B35 1.07795952 B36 1.07798544 B37 1.07794319 B38 1.07803016 B39 1.07806743 A1 6.31961440 A2 94.07195309 A3 95.30512753 A4 94.38669042 A5 94.58713162 A6 127.10303272 A7 109.89556522 A8 107.10994453 A9 126.54307980 A10 126.81263629 A11 109.85724161 A12 107.06302128 A13 126.91593763 A14 126.68048694 A15 109.77898998 A16 107.13028054 A17 126.98352930 A18 126.67933314 A19 109.65642922 A20 107.43086078 A21 126.85650291 A22 125.36309437 A23 117.60829677 A24 125.56378588 A25 117.60334742 A26 125.52442366 A27 117.63863540 A28 125.46104554
A29 117.64536447 A30 98.48331793 A31 124.49018481 A32 128.37464062 A33 124.51087771 A34 128.41640662 A35 124.47708155 A36 128.39197501 A37 124.32863598 A38 128.24305473 D1 -108.51711460 D2 161.17665430 D3 70.83221148 D4 -18.82295908 D5 85.79136868 D6 -174.21655652 D7 0.09638172 D8 -87.07481057 D9 86.40576040 D10 -175.60372334 D11 -0.05103700 D12 -88.28171080 D13 88.94769251 D14 -174.64817315 D15 0.35121394 D16 -84.34232126 D17 84.15970765 D18 -172.91552513 D19 0.72935930 D20 -85.63442533 D21 -5.66240036 D22 179.49688145 D23 -6.19734957 D24 -179.54995066 D25 4.93855843 D26 179.77279426 D27 6.47288831 D28 -179.56429342 D29 -89.13791474 D30 -179.79660459 D31 -179.88612542 D32 -179.83142492 D33 -179.77466949 D34 -179.60732447
D35 -179.84766166 D36 -179.39904616 D37 179.72163849
我也在做大分子的计算，希望能交流交流。我的电子邮件：chsun@imr.ac.cn 我的经验是先做PM3优化，然后再用精度高的方法和基组。
你的前三个原子，不是一条直线吗？如果是直线，那对称性就可以是C4v了。 如果不是直线，那至少要做成cs对称性吧？把前三个分子所在的平面作为对称面。 你现在的对称性是没有——虽然看上去很对称。 Fe用赝势基组不知道会怎么样，你可以看看别人的文章。
做优化的时候，可以固定一些外层的原子，只让和中间的原子相连的动。
大约有3M
另外，Fe、O、S分子是在一条直线上么？如果不是，你可以将吡啶环上的参数先减少了，然后再加Fe、O、S，肯定能够减少不少参数，但是不知道由没有实际意义了
1）我谈一下我做过的体系加虚原子的经验，希望对你有帮助：一般情况下，虚原子加在分子的对称中心，或键轴Cn轴上，如要加第二个，可放在垂直于Cn轴的方面上。 （2）你可以到基组库去寻找基组，可能要比单单的Lanl2Dz效果好
（3）你的体系高度对称，肯定有不少的CH键长键角是一样的（或者差别很小），一律设为一个参数，让它同步优化。
PS：你的输入文件的参数也太多了点
我看了一下你的结构，你说是高度对称？我想问你一下这是你手建的吗？从你的坐标看好象没有对称性。32号H怎么位置这么怪！我认为你的结果最高对称性是 C(4V)，最好手建一下为好（如不考虑外场作用）。
这种体系的计算现在非常热门也非常有难度，上次我们听的以色列教授介绍的也是这方面的内容，这种体系以我们的水平可能很难处理好。你可能要找到做这方面的专家指导才会有很好的效果，否则很难分析清楚计算出的结果。
上次以色列教授沙松过来我们这边就是介绍的有关P450的计算，虽然我没有听懂什么，不过以我自己计算含金属的小分子经验来看，关于这种体系很难处理的部分就在于如果控制中心金属原子铁和周围配体的配位情况，以及了解中心金属原子的电子占据，如果你仅是靠软件产生体系的输入坐标，可能计算很难达到你想要的结果。这些坐标的产生可以要靠自己的经验来步步设置，所以就像你所说的太细节的东西是很难说清楚的，因为他们在计算过程中拥有比你更多的经验。只能给你这些参考性建议了，没有什么实质性的帮助。
我的输入文件是：
%chk=zux3
#p hf/lanl2dz opt symm=loose scf=(maxcycle=200) test zux3 0 1
74 0.000000 0.000000 1.463625 6 -2.059783 0.151849 1.234523 8 -3.189332 0.235187 0.995024 6 0.153438 2.073727 1.558845 8 0.238770 3.225066 1.644723 6 0.000000 0.000000 3.544137 8 0.000000 0.000000 4.703699 6 0.000000 0.000000 -0.664692 6 0.000000 0.000000 -1.891052 6 0.000000 0.000000 -3.251669 7 -1.186271 -0.125670 -3.989741 6 -2.427060 -0.359632 -3.218024 6 -1.381218 0.710276 -5.201852 1 -2.263789 -1.158724 -2.497766 1 -2.749787 0.537294 -2.683916 1 -3.209142 -0.660105 -3.913517 1 -0.432645 1.124763 -5.529600 1 -1.810550 0.110340 -6.005013 1 -2.062473 1.534158 -4.980011 6 2.059783 -0.151849 1.234523 8 3.189332 -0.235187 0.995024 6 -0.153438 -2.073727 1.558845 8 -0.238770 -3.225066 1.644723 7 1.186271 0.125670 -3.989741 6 2.427060 0.359632 -3.218024 1 2.263789 1.158724 -2.497766 1 2.749787 -0.537294 -2.683916 1 3.209142 0.660105 -3.913517 6 1.381218 -0.710276 -5.201852 1 0.432645 -1.124763 -5.529600 1 1.810550 -0.110340 -6.005013 1 2.062473 -1.534158 -4.980011
分子的点群是c2，但是优化第二轮是就出错了， Omega: Change in point group or standard orientation. Error termination via Lnk1e in d:\\g98w\\l202.exe.
如果我不用关键词symm=loose,优化时分子点群就变成了c1。 请高手指教，这是何原因呢？若我想保持分子仍是c2，改如何处理？ 3x!
换内坐标吧 opt=Z-matrix
你用内坐标可能会好一点吧？而且也能保证你的对称性的。[em12][em12] Angle between quadratic step and forces= 29.72 degrees. Linear search not attempted -- first point.
Iteration 1 RMS(Cart)= 0.00036483 RMS(Int)= 0.00000015 Iteration 2 RMS(Cart)= 0.00000015 RMS(Int)= 0.00000006
Variable Old X -DE/DX Delta X Delta X Delta X New X (Linear) (Quad) (Total)
R1 2.05001 -0.00006 0.00000 -0.00017 -0.00017 2.04984 R2 2.83834 0.00000 0.00000 0.00009 0.00009 2.83844 R3 2.05240 -0.00007 0.00000 -0.00023 -0.00023 2.05218 R4 2.05240 -0.00007 0.00000 -0.00023 -0.00023 2.05218 R5 2.03466 -0.00003 0.00000 -0.00007 -0.00007 2.03459 R6 2.47511 0.00010 0.00000 0.00015 0.00015 2.47526 R7 2.02743 -0.00009 0.00000 -0.00024 -0.00024 2.02719 R8 2.51854 0.00015 0.00000 0.00038 0.00038 2.51892 A1 1.94852 -0.00009 0.00000 -0.00069 -0.00069 1.94783 A2 1.88382 0.00007 0.00000 0.00042 0.00042 1.88424 A3 1.88382 0.00007 0.00000 0.00042 0.00042 1.88424 A4 1.93636 -0.00006 0.00000 -0.00034 -0.00034 1.93602 A5 1.93636 -0.00006 0.00000 -0.00034 -0.00034 1.93602 A6 1.87204 0.00008 0.00000 0.00063 0.00063 1.87267 A7 2.06277 0.00007 0.00000 0.00058 0.00058 2.06335 A8 2.15661 0.00001 0.00000 0.00012 0.00012 2.15674 A9 2.06380 -0.00008 0.00000 -0.00071 -0.00071 2.06310 A10 2.19767 0.00004 0.00000 0.00041 0.00041 2.19808 A11 2.13077 -0.00008 0.00000 -0.00050 -0.00050 2.13027 A12 1.95474 0.00004 0.00000 0.00009 0.00009 1.95483 D1 3.14159 0.00000 0.00000 0.00000 0.00000 3.14159 D2 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 D3 1.03955 0.00001 0.00000 0.00017 0.00017 1.03972 D4 -2.10204 0.00001 0.00000 0.00017 0.00017 -2.10187 D5 -1.03955 -0.00001 0.00000 -0.00017 -0.00017 -1.03972 D6 2.10204 -0.00001 0.00000 -0.00017 -0.00017 2.10187 D7 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 D8 3.14159 0.00000 0.00000 0.00000 0.00000 3.14159 D9 3.14159 0.00000 0.00000 0.00000 0.00000 3.14159 D10 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 Item value Threshold Converged? Maximum Force 0.000148 0.000450 YES RMS Force 0.000059 0.000300 YES Maximum Displacement 0.000912 0.001800 YES RMS Displacement 0.000365 0.001200 YES
Optimization completed. -- Stationary point found.
详细解释一下好吗，这一段几乎全看不懂
Iteration 1 RMS(Cart)= 0.00036483 RMS(Int)= 0.00000015 Iteration 2 RMS(Cart)= 0.00000015 RMS(Int)= 0.00000006 这应该是按笛卡尔坐标和内坐标计算的均方误差，比较的标准不清楚！
Variable列是输入的变量；Old X列是上一次的变量值；-DE/DX是能量对坐标的一阶导数，是判断驻点的依据；Delta X是与上次相比变量的变化；New X是新的变量数值！
SCF Done: E(RHF) = -113.866331170 A.U. after 12 cycles Convg = 0.1529D-08 -V/T = 2.0020 S**2 = 0.0000
NROrb= 34 NOA= 8 NOB= 8 NVA= 26 NVB= 26 最后一行是什么东东？
总轨道数34，Alpha占居轨道8个，Beta占居轨道8个，Alpha虚轨道26个，Beta虚轨道26个！ NROrb＝number of orb
NOA= numberof alpha occupied orb NVA=number of virtual alpha orb
结构优化的时候有没有可能出现这种情况：
即构型在两个或几个结构循环变化，总也不能达到终点 可能出现么？如果这样的话该如何处理？非线形搜索？ 多谢！
“构型在两个或几个结构循环变化，总也不能达到终点，可能出现么？”这个问题是你实际遇到的，还是凭空想出来的？几个结构在势能面上对应着局部最小，自洽场的结果是会到达一个局部最小值，而不会出现在几个个结构循环变化。我觉得你想说的是收敛振荡的问题吧？
他指的可能是结构优化振荡，我曾经用过iop（1/8)里的一些选项，对控制结构优化有一定帮助
我现在遇到的问题是这样的，优化了很多圈，但每次的能量总是在很小的范围内起伏
只能在Force上达到两个YES，而Displacement上总是NO,我觉得就象是在几个结构中振荡一样，这是所谓的优化振荡么？
iop（1/8)里有设定优化步长的功能：
0 DXMAXT = 0.1 BOHR OR RADIAN (L103, Estm or UnitFC).
= 0.3 Bohr or Radian (L103, Read or CalcFC). = 0.2 Bohr or Radian (L105). = 0.3 Bohr or Radian (L113, L114). N DXMAXT = 0.01 * N
哪里有关于iop（1/8)的详细内容呢？？
在免费资源版有IOP的pdf文件
谢supi!使用了iop(1/8=1)结构终于优化完了，但这说明什么呢？体系的能量对于结构的微小变化很敏感？还是说势能面比较平坦？
版上的今天刚贴的Guassian tips这篇贴子中最后有提到收敛遇到的一种情况就是\"taking large jumps\".对于优化默认采用的步长为0.3Bohr对你优化的结构来说太大,所以你看到结果能量总是在很小的范围内起伏.因此你使用IOP(1/8=1)设置步长为0.01Bohr,而达到收敛的最终结果,说明这个结构的势能面是起伏不平的,对结构的变化非常敏感.
不过,这种情况也像是oscillating,建议用SCF(Vshift=300)试一下,看能不能收敛到同一结果?
试了好几次,总觉得使用二次收敛方法scf=qc,对收敛没有什么太大的效果,反而不如以前.
而且对于oscillating采用SCF(Vshift=N)来说,由于DFT已经使用N=200,即使增大N的值对收敛的结果也无明显效果,比不上减小优化步长.但是对于减少步长来说,机时消耗太大,虽然在能量上比前者来得低,但是仍可能达不到收敛的预期效果,因此这些都说明了一个问题:The initial guess isn't very good.关于如果处理这些问题,《计算化学》上给出了一些更好的建议。
看了关于学术腐败的讨论，很惭愧，学习高斯已经有好几个月，还停留在Opt阶段，而且还有不少问题，先问两个吧，请莫见笑。
1.我发现构型优化时，Opt(tight)和只有Opt时，机时差别若干倍，但结果好像差的不是很多，为什么？ 2.对于大分子，可以先用小基组优化，而后再用更高级基组优化，具体写输入文件时，小基组时的分子说明部分用初始分子构型的数据，而更高基组优化时的分子说明部分应该用小基组的计算结果，怎么写呢？
。opt=tight采用了更为严格的收敛标准，使用选项Opt=Tight的优化计算比使用默认截止点的优化计算要多算几步。所用的时间当然也大得多了。要注意我们使用的能量的原子单位是个很大的数值1au＝
627.51kcal/mol，也许用两种方法算出来的能量只差了0。05个单位，可换算一下就有30多的kcal/mol！ 2。你可以先用小基组算，然后用前面的chk文件来作大基组的计算。在控制行加入geom＝check，guess＝read，就不用分子构型输入了
结构优化不收敛，一般有什么OPT=ReadFc, opt=CalcFc, opt=calcAll,或者设置maxcycle,opt=restart等等，但有没有人知道那种最为有效？具体方法怎样？当然有示例是最好的了！呵呵：） 谈谈我的想法:
我感觉结构不收敛可能与你的初始猜测有关,也许先固定某些键长键角再优化,事后再放开优化会比较容易收敛.
采用低级算法先计算只是其中的方法之一,你为何一直要用这种方法?有什么依据? 如果你顺利解决了你的计算,那就谈谈你的经验如何?
我觉得采用OPT=CALCFC并不一定得到好的结果,有时候收不收敛实际和我们处理的化学问题相关,而不仅仅是用计算上的技巧进行解决.比如两个原子间电子占据没有sigma成键,但是存在强的sigma非键作用,那么它们之间的斥力很强,计算就很难收敛形成稳定的体系.因此在处理收敛的问题时,我们应该多分析不收敛的问题实质,是电子的什么运动状态引起的,而非看到不收敛就想利用一些key来解决,我自己也试过一些key,包括opt=calcfc,allcfc,NoDIIS,QC等等(包括一些强制性的手段),并不是很有效,因为它们都是从技术的角度来解决问题.如果一个体系能形成稳定的体系,那么它不依靠这些key也是能解决的.反而如果一体系不收敛只是利用了这些技术就收敛了,那么这种结果如何可信?前面的贴子中,有些技巧诸如调整原子间的距离和键角,这些手段才是比较有效的,因为它直接影响了原子间的电子运动状态,即波函数,我们有时候初始猜测不正确就是因为程序给出的电子的占据状态错误而直接影响了原子间的成键,这时就需要对电子占据的轨道进行调整,Guess=alter,这对于过渡金属的计算来说特别要注意,因为G98总不能给出正确的初始猜测.
谢谢诸位意见，收益非浅。我起初一直认为我的初始模型是正确的，但后来我试过4个模型，才发现初始结构非常重要。开始三个是MnO2理想的晶体结构，通过改变键长希望结构优化收敛，但是失败。我又重新定义采用MnO2的真实晶胞参数，现在结果有些眉目，能收敛只是没有正常结束，可能是我吸附分子的坐标出了点问题，出现
Error on Z-matrix card number 23
angle Alpha is outside the valid range of 0 to 180.
所以我深切体会到楼上师兄的话：调整原子间的距离和键角,这些手段才是比较有效的。我现在才发觉，G98真是一个费内存，费时间的东西，很多计算要的是经验，你计算的垃圾多了，你才知道你需要什么！这个过程真够悲惨
。
因篇幅问题不能全部显示，请点此查看更多更全内容