计算化学公社

标题: 水溶液中的计算，已经加了关键词scrf，还要不要再考虑水的传递 [打印本页]

作者
Author: qwe1832518773 时间: 2018-4-15 19:37
标题: 水溶液中的计算，已经加了关键词scrf，还要不要再考虑水的传递
起因是因为看到了一篇phys. chem. chem. phys. 的文章，文章本体在这儿
(, 下载次数 Times of downloads: 37)
里面提出了MEA的两性离子，在质子化另一个MEA时要先质子化水，再由质子化的水去质子化MEA，过程如下图所示（正文Fig 1）
(, 下载次数 Times of downloads: 32)
（因为自己只用过gaussian，而这篇文章中除了gaussian之外还用了VASP，所以我理所当然的认为上图的效果是VASP所呈现的，所以我在“第一性原理-VASP”板块发布了对于这一个疑问更详细的问题描述，各位老师各位同仁有知道的希望可以“凶狠”的指正我！！真心求教）

所以我也用gaussian算了一下MEA两性离子质子化MEA的这一步，按照文献Fig 1中的样子分别算了下面两个反应路径的活化能
1、MEA+[MEA+COO-]=[MEA+]+[MEACOO-]
2、MEA+[MEACOO-]+H3O+=[MEA+]+[MEACOO-]+H2O
两个路径的输入都一样，如下
# opt=(calcfc,ts,noeigen,gdiis) freq b3lyp/6-311++g(d,p) scrf geom=connectivity
可是最后得到的结果是：1的活化能比2低。还换成过6-31g的基组算了一次，得到的结果也是1的活化能比2低。。。分析这个结果得到的是，路径1比路径2更可能发生，那这样就和这篇文献中提出的路径不一致了。。。自己想了想有以下几个不确定的地方，望各位不吝赐教：
1、文献用了VASP直接找到如Fig 1所示的反应路径，与我用gaussian计算活化能来确定反应路径的方法有本质区别，所以我和他的结果无法比较
2、我在计算的时候加了scrf，已经考虑了水参与，所以我所做的路径2（先质子化水）根本就不存在

这只是我自己凭目前的浅薄水平能想到的，看起来可能很蠢，如果我哪里想错了或者想的不严谨，希望各位可以毫不留情的指正，非常感谢！！

作者
Author: mutron 时间: 2018-4-15 21:20
文章的SI有提到1的活化能比2低（Figure S3和Figure S4），但是原文也说到：

quote
"We also consider direct proton transfer from MEA+COO- to a nearby MEA, which turns out to be very facile with a negligible barrier provided they are in the right orientation (see Fig. S4, ESI†); however, solvation of MEA+COO- is more likely in the 30 wt% MEA solution(2MEA, 1CO2, 16H2O)"

even if path 1 is more favorable, proton from MEA+COO- is quite hard to reach a nearby MEA without bumping into a H2O molecule
提供点思路，不知道是否说得通

作者
Author: sobereva 时间: 2018-4-15 23:55
这个过程完全没必要用VASP。“还换成过6-31g的基组算了一次”这种做法毫无意义，看
谈谈量子化学中基组的选择
http://sobereva.com/336

这种涉及氢键问题用M06-2X或B3LYP-D3会好得多。

隐式水模型显式水模型没有矛盾之处，所以不要考虑“所以我所做的路径2（先质子化水）根本就不存在”

如果文献没有计算MEA+[MEA+COO-]=[MEA+]+[MEACOO-]的活化能，那么也没有理由说这个过程发生不了。我没仔细看这文章，不过如2L提示的，作者或许只是由于考虑到这种情况在溶液里出现的几率很小，所以才没有去特意提。

作者
Author: qwe1832518773 时间: 2018-4-16 11:21
十分感谢！两位老师的回复很有帮助！收益匪浅，自己一直钻在活化能低更容易发生的牛角尖里出不来，忽略了反应的实际情况，确实在这样一个体系中完全忽略水的存在是不合理的。而且对于自己的实验结果我也有些新的感悟，再次感谢二位！

作者
Author: qwe1832518773 时间: 2018-4-19 11:32

mutron 发表于 2018-4-15 21:20
文章的SI有提到1的活化能比2低（Figure S3和Figure S4），但是原文也说到：

quote

