计算化学公社

标题: 计算得到的吸附自由能比吸附能弱，请问是否正常？ [打印本页]

作者
Author: a9471163 时间: 2020-11-2 23:31
标题: 计算得到的吸附自由能比吸附能弱，请问是否正常？
体系：A,B两个分子在水溶液中的稳定吸附体系
吸附能计算方法：gromacs中rerun定义两个能量组，分别是A，B分子，然后通过gmx energy命令，得到静电相互作用能和非静电相互作用能，总能量约为-100kJ/mol
吸附自由能计算方法：gromacs伞形采样后计算PMF曲线，得到横坐标是质心距离，纵坐标是能量的图，最后用PMF曲线最高点的y值，减去最低点的y值，计算得到的值是40kJ/mol，明显小于100kJ/mol
我个人觉得结果不正常，更换体系后（换了两个不同的分子，在水溶液中的吸附体系）再次计算，仍然是吸附自由能比吸附能弱，不清楚问题出在哪，还请老师指教，谢谢各位~
伞形采样时用到的umbrella.mdp参数如下
title    = Umbrella  simulation
define    = -DPOSRES_B
; Run parameters
integrator  = md
dt       = 0.002
tinit    = 0
nsteps    = 5000000 ; 10 ns
nstcomm    = 10
; Output parameters
nstxout    = 50000    ; every 100 ps
nstvout    = 50000
nstfout    = 5000
nstxtcout = 5000    ; every 10 ps
nstenergy = 5000
; Bond parameters
constraint_algorithm = lincs
constraints          = all-bonds
continuation          = yes
; Single-range cutoff scheme
nstlist    = 5
ns_type    = grid
rlist    = 1.4
rcoulomb = 1.4
rvdw       = 1.4
; PME electrostatics parameters
coulombtype    = PME
fourierspacing  = 0.12
fourier_nx    = 0
fourier_ny    = 0
fourier_nz    = 0
pme_order    = 4
ewald_rtol    = 1e-5
optimize_fft = yes
; Berendsen temperature coupling is on in two groups
Tcoupl    = Nose-Hoover
tc_grps    = Other Water_and_ions
tau_t    = 0.5    0.5
ref_t    = 298.15    298.15
; Pressure coupling is on
Pcoupl       = Parrinello-Rahman
pcoupltype    = isotropic
tau_p          = 1.0
compressibility = 4.5e-5
ref_p          = 1.0
refcoord_scaling = com
; Generate velocities is off
gen_vel    = no
; Periodic boundary conditions are on in all directions
pbc    = xyz
; Long-range dispersion correction
DispCorr = EnerPres
; Pull code
pull                   = yes
pull_ngroups          = 2
pull_ncoords          = 1
pull_group1_name       = MOL
pull_group2_name       = SUB
pull_coord1_type       = umbrella    ; harmonic biasing force
pull_coord1_geometry = distance    ; simple distance increase
pull_coord1_groups    = 1 2
pull_coord1_dim       = Y N N
pull_coord1_rate       = 0.0
pull_coord1_k          = 1000
pull-coord1-start             = yes          ; define initial COM distance > 0

作者
Author: a9471163 时间: 2020-11-2 23:33
分子有明显的疏水作用，所以我预期结果应该是吸附自由能比吸附能（即静电能+非静电能）更强，计算结果恰恰相反，我个人觉得可能是PMF计算有误，因为我对PMF计算所需的合理参数，还不太熟悉

作者
Author: wzkchem5 时间: 2020-11-3 08:21

a9471163 发表于 2020-11-2 23:33
分子有明显的疏水作用，所以我预期结果应该是吸附自由能比吸附能（即静电能+非静电能）更强，计算结果恰恰 ...

可以这样去理解这个问题：当分子没有疏水作用时，吸附自由能应当比吸附能弱，因为吸附过程熵减。所以当分子虽然有疏水作用但是不是很强的时候，吸附自由能比吸附能弱也是有可能的。
你觉得疏水作用明显可能是跟你研究的其他体系比的，不一定意味着站在自由能的角度看也可以说疏水作用明显。

作者
Author: a9471163 时间: 2020-11-3 09:20

wzkchem5 发表于 2020-11-3 08:21
可以这样去理解这个问题：当分子没有疏水作用时，吸附自由能应当比吸附能弱，因为吸附过程熵减。所以当分 ...

你好，我是这样理解这个问题的：
ΔG=ΔU+PΔV-TΔS，目前ΔG=-40，ΔU=-100，PΔV项可以忽略不计，因为NPT后盒子体积的变化不大，代入公式得到：
-40=-100-TΔS
熵是表示系统混乱度的一个量,混乱程度越大,熵越大。所以在自发吸附过程中,只要是自发过程,系统的熵是必然增加的（网上查到的说法）
所以ΔS项是正值
最终等式左边≠右边，说明我还是觉得自己计算结果有问题，大概率是PMF计算出了问题，因为gromacs的rerun计算吸附能并没有什么难度

