求助：高压条件下的同位素动力学效应

dnlx

各位，我一直在模拟乙烷裂解这个地质反应。气体压力大约是500大气压，气体中95%是甲烷，5%是乙烷。温度是低温、长时间（200度，百万年尺度）。存在类似活性炭这样的催化剂（多孔材料，纳米级孔隙，存在自由基）。反应第一步是乙烷和自由基反应变成乙基（可逆），部分乙基转变成为不可逆的产物。

实际观察到的数据是，带有13C的乙烷反应得更快(inverse KIE)，12k/13k大概是0.98。如果我从reaction complex（乙烷+一个芳香基，不论是苯基还是芘基，代表活性炭的活性位）出发，DFT模拟可以得到Inverse KIE。但是如果从单独的反应物开始计算，总是normal KIE。我试验了其它催化途径，比如乙烷质子化，比如乙烷和黄铁矿（FeS2）反应，都是normal KIE。

现在想到的原因是1）因为这是个多孔材料，可能是传质控制；2）高压下的steric isotope effect；3) 高压下热力学上的非理想性。

关于1），传质过程是有利于轻同位素分子参加反应，而不是重同位素分子，所以想不出如何解释后者反应更快。

关于2），有一些文章（尤其是酶催化反应的文章）提到这个概念，但是没发现这样的量子化学模拟。它的根据是，重同位素分子及过渡态的相关键长更短，所以在体积局限的环境下（高压，或者多孔材料），反应更快。也有人做出了这样的体系（https://www.osti.gov/servlets/purl/983282）。
但是读社长的http://bbs.keinsci.com/thread-10159-1-1.html，一来极难模拟同位素对键长的影响，二来社长举例说：“对于BH（氢化硼），DBOC效应对平衡键长的影响才区区0.0007埃……对于更重的原子影响会更微乎其微”，所以很难想象高压和泛泛的多孔材料会带来明显的同位素效应。

关于3），这篇文章讲到同位素混合物的非理想性（https://www.degruyter.com/document/doi/10.1351/pac200476010011/pdf），这种效应应该很微弱，否则不会长期经常被忽略（见文章introduction）。但是我不知道它会不会带来百分之零点几的差别，以及能否进行理论模拟。

总之现在对模拟观察到的数据依然没有头绪，不知道是不是有哪里没考虑到或者没考虑对，请各位指教！

wzkchem5

我觉得更可能的原因是机理猜错了，比如可能是乙烷通过协同机理往活性炭（或者FeS2等）上面转移两个氢变成乙烯。百万年时间尺度下，常常会有实验室条件下观察不到的新机理出现，尤其对于乙烷这样的实验室条件下比较惰性的分子，可能完全走的是实验室条件下极难看到的非传统途径也说不定。

dnlx

wzkchem5 发表于 2022-6-25 21:36
我觉得更可能的原因是机理猜错了，比如可能是乙烷通过协同机理往活性炭（或者FeS2等）上面转移两个氢变成乙 ...

虽然有可能存在漏掉的机理，但是可能性不大，自由基抽氢应该是最优势的反应，自由能垒很低。其速度只是局限于自由基的低浓度。
乙烷参与的协同反应的能垒远高于抽氢反应，合适活性位的浓度更低，而且不能保持活性位的浓度（这个和自由基机理不一样）。
我估计了一下高压对KIE的影响，如果下面的公式（Basukevsky et al. 1985）可以用，



那么高压（恒压）下12k/13k = exp(p ΔΔV/RT)，ΔΔV是同位素替换引起的活化体积（volume of activation）的变化。即使用非常夸大的键长变化（同位素替换之后，过渡态C-H键长降低0.005埃）和夸大的分子动力学直径（1 nm），12k/13k 是0.996，和1差得不远，不能弥补S带来的KIE（12k/13k = 1.03）。

不知道上面的计算是否合适。另外，由于长链烷烃也观察到类似的现象，所以这个反应不会是在非常局限的小孔中发生的，所以孔隙的位阻可以不考虑。

现在不知道非理想性的影响该怎么估算一下。

wzkchem5

dnlx 发表于 2022-6-26 13:13
虽然有可能存在漏掉的机理，但是可能性不大，自由基抽氢应该是最优势的反应，自由能垒很低。其速度只是局 ...

你计算的时候有没有考虑13C取代的乙烷的symmetry number的影响？单个13C取代的乙烷的symmetry number和普通乙烷是不一样的，有的程序可能不能正确处理

dnlx

wzkchem5 发表于 2022-6-26 20:59
你计算的时候有没有考虑13C取代的乙烷的symmetry number的影响？单个13C取代的乙烷的symmetry number和普 ...

考虑了。如果对称数不对，结果会错得离谱。

wzkchem5

dnlx 发表于 2022-6-26 14:04
考虑了。如果对称数不对，结果会错得离谱。

那我感觉还是挺奇怪的，如果从reaction complex算是inverse KIE，从单独反应物算是normal KIE，那就意味着从单独反应物生成reaction complex这一步有很大的equilibrium isotope effect。但是这一步就是一个非共价吸附而已，不应该有很大的同位素效应才对。搞不好可能和低频振动的处理有关？

dnlx

wzkchem5 发表于 2022-6-26 21:09
那我感觉还是挺奇怪的，如果从reaction complex算是inverse KIE，从单独反应物算是normal KIE，那就意味 ...