老师，我还有一些疑问，已经纠结我很长时间了。您也看过我发在一楼的那篇文章了，文章中提出有水参与时能垒为0.08eV，而无水参与只有0.005eV。但是我自己的计算结果显示有水参与（我不知道文章中究竟是怎么计算有水分子网络传递质子这个过程的能垒的，我自己只是单纯的在无水参与的基础上，在两个MEA中间加了一个水分子）的能垒为0.27eV，而无水（无水的反应完全与文献中一样）的能垒达到了0.16eV。并且我使用的基组也是b3lyp/6-311++g(d,p)，和文献中提到的一样。
我的疑问有三点：
1、为什么我的无水反应得到的能垒与文献中的差别如此之大？（我用同样的方法验算过文献正文Fig 5的能垒，是与文献结果一致的）
2、文献中所提出的经水分子网络传递质子的过程我究竟该如何做？（我应该按照文献中所述的把2个MEA，1个CO2和16个水全部放进高斯的TS计算里吗？那这么多分子我又该如何排列呢？）
3、第三个问题倒不是针对这篇文献，是我自己想问的。在做TS计算的时候，初始的猜想结构如果摆向不同的化，即初始的几何结构不同，最后得到的过渡态会是一样的吗？

作者
Author: mutron 时间: 2018-4-20 22:18
本帖最后由 mutron 于 2018-4-23 20:18 编辑

1、由于原文没有给出任何结构坐标，得自己搜索出他们对应的结构再优化、振动分析，需要点时间。如果能提供你计算的坐标，或许能看出来差别在哪（SI里面哪怕没水传递质子的过程他们也加了一个显式水，而且反应物构象也不少，构象不同计算结果也很有可能有不小差别）

在尝试重复他们的计算结果过程中，留意到原文的计算水平是B3LYP/6-311++G（他们文章是这么写的），而不是6-311++G**，这有点奇怪

2、原文是用AIMD做的动力学模拟，这个你得问用过VASP的老师

而且过程也要看是什么内容，如果讨论两性离子传递质子给水分子，原文通过截取动力学模拟的快照说明有这个过程，然后通过对比加两个显式水分子和加一个显式水对能垒降低的趋势，说明水分子（排布）对质子转移过程有重要作用；最后文章用RDF讨论随着反应的进行，两性离子向水分子转移质子的可能性不断减少

3、如果变化的趋势是另一个过渡态的构象，则有可能会不一样（比如本例几个不同的反应物构象他们的质子转移过渡态也不一样）；如果变动的只是反应路径上的原子，很有可能得到一样的过渡态

作者
Author: mutron 时间: 2018-4-23 20:50
本帖最后由 mutron 于 2018-4-24 00:18 编辑

qwe1832518773 发表于 2018-4-19 11:32
老师，我还有一些疑问，已经纠结我很长时间了。您也看过我发在一楼的那篇文章了，文章中提出有水参与时能 ...

如果使用6-311++G**基组在PCM下算自由能，则SI里面的Figure S4的能垒为0.018 eV（int=fine）。6-311++G基组在此条件下甚至得不到过渡态，所以也就重复不了文章的结果，不过相差也才0.3kcal也基本差不多了

你可以把你优化好的结构坐标贴出来看一下，相差几个kcal我估计原因仅仅是你计算的构象并不是文章里面对应的结构（其实Figure S4反应物的结构甚至不是能量最低的构象）

顺便这是我优化好的结构，你也可以对比你自己的计算结果：
MEA+COO- MEA H2O
0 1
O                2.46450300 1.44052900 -1.21067500
C                1.99400600 2.02613100 -0.00010700
C                1.37912300 0.97429900 0.92291500
N                0.39904900 0.12891900 0.17103800
C                0.46205300 -1.43445600 0.42151000
O                1.58431800 -1.85775700 0.68722200
O                -0.63779100 -1.96598600 0.27986300
H                -0.59476800 0.45924900 0.35186200
H                0.59669100 0.24806500 -0.82915400
H                3.12244000 0.74862600 -0.99452200
H                2.79550700 2.53905500 0.54255300
H                1.24991100 2.77113200 -0.28911900
H                2.14737000 0.31100000 1.31512700
H                0.86079800 1.44019600 1.76163000
O                3.96470700 -0.79242300 -0.38833400
H                3.21744400 -1.31081900 -0.03049000
H                4.40666300 -1.36462700 -1.02529200
O                -2.81403600 -0.99302900 -1.18815100
C                -3.76035900 -0.27691500 -0.40064900
C                -3.32768100 1.14338100 -0.05612100
N                -2.17238500 1.16534000 0.86915700
H                -1.94349900 2.13095300 1.09337300
H                -2.43743000 0.73538200 1.75362200
H                -2.08147200 -1.29547800 -0.62016500
H                -3.99205300 -0.82873200 0.52096000
H                -4.68065100 -0.22276700 -0.98769000
H                -3.03618800 1.66238400 -0.97317400
H                -4.19386700 1.67124600 0.36442500