作者
Author: a9471163 时间: 2020-11-3 09:54

wzkchem5 发表于 2020-11-3 08:21
可以这样去理解这个问题：当分子没有疏水作用时，吸附自由能应当比吸附能弱，因为吸附过程熵减。所以当分 ...

请问有没有报道过吸附过程熵减小的文献呢？假如我的体系，吸附过程确实是熵减小的，也就能解释为什么吸附自由能比吸附能弱了

作者
Author: wzkchem5 时间: 2020-11-3 10:45

a9471163 发表于 2020-11-3 09:20
你好，我是这样理解这个问题的：
ΔG=ΔU+PΔV-TΔS，目前ΔG=-40，ΔU=-100，PΔV项可以忽略不计，因为 ...

不是不是，你搞混了
ΔG=ΔH-TΔS的ΔS是等温条件下的熵变，自发过程ΔS必然大于0的ΔS是绝热条件下的熵变。所谓ΔS必然大于0的说法只有对于孤立体系（NVE/NPE系综）才成立，而你用的是NPT系综，NPT系综的自发过程只有ΔG<0的结论，没有ΔS>0的结论。

作者
Author: wzkchem5 时间: 2020-11-3 10:59

a9471163 发表于 2020-11-3 09:54
请问有没有报道过吸附过程熵减小的文献呢？假如我的体系，吸附过程确实是熵减小的，也就能解释为什么吸附 ...

多的是，没准比熵增大的例子还多。
给个最简单的例子，DNA单链结合成双链，除了端基的AT碱基对以外，所有碱基对对于熵变的贡献都是负的。参见https://en.wikipedia.org/wiki/Nucleic_acid_thermodynamics，Table 1

作者
Author: sobereva 时间: 2020-11-4 04:17

a9471163 发表于 2020-11-3 09:54
请问有没有报道过吸附过程熵减小的文献呢？假如我的体系，吸附过程确实是熵减小的，也就能解释为什么吸附 ...

对于溶质来说，熵减小是普遍现象。吸附前单体可以自由运动，吸附后运动受限，因此熵减小才是常见的。
虽然溶剂下有疏水作用，但也别忘了溶剂也等效地极大地削弱了溶质间的静电作用。

作者
Author: a9471163 时间: 2020-11-4 09:07

wzkchem5 发表于 2020-11-3 10:59
多的是，没准比熵增大的例子还多。
给个最简单的例子，DNA单链结合成双链，除了端基的AT碱基对以外，所 ...

非常感谢老师，我明白了

作者
Author: a9471163 时间: 2020-11-4 09:18

sobereva 发表于 2020-11-4 04:17
对于溶质来说，熵减小是普遍现象。吸附前单体可以自由运动，吸附后运动受限，因此熵减小才是常见的。
虽 ...

谢谢sob老师，请问静电作用一般如何对熵作出贡献？
我的理解是，假如两个物质之间因为有疏水作用、静电作用而吸附在一起了，物质A和B表面的水化膜会解缚成水分子回到溶液中，这部分熵是增加的，而两个分子的吸附是熵减小的，所以也很难说疏水作用、静电作用对熵的总体贡献是增大熵还是减小熵，不知是否正确？

作者
Author: a9471163 时间: 2020-11-4 09:20

sobereva 发表于 2020-11-4 04:17
对于溶质来说，熵减小是普遍现象。吸附前单体可以自由运动，吸附后运动受限，因此熵减小才是常见的。
虽 ...

假如物质疏水作用很强，物质表面的水会减少，所以能不能说，疏水作用强更有利于熵减小？

作者
Author: k64_cc 时间: 2020-11-4 12:23
本帖最后由 k64_cc 于 2020-11-4 12:51 编辑

a9471163 发表于 2020-11-3 09:20
你好，我是这样理解这个问题的：
ΔG=ΔU+PΔV-TΔS，目前ΔG=-40，ΔU=-100，PΔV项可以忽略不计，因为 ...

两个Energy Group之间的相互作用能不是dU，PMF曲线上两个state之间<U>的差值才是能和dG对比分析的dU。能量分解得到的相互作用能可能对自由能变化起主导作用，但是不能拿来和自由能曲线一起做定量分析。

作者
Author: sobereva 时间: 2020-11-4 14:34

a9471163 发表于 2020-11-4 09:20
假如物质疏水作用很强，物质表面的水会减少，所以能不能说，疏水作用强更有利于熵减小？

看是对谁来说，对溶质的话，从结果上是如此
溶质之间的静电作用体现在结合的焓变的层面上，不是熵的层面上

作者
Author: a9471163 时间: 2020-11-4 18:35