是很奇怪，我可以保证数据处理没错，因为计算了不同的同位素取代（两个位置的13C，以及各个位置的氘），互相可以印证。

这个反应复合物的单点能低于反应物的加和，但是自由能高于反应物加和（几个kcal/mol）。确实如你所说，EIE不应该这么大。如果是低频振动处理带来的问题，请问哪个误差更小，是从反应复合物开始算，还是按照单独反应物加和？

=====================
补充：有可能又回到这个问题了：
http://bbs.keinsci.com/forum.php ... 2&fromuid=32447

wzkchem5

dnlx 发表于 2022-6-26 14:28
是很奇怪，我可以保证数据处理没错，因为计算了不同的同位素取代（两个位置的13C，以及各个位置的氘）， ...

该不该从反应复合物开始算和有没有低频没关系，应该总是取反应物和反应复合物里面自由能较低的那个。低频问题可以用Shermo算一下试试

dnlx

wzkchem5 发表于 2022-6-27 02:54
该不该从反应复合物开始算和有没有低频没关系，应该总是取反应物和反应复合物里面自由能较低的那个。低频 ...

Shermo（ilowfreq设为2，Grimme's entropy interpolation）跟ORCA的结果是一样的。比如过渡态，熵只差了0.004 cal/mol/K，最后一位有效数字的区别。
我现在的疑问在这里：如你前面指出的，“单独反应物加和”与“反应复合物”之间是个很弱的吸附，不应该有明显的同位素分馏。所以我得到的二者的熵，应该至少有一个不对。在排除数据处理出错的情况下，是不是因为模拟精度问题或者 计算振动熵的模型的问题。
附件是单独反应物、反应复合物、过渡态的输出文件。如果需要，我可以上传同位素替换结果（这个简单而罗嗦）。

【补充：明白你的意思了，是用Shermo里面Truhlar的方法去计算频率100以下的部分】

CompareReactantsComplex.zip (2.12 MB, 下载次数 Times of downloads: 2)

2022-6-27 04:15 上传 Uploaded

点击下载
Click to download
阅读权限: 40

wzkchem5

dnlx 发表于 2022-6-26 21:11
Shermo（ilowfreq设为2，Grimme's entropy interpolation）跟ORCA的结果是一样的。比如过渡态，熵只差了0 ...

哦，你本来用的就是orca是吧，那就不是Shermo能解决的问题了，Grimme的QRRHO一般比Truhlar的方法好。
可以试一下构象搜索。有一个可能性是你观察到的EIE是因为你只把乙烷里的一个C换成了13C，如果考虑到两个C都可能是13C、概率基本相等，可能EIE就抵消了。

dnlx

问题解决了。计算KIE需要知道熵变，但是熵既不能用ORCA默认值（Grimme方法），更不能用Truhlar方法。二者都存在近似，精度不够。应该自己用频率去计算振动熵，同位素异构体的频率都是一一对应的。这个本来是统计热力学的基础知识，但是之前一切交给软件，就把这些细节忽略了。

wzkchem5

dnlx 发表于 2022-7-1 14:26
问题解决了。计算KIE需要知道熵变，但是熵既不能用ORCA默认值（Grimme方法），更不能用Truhlar方法。二者都 ...

不对吧，你是怎么用频率计算熵的，套的是什么公式？如果你手算的方法确实比Grimme、Truhlar的方法好，那问题就来了，为什么他们没有把你的手算方法写到程序里？

dnlx

wzkchem5 发表于 2022-7-1 21:43
不对吧，你是怎么用频率计算熵的，套的是什么公式？如果你手算的方法确实比Grimme、Truhlar的方法好，那 ...

同样的基本公式（Grimme文章里的公式3）。问题在于低频处理，Truhlar的方法是把低频设成100，这个当然是经验性的。Grimme的方法包括了经验拟合和近似，见他文章公式5之后的讨论。这些近似在计算反应熵的时候或许足够准确，但是计算异构体的熵变之差（Delta Delta S）的时候，误差很可能就被放大了。

实际情况就是前面说的，从Grimme的方法出发，反应物到反应复合物之间有很大的同位素分馏，这个肯定不对。

所以还不如从基本公式开始计算。如果认为低频不应该用公式3，那就全部改成100，也就回到Truhlar的方法了（同位素异构体的频率一一对应，Truhlar方法的弊端，如你说的http://bbs.keinsci.com/forum.php ... 92259&fromuid=32447，在此没有影响）。
【补充：用量子化学的结果去计算KIE的时候，采用的方法 可能还得回到Bigeleisen 1950年前后的工作。一个例子是10.1016/j.gca.2008.12.022，公式7。】

wzkchem5

dnlx 发表于 2022-7-1 14:59
同样的基本公式（Grimme文章里的公式3）。问题在于低频处理，Truhlar的方法是把低频设成100，这个当然是 ...

没明白你是什么意思，所以你的方法和Truhlar的方法是等价的，只是后者是程序自己算的，前者是你手算的？那你上一个帖子为什么说“更不能用Truhlar方法”？

dnlx

wzkchem5 发表于 2022-7-1 23:06
没明白你是什么意思，所以你的方法和Truhlar的方法是等价的，只是后者是程序自己算的，前者是你手算的？ ...

根据得到的频率，用EXCEL算，公式推导见10.1002/9780470143476.ch2，实例见10.1016/j.gca.2008.12.022。