TS
0 1
O                -2.73857200 0.77235900 1.31610100
C                -2.06450900 1.81510700 0.62465800
C                -1.13963300 1.32767400 -0.48972500
N                -0.10503500 0.36938400 -0.02383000
C                -0.39555100 -1.11022700 -0.09304400
O                -1.30333800 -1.46719400 -0.85505500
O                0.39569400 -1.78440500 0.59296700
H                0.96229400 0.51446500 -0.67740300
H                0.17991900 0.57758000 0.93306200
H                -3.30126300 0.28164600 0.68601000
H                -2.77346100 2.52811900 0.18478900
H                -1.48291400 2.35435800 1.37783400
H                -1.69817800 0.83438000 -1.28391600
H                -0.63261300 2.19471600 -0.92267000
O                -3.96558900 -0.82448300 -0.65158600
H                -3.07260000 -1.19382000 -0.81331500
H                -4.52116100 -1.56908300 -0.39651300
O                2.66432100 -0.48144900 1.43824100
C                3.60533200 -0.35503200 0.37865800
C                3.33811500 0.83481600 -0.53768200
N                2.11654600 0.65412400 -1.37345200
H                2.00270000 1.45610200 -1.99098600
H                2.22552700 -0.15522500 -1.98408900
H                1.89966400 -1.02121800 1.14529400
H                3.65754800 -1.27987500 -0.21004200
H                4.58802100 -0.20378200 0.83282300
H                3.19618300 1.73609200 0.06191500
H                4.20709600 0.98834000 -1.18433800

作者
Author: qwe1832518773 时间: 2018-7-2 11:58
本帖最后由 qwe1832518773 于 2018-7-2 12:00 编辑

mutron 发表于 2018-4-23 20:50
如果使用6-311++G**基组在PCM下算自由能，则SI里面的Figure S4的能垒为0.018 eV（int=fine）。6-311++G基 ...

很感谢你认真的解答！~
我用你优化完的坐标算了一下，用的当然也是B2LYP/6-311++G**，活化能是0.028eV，我对比了一下我之前的结果（0.165eV）发现差别就在于：
a.你和那篇文献所用的结构在两性离子的羟基处都连有一个水分子
b.我之前的结果中没有连水分子（因为我自己当初认为这个水分子只是作者画上去证实是在水环境下，并没有参与反应的）

没有想到就是这个水分子使得结果差了这么多。。所以我还是想向你请教2个问题：
1.为什么只有两性离子的羟基连有水分子而边上MEA的羟基没有连，如果是真实情况的话不应该两者的羟基都会连水分子吗？
2.这个反应（只讨论无水传递，即文献FIG S4）是如何确定的呢？是因为在不同情况下（指的是与羟基相连的水分子的个数），这个反应的能量最低吗？还是用VASP先做的模拟后直接确定了这种反应方式呢？

十分不好意思。问题隔了很长时间，可能你都忘了，但还是希望能为我解答一下疑惑，感谢！
附图：
(, 下载次数 Times of downloads: 36) 图a

(, 下载次数 Times of downloads: 41)
E:\bb.png

图b

作者
Author: sobereva 时间: 2018-7-2 22:16
B2LYP根本没法用
如果指的B2PLYP，基组起码def2-TZVP，6-311++G**这么low的只能给非双杂化泛函用

作者
Author: mutron 时间: 2018-7-4 20:27
本帖最后由 mutron 于 2018-7-4 20:28 编辑