sobereva 发表于 2020-11-4 14:34
看是对谁来说，对溶质的话，从结果上是如此
溶质之间的静电作用体现在结合的焓变的层面上，不是熵的层面 ...

明白了，谢谢老师

作者
Author: a9471163 时间: 2020-11-4 18:44
本帖最后由 a9471163 于 2020-11-4 18:45 编辑

k64_cc 发表于 2020-11-4 12:23
两个Energy Group之间的相互作用能不是dU，PMF曲线上两个state之间的差值才是能和dG对比分析的dU。能量分 ...

谢谢您的提醒~
请问PMF曲线上两个state之间<U>的差值如何求出来，我目前想的办法是把state1和state2的分子构型单独拿出来，用gromacs定义能量组rerun再求一遍静电能+非静电能，然后两个state的能量作差，不知道这种算法是不是合适啊

作者
Author: k64_cc 时间: 2020-11-5 11:58
本帖最后由 k64_cc 于 2020-11-5 17:09 编辑

a9471163 发表于 2020-11-4 18:44
谢谢您的提醒~
请问PMF曲线上两个state之间的差值如何求出来，我目前想的办法是把state1和state2的分子 ...

不行，能量收敛比PMF慢，别想rerun了。而且G=H-TS的H是总势能，不是分解的，不用定义组，看总的Potential Energy就行。另外，水和溶质的相互作用能也会变化，中和一部分AB间作用能的变化，所以只分析这AB间相互作用也不是很充分。

饭得一口一口吃，你得先确定H主导了自由能变化，再分析什么因素主导了H变化。对于<H>的计算，你就在CV曲线上选两个点，用非常强的restraint把CV控制在那个点附近，然后跑三组长时间的模拟，算每组的平均势能，并用三组数据算variance。
Umbrella Sampling的原理是打格子数sample。你可以用一个模拟给多个格子提供数据。但是算PMF上某点平均势能的时候，你需要要尽量让sample处于你选定的那个格子里，所以要用一个很强的restraint来避免CV的波动。也是由于同个原因，US时的数据不能拿来直接用。

作者
Author: a9471163 时间: 2020-11-6 15:35

k64_cc 发表于 2020-11-5 11:58
不行，能量收敛比PMF慢，别想rerun了。而且G=H-TS的H是总势能，不是分解的，不用定义组，看总的Potential ...

我明白了，您讲得很清楚，非常感谢老师

作者
Author: a9471163 时间: 2020-11-6 15:36

k64_cc 发表于 2020-11-5 11:58
不行，能量收敛比PMF慢，别想rerun了。而且G=H-TS的H是总势能，不是分解的，不用定义组，看总的Potential ...

我想用您说的这个方法计算熵变，请问有没有用过这种方法的参考文献可以引用

作者
Author: k64_cc 时间: 2020-11-7 17:52

a9471163 发表于 2020-11-6 15:36
我想用您说的这个方法计算熵变，请问有没有用过这种方法的参考文献可以引用

随便整一篇我老板从前的工作：

https://pubs.acs.org/doi/abs/10.1021/jacs.5b07232

作者
Author: Tifo 时间: 2022-1-10 21:10

k64_cc 发表于 2020-11-5 11:58
不行，能量收敛比PMF慢，别想rerun了。而且G=H-TS的H是总势能，不是分解的，不用定义组，看总的Potential ...

老师您好，我也用和您说的一样的方法计算系统的总势能，就是在多个吸附离子中选择一个使用restraint将其控制在CV路径上的某一点3ns，再将获得的系统势能平均，但这样最后获得的势能随距离的变化曲线波动很大，与该离子的自由能变化对不上，请问是还有其他措施我没做吗？

作者
Author: k64_cc 时间: 2022-1-12 13:17

Tifo 发表于 2022-1-10 21:10
老师您好，我也用和您说的一样的方法计算系统的总势能，就是在多个吸附离子中选择一个使用restraint将其 ...

Variance大是一码事，和自由能曲线对不上是另一码事。

3 ns模拟，variance大是很正常的。一般来说，平均势能的收敛需要更长的时间。我之前做水的PMF，已经很好弛豫了，还是需要几十ns的轨迹才能确保它收敛。如果收敛之后势能曲线还是和自由能对不上，那就说明是熵驱动的，不是焓驱动。这俩玩意也不一定非要对上。

作者
Author: Tifo 时间: 2022-1-13 15:44

k64_cc 发表于 2022-1-12 13:17
Variance大是一码事，和自由能曲线对不上是另一码事。

3 ns模拟，variance大是很正常的。一般来说，平 ...

感谢老师回复，我还有个问题想请教一下。就是我现在以一个离子为探针，把它从吸附稳定位置拉到体相位置，那等这个探针离子到体相位置的时候，其他离子难道不会重新吸附到之前探针离子所在的位点上吗，那这时候这个系统的总势能为什么还会和我拉那个探针之前的时候不一样？

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