求助：温度320 K在乙醇中[H3O]+的自由能计算求审核

liaorongbao

温度320K在乙醇中[H3O]+的自由能计算(求审核)

现假设[H3O]+离子在真空中的结构和在乙醇溶剂中的结构有显著差别。那么如何计算[H3O]+离子在乙醇溶液中的自由能呢？并且要求气相自由能的计算采用http://sobereva.com/327中提到的“讲究的人”的算法。

1 计算气相(320 K, 1 atm)的自由能
1.1 先做气相(320 K, 1 atm)的自由能计算。

考虑到这个离子最终的结构是溶剂中的结构。因此采用B3LYP/6-311G* em=GD3BJ opt freq scrf(SMD,solvent=ethanol) temperature=320获得溶剂中的结构，输入文件记为111.gjf。Sob说，对于多数中性且没有局部显离子性的体系，带不带隐式溶剂模型对结构和自由能的热校正量影响甚微。这里刻意带了溶剂模型，是因为假设溶剂中的结构显著不同于气相结构。不带溶剂模型优化出的结构与溶剂环境下的实际状态可能定性不符而无法接受，索性就带着溶剂模型。即对应SMD溶剂模型表现乙醇环境的情况下用B3LYP-D3(BJ)/6-311G*对体系进行优化和振动分析。这个级别不昂贵，而用于opt freq目的足够给出可靠的气相下标况下(320 K, 1 atm)的自由能。在111.out中搜索Gibbs可找到下面这行Thermalcorrection to Gibbs Free Energy=  0.015633 Hartree。

我猜想：把111.out结构取出后，再做一次不带溶剂的B3LYP/6-311G* em=GD3BJ freq temperature=320计算，把计算文件叫112.gjf，结果文件叫112.out。得到的Thermal correction to Gibbs FreeEnergy= 0.016065数据(这个真空计算是带着色散效应，不知是否应该带色散效应，我认为应该带着色散。不带色散效应的113.gjf的B3LYP/6-311G* freq temperature=320得113.out中结果是0.016066)。

问题1：三个数据0.015633，0.016065，0.016066接下来用哪个最好(我接下来用的是0.015633，但不知是否正确)？

http://bbs.keinsci.com/forum.php?mod=viewthread&tid=3805给出B3LYP/6-311G(d,p)的ZPE矫正因子是0.9888(来自文献10.1021/jp073974n)，没查到B3LYP/6-311G*的ZPE的矫正因子，只好采用B3LYP/6-311G(d,p)的0.9888作为近似替代。http://sobereva.com/221给出了各种计算级别对应的各种热力学数据的矫正因子。

1.2 采用气相结构计算气相(320 K, 1 atm)下电子能

采用以上111.out的输出结构制作121.gjf文件，关键词是B2PLYP/def2TZVP em=GD3BJ(Sob博文http://sobereva.com/327中没说计算电子能这一步要加上温度，我理解成低温下电子不发生跃迁)使用较好的B2PLYP-D3(BJ)/def2-TZVP级别计算真空下较高精度电子能量。结果文件为121.out。

问题2：如果想得到更精确的电子能，也可以采用CCSD(T)/cc-pVTZ em=GD3BJ格式计算电子能么？

然后在121.out中读取电子能量MP2=-76.6814719 Hartree。Sob说读双杂化泛函算的电子能量的时候千万别读错地方，不懂的话看《谈谈该从Gaussian输出文件中的什么地方读电子能量》（http://sobereva.com/488）。Sob说Gaussian算大体系的双杂化泛函巨慢还耗硬盘，用ORCA算省时得多，这里专门提到了《简谈量子化学计算中DFT泛函的选择》（http://sobereva.com/272）。如果只能用Gaussian，当体系超过50原子的话，建议用M06-2X代替B2PLYP-D3(BJ)，便宜多而精度通常下降不多（对主族体系而言）

1.3 气相标况(320 K, 1 atm)自由能讲究的人的计算
1.3.1settings.ini文件的修改

121.out中B2PLYPD3/def2TZVP得1 atm气相电子能MP2=-76.6814719 Hartree。用文本编辑器编辑Shermo目录下的settings.ini文件，将里面的E=后面的值改为-76.6814719，注意等号和数值之间必须有个空格。再将T=后面的值设为320，代表在320 K下分析。B3LYP/6-311G*级别下(采用B3LYP/6-311G(d,p)代替)的ZPE校正因子为0.9888，因此将sclZPE=后的值设为0.9888。ilowfreq设2使用Grimme的准RRHO模型计算热力学量。

问题3：这个0.9888相关的问题有没有更好的处理方法？

1.3.2 气相标况(320 K, 1 atm)情况下的自由能的获得

打开Shermo.exe，然后将111.out的图标拖到Shermo窗口里，此时路径就自动产生了。按回车后首先看到分子的信息。最后得到的Sum of electronic energy and thermal correction to G:   -76.6684846 a.u.就是精确的(320 K, 1 atm)气相吉布斯自由能(若111.gjf中温度是298.15 K做完这一步后也得同样结果)。

2 计算溶解自由能

2.1 采用111.out的输出结构和关键词M052X/6-31G*建立文件21.gjf，在真空下计算电子能量。21.out中查到SCF Done:  E(RM052X) = -76.6720728591

2.2采用111.out的输出结构和关键词M052X/6-31G*scrf(SMD,solvent=ethanol)构建文件22.gjf计算溶剂模型乙醇环境电子能，22.out中查到E(RM052X) =  -76.8228157375  。

2.3计算在乙醇中的溶解自由能

从2.1步和2.2步前后都是1M标准态。因此溶解自由能为-76.8228157375 + 76.6720728591=-0.150742878 Hartree (折合-94.592588 kcal/mol)。只能得到298.15 K的数据，无法得到320 K的数据(见http://sobereva.com/327)。只能用298.15 K的数据近似代替320 K的数据。

3 计算溶质在1 M的溶液中的自由能

1atm标态的气相自由能-76.6684846 au。从1 atm气态压缩到1M(即1 mol/L)气体过程中吉布斯自由能增加+1.89 kcal/mol，1.89/627.51 = 0.003011904 au。溶剂模型算的1 M气体变成1M溶液过程的溶解自由能-0.150742878 au。因此320 K时1M标准态的[H3O]+在乙醇中的自由能如下：

-76.6684846+ 0.003011904 - 0.150742878 = -76.81621557 au

问题4：以上过程有没有其他错误？

问题5：针对溶液中化学反应的吉布斯自由能的计算；用以上方法分别计算产物和反应物，然后用产物自由能减反应物即可？

北大-陶豫

别的我不太懂，但是 CCSD(T)/cc-pVTZ 肯定不需要 em=GD3BJ

liaorongbao

北大-陶豫 发表于 2024-7-1 15:49
别的我不太懂，但是 CCSD(T)/cc-pVTZ 肯定不需要 em=GD3BJ

谢谢交流指导