分子吸附材料表面，可能的构象太多，求解！

MrD

研究分子吸附一个材料表面问题，表面是周期性的，分子是不规则的，因此可能的吸附形态太多。

而几何优化的时候，分子不会大幅度翻转或平移，只是在给定初始结构附近，在Z方向上有变化，这样就只能人为做很多种初始构型。

请问有没有更好的方法或算法，能让分子在表面驰豫幅度大一些，尽量找到一个能量最小点。

尝试过第一性原理的动力学，但由于体系大，80多原子，即使降低到最差精度，一天也只能几百步（每步一飞秒），4核8GB的win7机器表示很骄傲。

另外想请教，如果想做第一原理动力学的话，想对这种结构跑平衡，恒温器有没有太大影响，或者使用恒温器应根据什么样的经验来选择呢？时间尺度多少事宜呢？是否还要模拟退火？

我不太懂动力学，感觉在跑平衡的过程中，势能应该是降低然后趋于平衡，但我跑出的势能为什么是波浪线呢？（我用的NVT系宗，300K，simple NH恒温器）


PS：我用的是MS-Dmol, 有人说dmol优化周期性体系不给力，让转VASP，但我始终觉得只要算法完善，不应该有程序不给力的一说。请各位指点。

Aristotler

sobereva 发表于 2017-10-28 12:32
原理上会，但是你需要多跑的时间以数量级来论

谢谢老师

sobereva

Aristotler 发表于 2017-10-28 07:53
谢谢老师，请问老师这个加速温度如何选取呢？

看具体情况，并且需要测试，比如总不能高到整个体系都完全散架，但又不能太低，否则起不到效果

sobereva

Aristotler 发表于 2017-10-28 08:01
假如只加热300K，但时间够长，是不是也会和高温下达到的效果差不多？得到一些高温下才容易产生的结构？谢 ...

原理上会，但是你需要多跑的时间以数量级来论

Aristotler

Aristotler 发表于 2017-10-28 07:53
谢谢老师，请问老师这个加速温度如何选取呢？

假如只加热300K，但时间够长，是不是也会和高温下达到的效果差不多？得到一些高温下才容易产生的结构？谢谢

Aristotler

sobereva 发表于 2017-10-28 01:32
你的具体情况我不清楚。不同的外环境都会影响柔性分子的最稳定构象，真空下和在晶体环境下的情况有区别 ...

谢谢老师，请问老师这个加速温度如何选取呢？

sobereva

Aristotler 发表于 2017-10-27 22:47
老师，看到这个帖子，学生有一个迷惑的事情。为了得到合理的能量最低的构型，我先采用淬火的方法在1000K ...

你的具体情况我不清楚。不同的外环境都会影响柔性分子的最稳定构象，真空下和在晶体环境下的情况有区别很正常。
虽然显示世界中不会加热到1000K去，但是你做模拟升温到1000K再退火，这是个加速采样的过程，用比较少的时间模拟就可以找到可能需要非常长时间的低温下MD（通过翻越势垒）才能得到的真正低能的构象。

Aristotler

sobereva 发表于 2015-9-2 23:41
很可能这样，所以不能光取构象搜索时能量最低的那个，而是最低的一批都应当考虑，看谁吸附后能量最低。

老师，看到这个帖子，学生有一个迷惑的事情。为了得到合理的能量最低的构型，我先采用淬火的方法在1000K下模拟100ps, 得到了三种能量梯度的构型（这个可能有些不专业，就是用forcity 淬火后优化，100000步中三种大概率的结构，分别是583kcal/mol, 587kcal/mol 和591kcal/mol。之后拿去优化，还是591kca/mol的能量最低。但是我最近找到了这个物质（C4Q)的晶体结构，单个分子的结构和587kcal/mol类似。学生疑惑的是：分子动力学的意义是什么？给了一个很高的温度，为了分子有一个初始的动能，可以越过一些能垒，找到能量最低（最稳定？）的构型，但实际中也许这个物质就不会受到这个动能，不会变到类似591kcal/mol 的结构。那么如何搜索构想？如果不给一定高的温度，结构的多样性又会减少，如何判断加热温度？学生刚刚接触动力学，有不对的地方请您指出，感谢。

卡开发发

鹃子 发表于 2016-11-28 16:49
非常清晰，谢谢大神！
1. 对于4还有一事不明，既然Mulitipole expansion的值在Dmol3中最大是16，一般默认 ...

1、有机系统一般原子的角动量都不会很高，所以差异不大，所以我个人觉得不必重新算。但过渡金属就得小心了，能量有时候影响不小，力也会有影响。
2、增加精度之后并非会更接近实验值，随着精度增加到一定程度会接近理论极限值，在理论框架确定的情况下要得到的也是这个极限值才是，如果理论极限值和实验值差异太大，那就是理论方法选取的不好。单纯的提高精度无意义，必然伴随着收敛测试及其标准。
3、截断半径和多极展开并非像基组尺寸那样是因为变分原理限制导致的精度差异，因此这两者的精度提升能量不见得是下降的，或者说之前能量低仅仅只是数值误差。

teller3531

精彩的帖子，学习了。我做了点石墨烯吸附苯，不过没用周期性，我就去了一片，用packmol得到的初始结构都是π-stacking，不论初始结构怎么放，packmol跑了之后都是π吸附

鹃子

