《计算化学-10.计算类型及高斯实现方法.ppt》由会员分享，可在线阅读，更多相关《计算化学-10.计算类型及高斯实现方法.ppt（49页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、袁淑萍袁淑萍国家重点实验室培育基地国家重点实验室培育基地408408室室E-mail:计算化学计算化学计算作业类型计算作业类型势能面扫描势能面扫描单点能单点能几何优化几何优化反应过渡态寻找反应过渡态寻找反应路径反应路径振动频率分析振动频率分析几何优化的目的几何优化的目的寻找势能面上的极小点寻找势能面上的极小点,确定分子可能的稳定结构确定分子可能的稳定结构几何优化得到的仅仅是势能面上的局部极小点!能量OPT的更多选项的更多选项Maxcycle=n:设定优化的最大步数设定优化的最大步数NoFreeze:活化所有冻结的坐标活化所有冻结的坐标(Geom=Check)Expert:在力常数和步长上放松限
2、制在力常数和步长上放松限制Tight,VeryTight:提高收敛标准提高收敛标准(力和步长力和步长)Loose:在工作初期使用在工作初期使用MaxStep=m:最大步长最大步长=0.01*m；IOP(1/8=m):改变步长改变步长 0.01*mRestart:优化作业的重启优化作业的重启GEOM的选项的选项Checkpoint:从Chekpoint文件中读入分子数据(通常还用到Guess=Read).Modify:读入并修改Checkpoint文件中的分子数据.NoCrowd:允许原子间距小于0.5A.NoKeep:丢弃关于冻结的信息.Step=n:在前一次优化失败后,把它的第n步作为起点进
3、行再次优化（Geom=(Check,Step=n)）Handling Difficult Optimization Cases复杂体系的优化复杂体系的优化q对于大体系，先先用用简简单单方方法法作作一一粗粗略略优优化化、后后用用高高级级方方法和高级基组优化法和高级基组优化，有时可能更快地达到目的。q后续作业可以从前步作业的结果中取得几何构型和波函数作为继续优化的起点，此时，后续作业再不需要输入分子几何（但仍需要分子电荷和自旋多重度)。q也可先用分子力学优化，作为SP的输入构型或进一步作量子化学优化的初始构型。分子力学优化相当快，机时远少于量子化学优化。几何优化的策略初学者有时急于求成，希望一步成
4、功，一开始就初学者有时急于求成，希望一步成功，一开始就采用高级方法和高级基组而忘了欲速则不达的道采用高级方法和高级基组而忘了欲速则不达的道理。理。q优化中途变量出错。最常见的现象是键角超越了最大限度180o：q可以减小步长OPT(stepsize=N)或在进程行增加Opt=calcfcq先用较低水平的方法和小基组作频率计算，再在后续的优化作业中使用Opt=Readfc。但这种计算往往是费时的。q初始构型或优化过程中的某种构型在势能面上处于较为平坦的区域，容易造成振荡,甚至永远不收敛。q解决的办法是加一个哑原子,重新定义键角。q步数超出q检查非常容易改变的坐标和/或强烈耦合的坐标q增加几何优化的
5、最大步数 OPT(Restart,Maxcyc=N)q在优化过程中点群改变了q检查结构和/或不使用对称性(NoSymm)批处理任务批处理任务任务比较多，互相之间没有关联，不能用多步处理的方法写在同一个输入文件里。编辑批处理文件，程序按次序调用并逐步计算计算作业类型计算作业类型势能面扫描势能面扫描单点能单点能几何优化几何优化反应过渡态搜索反应过渡态搜索反应路径反应路径振动频率分析振动频率分析搜索过渡态的目的寻找势能面上的一级鞍点，确定分子可能的稳定寻找势能面上的一级鞍点，确定分子可能的稳定结构之间的变换可能经过的态结构之间的变换可能经过的态一级鞍点满足的条件一级鞍点满足的条件:过渡态优化算法过渡
6、态优化的算法:TS，QST2，QST3寻找过渡态的实用建议关键词:Opt=QST2Optimization of Transition StatesOPT=TS -输入初始猜测的过渡态结构 -确认初始Hessian矩阵有一个本征向量有负的本征值 -可能的话计算力常数矩阵(CALCFC 或READFC)寻找过渡态用OPT=QST2 或或 OPT=QST3在冗余内坐标下寻找在冗余内坐标下寻找QST2:输入输入反应物附近反应物附近和和产物附近产物附近的两个结构的两个结构(通过对冗余坐标的线性插值估计过渡态通过对冗余坐标的线性插值估计过渡态)QST3:输入反应物、产物和估计的过渡态输入反应物、产物和估
7、计的过渡态Opt=QST2 的输入#B3LYP/6-31G(d,p)OPT=QST2H3CO-Title 10 2C1 0.0 0.0 0.002 0.0 0.0 1.3H3 0.0 0.9 -.3H4 0.8 -.2 -.6H5 -.8 -.2 -.6CH2OH-Title 20 2C1 0.0 0.0 0.002 0.0 0.0 1.4H3 0.0 0.92 1.7H4 0.7 -.1 -.7H5 -.7 -.1 -.7每个结构找到原子每个结构找到原子顺序应该一致顺序应该一致输入结构不需要是输入结构不需要是优化的稳定结构优化的稳定结构QST3 增加了第增加了第3个个标题和过渡态结构标题和过
8、渡态结构的估计的估计为过渡态估算Hessian矩阵初始的初始的Hessian矩阵必须有一个负的本征值矩阵必须有一个负的本征值,及及其合适的本征向量其合适的本征向量从低等级计算得到近似的Hessian矩阵从低等级的完整Hessian矩阵来计算(READFC:从频率计算中读入)在同等级计算完整的Hessian矩阵(CALCFC)在优化的每一步都重新计算完整的Hessian矩阵(CALCALL)，非常耗时，非常耗时过渡结构算例过渡结构算例#T UHF/6-31G(d)Opt=QST2 H3CO-H2COH Reactants 0,2CO 1 1.48H 1 R 2 AH 1 1.08 2 110.3
9、 120.H 1 1.08 2 110.3-120.R=1.08A=110.H3CO-H2COH products 0,2CO 1 1.48H 1 R 2 AH 1 1.08 2 110.3 120.H 1 1.08 2 110.3-120.R=1.9A=30.过渡结构的确认过渡结构的确认计算完整的计算完整的Hessian矩阵，确保有且仅有一个负本征值矩阵，确保有且仅有一个负本征值通过振动分析计算完成，可以使用连续作业通过振动分析计算完成，可以使用连续作业#B3LYP/6-31g opt=(ts,calcfc)freq确认本征矢量的方向正确连接到反应物和产物确认本征矢量的方向正确连接到反应物和
10、产物通过计算反应途径进行确认（通过计算反应途径进行确认（IRC）寻找过渡态失败怎么办1)有多余的负本征值沿着不是过渡向量的那些负本征值向量的方向寻找2)没有负本征值在反应坐标上扩大扫描,寻找最高能量对应的坐标点，沿最小本征值方向优化需要注意的问题需要注意的问题一般过渡态都涉及到键的断裂与生成，因此即使是自旋多重度为1的体系，也应该用非限制性的理论方法来寻找。设想一个键的均裂反应，原来成对的两个电子在反应过程中是要分别属于两个自由基的，不会占据同一个空间轨道了。计算作业类型计算作业类型势能面扫描势能面扫描单点能单点能几何优化几何优化反应过渡态搜索反应过渡态搜索振动频率分析振动频率分析反应路径反应
11、路径振振动动频频率率分分析析计计算算可可以以作作为为一一个个单单独独的的作作业业，但但更更常常见见的的是是附附在在分分子子构构型型优优化化之之后后成成为为多多步步作作业业，以验证优化结果是否稳定构型或过渡态。以验证优化结果是否稳定构型或过渡态。频频率率计计算算与与构构型型优优化化用用同同样样的的计计算算方方法法和和基基组组进进行行，关关键键词词是是Freq。下下面面是是频频率率计计算算作作为为后后继继作作业的输入实例：业的输入实例：%chk=H2O.chk#HF/6-31G(d)Freq Geom=Checkpoint Test H2O Frequencies 0 1频率计算的目的频率计算的目
12、的预测分子的红外和拉曼光谱预测分子的红外和拉曼光谱(频率和强度频率和强度)在在HF水平计算时，由于忽略相关能，计算的频率高于实水平计算时，由于忽略相关能，计算的频率高于实验值。通常在计算前或计算后利用校正因子来校正。验值。通常在计算前或计算后利用校正因子来校正。影响频率计算准确性的因素还有：确定力常数矩阵所需影响频率计算准确性的因素还有：确定力常数矩阵所需的微分技术的不准确性和平衡几何构型的不确定性等。的微分技术的不准确性和平衡几何构型的不确定性等。几何构型优化时几何构型优化时为难进行的优化计算力为难进行的优化计算力常数常数矩阵矩阵力常数是能量对坐标的导数。从数学的角度，导数显示力常数是能量对
13、坐标的导数。从数学的角度，导数显示着函数的变化趋势（可以根据导数值确定函数在某个范着函数的变化趋势（可以根据导数值确定函数在某个范围内是变大还是变小）。围内是变大还是变小）。优化的目的是寻找对应能量最低的构型。所以对难进行优化的目的是寻找对应能量最低的构型。所以对难进行的优化可以通过提供力常数值来指引优化的方向，使优的优化可以通过提供力常数值来指引优化的方向，使优化能够顺利完成。化能够顺利完成。确定势能面上驻点性质（极小点？鞍点确定势能面上驻点性质（极小点？鞍点？）？）虚频率的有无及其数目可以帮助我们判断优化结果是鞍虚频率的有无及其数目可以帮助我们判断优化结果是鞍点或能量极小构型，分子稳定构型
14、不应有虚频率，而过点或能量极小构型，分子稳定构型不应有虚频率，而过渡态则应该且只有一个虚频率。渡态则应该且只有一个虚频率。计算零点能和热化学数据，计算零点能和热化学数据，对体系总能对体系总能量、量、焓、熵进行校正焓、熵进行校正极化度与超极化度。极化度与超极化度。Hessian矩阵本征值与振动模式的关系矩阵本征值与振动模式的关系1，Hessian矩阵的本征值对应于简正振动的力常数矩阵的本征值对应于简正振动的力常数2，简正振动的，简正振动的频率频率等于力常数的等于力常数的平方根平方根3，简正振动的振动模式是相应本征值的本征矢，简正振动的振动模式是相应本征值的本征矢振动强度振动强度用于光谱指认IR光
15、谱的振动带的强度由偶极矩对简正模式的微商确定Raman光谱的振动带强度由极化率对简正模式的微商的平方确定振动分析算例振动分析算例%chk=freqhessian.chk#p UHF/6-31G(d)freq=(noraman)H3CO-H2COH TRANSITION STATE0 2 C 0.03104 0.63055 0.0 O 0.03104 -0.73696 0.0 H-0.99138 -0.13555 0.0 H 0.27839 1.12398 0.92614 H 0.27839 1.12398 -0.92614 注意注意:频率计算只能在优化频率计算只能在优化好的结构上进行，只能使用
16、好的结构上进行，只能使用与优化方法同样的方法来做与优化方法同样的方法来做频率计算。频率计算。振动频率的校正因子振动频率的校正因子计算得到的简正频率比实验值一般高约计算得到的简正频率比实验值一般高约10%这是由于谐振子近似和理论的近似而产生的这是由于谐振子近似和理论的近似而产生的热化学数据热化学数据所有频率计算都包含有体系的热化学分析。默认是在所有频率计算都包含有体系的热化学分析。默认是在298.15K一个大气压条件下计算得到。对每一种元素都使一个大气压条件下计算得到。对每一种元素都使用含量最高的同位素。下面是热化学的输出部分用含量最高的同位素。下面是热化学的输出部分:-Thermochemis
17、try-Temperature 298.150 Kelvin.Pressure 1.00000 Atm(计计算算的具体条件：温度，气压的具体条件：温度，气压)Atom 1 has atomic number 35 and mass 78.91834(计计算算所所用元素的确切质量用元素的确切质量)Atom 2 has atomic number 16 and mass 31.97207。可可以以通通过过参参数数ReadIsotopes根根据据自自己己的的具具体体需需要要设设定定计计算算条条件件Freq(ReadIsotopes)选项：选项：Specify alternate temperatur
18、e,pressure,and/or isotopesThis information appears in a separate input section having the format:temppressurescaleMustberealnumbersisotopemassforatom1isotopemassforatom2.isotopemassforatomn where temp,pressure,and scale are the desired temperature,pressure,and an optional scale factor for frequency
19、data when used for thermochemical analysis(the default is unscaled).%Chk=freq#HF/6-31G(d,p)Freq Test Frequencies at STP moleculespecification -Link1%Chk=freq%NoSave#HF/6-31G(d,p)Freq(ReadIso,ReadFC)Geom=Check Test Repeat at 300 K 0,1 300.0 1.0 16 2 3.The isotope masses for the various atoms in the m
20、olecule,arranged in the same order as they appeared in the molecule specification section.If integers are used to specify the atomic masses,the program will automatically use the corresponding actual exact isotopic mass(e.g.,18 specifies O18,and Gaussianuses the value 17.99916).后面四行的四个能量分别为E0，E，H，G。
21、计算关系为E0=E(elec)+ZPE(表示零K时体系的能量等于Eelec（就是优化计算中所得的SCF能量）和零点能之和)E=E0+E(vib)+E(rot)+E(trans)(在某一特定温度和压力下体系的能量等于零K时体系的能量加上内能值)H=E+RT(在某一特定温度和压力下体系的焓等于在这一条件下的能量和RT之和)G=H-TS(在某一特定温度和压力下体系的自由能等于在这一条件下的焓和TS的差。S表示熵)零点能是对分子的电子能量的矫正，表明了在零点能是对分子的电子能量的矫正，表明了在0K温度下的分子的振动能温度下的分子的振动能量对电子能量的影响。需要根据所用的算法采用不同的因子进行矫正。量对
22、电子能量的影响。需要根据所用的算法采用不同的因子进行矫正。在较高的温度下预测体系能量时，要给总能加上内能项。内能包括特定温度下的平动能，转动能和振动能，这些因素都受到具体条件的影响。E(Thermal)E(Thermal)ContributionstothethermalenergycorrectionCVCV Constantvolumemolarheatcapacity S S Entropy计算反应焓变计算反应焓变水合反应的标准反应焓变水合反应的标准反应焓变H298H+H2O-H3O+Example 8_02其中：其中：Ee0:0K时产物与反应物的能量差时产物与反应物的能量差;(Ee)2
23、98:0K到到298K电子能量的变化，此例中可忽略电子能量的变化，此例中可忽略;Ev0:0K时反应物和产物的零点能之差时反应物和产物的零点能之差;(Ev)298:0K到到298K振动能量的变化振动能量的变化;Er298:产物和反应物的旋转能之差产物和反应物的旋转能之差;Et298:产物和反应物的平动能之差产物和反应物的平动能之差;(PV):由于有一摩尔气态分子消失，由于有一摩尔气态分子消失，PV=-RT。其计算方法可以表示为其计算方法可以表示为H298=E298+(PV)E298=Ee0+(Ee)298+Ev0+(Ev)298+Er298+Et298各项能量可通过频率分析得到。各项能量可通过频
24、率分析得到。这这样样，所所要要做做的的工工作作就就是是分分别别对对反反应应物物和和产产物物进进行行优优化化然然后后进进行行频频率率分分析析得得到到所所需需数数值值。采采用用B3LYP/6-31G(d)就就能够得到足够精确的结果。能够得到足够精确的结果。(B3LYP/6-31G(d)opt,freq)这这里里注注意意我我们们不不用用计计算算H+，由由于于没没有有电电子子，它它的的电电子子能能量量显显然然是是0；由由于于只只有有一一个个原原子子，其其振振动动，转转动动能能显显然然也也是是零零，这这样样，其其只只有有平平动动能能，其其值值为为1.5RT=0.889 kcal/mol。(详见统计热力学详见统计热力学)。最最 终终 计计 算算 得得 到到 H298=-163.3kcal/mol,实实 验验 值值 为为-165.3 1.8kcal/mol，两者符合的相当好。，两者符合的相当好。