qwe1832518773 发表于 2018-7-2 11:58
很感谢你认真的解答！~
我用你优化完的坐标算了一下，用的当然也是B2LYP/6-311++G**，活化能是0.028eV， ...

因为当初的目的是重复文献结果，我尽量摆成文献的结构进行计算，所以有个水分子。

1.真实情况是每个可以成氢键的基团都很可能会与溶剂水分子成氢键，包括N-H、O-H等。但是由于这是原文计算MEA+COO-直接向MEA转移质子的反应，水分子并不参与到反应当中，可以预见加不加显式水分子对结果的影响不会太大。我把Fig.4结构里的水分子删掉重算了一遍能垒为 0.020 eV（int=fine），基本和之前的一样。至于原文为什么要（在那个地方以那种方式只）加一个显式水，它并没有相应的说明，由于对结果影响不大可以不用管它。而且哪怕是你计算的结果0.165 eV（3.8kcal），能垒也太小了，室温瞬间反应完毕，也不影响原文结论
减掉H2O的TS：
0 1
O                -3.77113200 -0.08640500 -0.07809200
C                -2.88294300 -0.92531200 -0.80991700
C                -1.58365900 -1.25292800 -0.07176800
N                -0.48806700 -0.25155500 -0.20414300
C                -0.68046200 1.14910200 0.32959000
O                -1.71461600 1.36459000 0.97519100
O                0.28156500 1.89459300 0.06974200
H                0.54575800 -0.68492600 0.37278700
H                -0.19620900 -0.16913700 -1.17819600
H                -3.21410500 0.58144500 0.36567800
H                -3.41558200 -1.86150000 -0.99247200
H                -2.64715100 -0.48643900 -1.79001400
H                -1.78164000 -1.39527200 0.99165000
H                -1.17856300 -2.18856900 -0.46466500
O                2.47111500 0.82450000 -1.18524100
C                3.33951900 0.21016600 -0.24156700
C                2.89123600 -1.18293400 0.18703700
N                1.64816400 -1.15319100 1.00918500
H                1.42542500 -2.09526000 1.32635900
H                1.79974800 -0.60232700 1.85365400
H                1.72868000 1.26158500 -0.71567200
H                3.47171300 0.84912600 0.64146500
H                4.31762600 0.10936300 -0.71861100
H                2.68876600 -1.79403500 -0.69476200
H                3.69725500 -1.65762500 0.75445000

2.找出反应物（结合显式水或没有显式水）的最稳定构象和相应的过渡态。你的计算结果和原文差别也只有3kcal，很有可能就是构象差异造成能量差别而已

作者
Author: qwe1832518773 时间: 2018-7-12 13:00

sobereva 发表于 2018-7-2 22:16
B2LYP根本没法用
如果指的B2PLYP，基组起码def2-TZVP，6-311++G**这么low的只能给非双杂化泛函用

谢谢sob老师解答！

作者
Author: qwe1832518773 时间: 2018-7-12 13:02

mutron 发表于 2018-7-4 20:27
因为当初的目的是重复文献结果，我尽量摆成文献的结构进行计算，所以有个水分子。

这几天在科音听sob老师系统的讲了量化，我想极大可能是因为我用了加弥散函数的基组（自己一开始都没意识到。。。）所以结果和文献的差别这么大

此外，非常感谢您的解答！！！

作者
Author: Norman 时间: 2021-3-26 03:58

sobereva 发表于 2018-4-15 23:55
这个过程完全没必要用VASP。“还换成过6-31g的基组算了一次”这种做法毫无意义，看
谈谈量子化学中基组的 ...

请问这个反应计算过渡态和能垒，用6-311+G**应该和6-311++G**差距不大吧？用def-TZVP是不是差点呢？

作者
Author: sobereva 时间: 2021-3-26 09:14

Norman 发表于 2021-3-26 03:58
请问这个反应计算过渡态和能垒，用6-311+G**应该和6-311++G**差距不大吧？用def-TZVP是不是差点呢？

此文明确说了别随便给氢加弥散
谈谈弥散函数和“月份”基组
http://sobereva.com/119

该用什么这里明确说了
谈谈量子化学中基组的选择
http://sobereva.com/336（http://bbs.keinsci.com/thread-3545-1-1.html）

欢迎光临计算化学公社 (http://bbs.keinsci.com/)