非常清晰，谢谢大神！
1. 对于4还有一事不明，既然Mulitipole expansion的值在Dmol3中最大是16，一般默认设置取到8，记得之前对于精度的设置Sob老师提过一般精度的设置对于几何优化和频率计算的结果影响很小，是否可以在opt 和freq计算中采用默认，只对单点进行优化设置？
2. 依然采用TBHP做测试，Mulitipole expasion 设置为16，对于计算量的影响确实很小，几何优化一步即可完成，但是能量相对于8的增高约0.01eV，如前面cutoff从4.3-4.6（能量增大0.096eV），增大设置精度之后能量都有所增大，我的理解是计算精度越高，能量越接近真实值，其绝对值越大，即能量越负才对。如果是这样，单纯地提高计算精度是否有意义？
PS：对于相对提高计算精度能量的增大有点迷糊了，还请赐教!

鹃子

非常清晰，谢谢大神！
1. 对于4还有一事不明，既然Mulitipole expansion的值在Dmol3中最大是16，一般默认设置取到8，记得之前对于精度的设置Sob老师提过一般精度的设置对于几何优化和频率计算的结果影响很小，是否可以在opt 和freq计算中采用默认，只对单点进行优化设置？
2. 依然采用TBHP做测试，Mulitipole expasion 设置为16，对于计算量的影响确实很小，几何优化一步即可完成，但是能量相对于8的增高约0.01eV，如前面cutoff从4.3-4.6（能量增大0.096eV），增大设置精度之后能量都有所增大，我的理解是计算精度越高，能量越接近真实值，其绝对值越大，即能量越负才对。如果是这样，单纯地提高计算精度是否有意义？
PS：对于相对提高计算精度能量的增大有点迷糊了，还请赐教!

卡开发发

鹃子 发表于 2016-11-28 11:04
楼楼们的讨论确实很精彩，小菜也是一直在计算中只将cutoff按照Coarse、Medium、Fine等设置，今天看了@卡开 ...

1、一般Fine以上的级别，能量对截断半径的偏差是0.1eV/Ang，如果是比较相对能量的话，这点误差估计就很小了，所以大部分情况选取Fine的精度就能够很好。
2、1中提到的误差必然存在，如果要尽可能抵消那个误差，截断半径最好是一致的，这也是最合理的做法，否则就会引入新的误差。不少人确实只知道设置精度，不知道手工修改截断半径来保证同一个体系中的截断半径一致，所以算出来相互作用能变成正的事情还是时常发生。
3、那么，既然每个原子都有最佳的截断半径，是否可以在同一个计算中单独给每一个原子设定一个合理的截断半径以节约计算量呢？答案是可以，但是我并不推荐，我自己测试的经验来看（这当中其实做了不少测试，尤其是以前几个算CeO2的朋友其实也都希望能只对Ce设置那么大的截断半径，其他的设置小点），这种做法给SCF带来很大的数值不稳定，不是SCF难收敛就是结构优化不好做，具体原因我不是很清楚。
4、Multipole expansion是一种计算Coulomb势的方案。这个过程会把单中心电荷密度做多极展开（也就是按照不同l、m的球谐函数的级数），然后求得相应的Coulomb势的贡献，然后再对这些势进行求和。m只能取到-l~l，而l其实原则上是0~∞，所以原则上应该是要测试决定l一般是多少，必然也是越大越合理。DMol3最大做到16极，不过即便做到16极计算量也不会比8极大太多。尤其是包含角量子数比较高的原子的系统，我建议还是设置大点好，比如d、f电子系统。

鹃子

楼楼们的讨论确实很精彩，小菜也是一直在计算中只将cutoff按照Coarse、Medium、Fine等设置，今天看了@卡开发发的解答，查看我的体系，发现不同原子种类的结构中cutoff确实差异不小。于是乎，将一个叔丁基过氧化氢（TBHP）的分子按照体系中最大的global cutoff重新优化，发现cutoff从3.3-4.3，TBHP能量降低约0.381eV，但是再按照Fine将cutoff设置为4.6之后，TBHP能量相对于cutoff4.3升高约0.096eV（只做几何优化，并未计算单点能）。
这样的结果是不是意味着不同的原子种类确实对应有相应的最适宜的cutoff值，但是为了计算相互作用，我们必须采用相同的参考标准结果才有意义？
PS:看到卡开发发大神在一层提到Multipole expansion的设置，MrD后续的回答中也提到自己的多级展开是Octupole，对比我自己的体系，发现软件给的参数也是Octupole，不知道Dmol3根据什么设置这个参数?(Dmol3说明文件中的Multipolar expansion定义为Specify the maximum angular momentum function used in the multipolar representation of the charge density. 不知道大神的意思是不是这个参数也因该保持一致？）

卡开发发

jiewei 发表于 2015-9-4 13:48
可能是我表达有误， 能量最低肯定是有原因， 不过在某些条件下的分子构象未必是以能量最低的形态存在。  ...

一般来说应该很多情况都是polycrystal，如果是amorphous计算就真的没办法按照表面去做了。polycrystal就可以按照上次的想法，实际上每个晶粒及各个暴露面还是很大的，足以使用Slab模型进行近似。原则来说确实像兄台所说的，甚至不同的面的不同的吸附方式也得考虑。amorphous有时候用一些disorder处理（比如VCA）还凑合，但也很难说可不可靠。

当然，确实有不少相当复杂的情况其实理论非常难于研究，理论不是万能的。顺带吐槽一个：有时候这样的体系老师不了解计算还想让学生做，学生自己也没判断能力来判断这样的体系能不能通过ab initio的方式实现，这种情况非常让人头疼。