请选择 进入手机版 | 继续访问电脑版

计算化学公社

 找回密码
 现在注册!
查看: 17572|回复: 66

[molclus] genmer:生成团簇初始构型结合molclus做团簇结构搜索的超便捷工具

  [复制链接]

1万

帖子

25

威望

2万

eV
积分
41955

管理员

公社社长

发表于 2015-12-15 12:18:02 | 显示全部楼层 |阅读模式
注:本文内容对应molclus官网上最新版molclus程序包里的genmer组件

genmer:生成团簇初始构型结合molclus做团簇结构搜索的超便捷工具

文/Sobereva@北京科音
First release:2015-Dec-15   Last update: 2019-May-26

1 前言

笔者开发过免费的molclus程序用于团簇构型搜索和分子构象搜索,不熟悉的话应当先看此文《使用molclus程序做团簇构型搜索和分子构象搜索》(http://bbs.keinsci.com/thread-577-1-1.html)。简单来说,molclus的主要执行过程就是读取含有多帧结构信息的traj.xyz文件,依次对每个结构调用Gaussian/MOPAC/ORCA/xtb/OpenBabel进行优化,然后对结果用isostat子程序进行统计分析。

记录初始构型的traj.xyz可以通过多种方式产生,此文介绍用笔者开发的genmer程序来快速方便地生成molclus所需要的初始构型,读者会发现基于genmer来使用molclus做原子和分子团簇构型搜索超级方便和高效!!!获得genmer的方法就是到molclus主页(http://www.keinsci.com/research/molclus.html)下载最新版本molclus程序,预编译的genmer就在其中。使用genmer/molclus发表文章时请按照molclus主页的说明进行引用。

目前已经有不少采用genmer+molclus的团簇构型搜索文章发表了,比如J. Phys. Chem. A, 120, 4560 (2016)、J. Theor. Comput. Chem., 16, 1750058 (2017)、Phys. Chem. Chem. Phys., 19, 28434 (2017)、Chem. Phys. Lett., 663, 128 (2016)等。


2 原理

genmer可以把任意种、任意数目的单体(可以是原子或者分子)放在一起通过随机方式生成指定数目的团簇初始构型。

genmer生成团簇初始构型的算法是洒家想出来的,又简单又实用,耗时又极低。就是让被加入的单体的中心从体系的原点出发向着一个随机的方向“发射”,发射过程中单体不断移动,当这个单体被移动到和体系已有的原子没有不合理接触(即被加入单体里的原子和当前团簇的原子间的距离都大于相应两个原子的半径和),这个单体就被加入当前团簇;接下来再这样加下一个单体,反复如此。最开始时的团簇就是第一个单体,随着单体不断加入,团簇就会不断长大,当指定数目的单体全都被加进去时这个团簇的生成就宣告结束。由于每个单体发射的方向都是随机的,还可以允许单体自身随机旋转,所以每次用这个算法生成的团簇构型都可以不一样。在生成团簇的时候可以通过设置一些约束条件来控制产生的团簇的大致形状,比如可以设置单体原子允许的Z坐标范围(通过zmax、zmin指定),再比如可以设置单体中心距离体系中心的最大和最小距离(通过rmax、rmin指定),等等。如果沿着某个方向发射单体过程中,在单体移动到合适的加入位置之前就已经违背了约束条件,那么这次发射就失败了,程序会再尝试向着另一个随机朝向发射,反复如此,最大尝试的次数可以通过maxemit指定。如果超过这个尝试次数还没能成功地将单体加入体系,那么这个团簇的产生就宣告失败。

将上述团簇产生过程画图示意的话就是下面这样。青色的椭圆是已经加入团簇的单体,三个箭头示意了三种随机方向,当要被加入的单体按照箭头所示方向运动到图中黄色椭圆的位置时,由于和团簇里现有的原子已经没有不合理的接触了,此时这个单体就成功加入了,其坐标就被写入到团簇坐标里了。
emit.png

genmer的这个算法既可以生成分子团簇的初始构型,也可以生成原子团簇的初始构型。对于前者,单体是分子,原子半径使用范德华半径;对于后者,单体是原子,原子半径使用共价半径。

假设设定了团簇由10个A、5个B单体构成,程序会先把10个A都加进去,再把5个B加进去,根据如上算法,得到的团簇中10个A就会倾向于在团簇内层,5个B会倾向于在团簇外层。如果想让它们的出现区域不受设定顺序的影响而是完全随机的,genmer也可以做到自动打乱顺序,这称为shuffle。

在产生团簇时每个分子都被当成是刚性的。但如果一种分子有多种可能的构象都需要考虑,可以提供包含多个构象结构的分子的.xyz文件,并且设定每种构象出现的概率。这样每次加入分子时就会按概率抽取一种构象。

虽然genmer产生的团簇初始构型可能看起来略微松散,但这完全没关系,一做优化就自动收缩变得致密了。molclus用的团簇构型搜索的关键是让初始构型分布得足够随机且全面,这样才可能经过优化后获得尽可能多的极小点构型,特别是采样到能量最低结构。只要让genmer产生的构型数目足够多,是肯定可以达到这个目的的。


3 使用方法

genmer的使用非常容易,只需要提供每种分子的几何结构.xyz文件,设定出现数目就行了。如果你不懂什么叫.xyz格式,看此文:《谈谈记录化学体系结构的xyz文件》(http://sobereva.com/477)。运行参数和分子设定都记录在genmer.ini文件中,直接启动genmer后程序会读取当前目录下的genmer.ini,然后把生成的指定数目的初始构型都输出到当前目录下的traj.xyz当中。之后就可以直接运行molclus来自动调用其它程序去优化其中的每个结构了。对于Windows版,双击genmer.exe图标即可启动,对于Linux版,通过输入genmer可执行文件的路径即可启动。

参数文件genmer.ini的每个选项在此文件里都有注释,一看就懂,下面介绍一下。通常大家就在程序自带的genmer.ini基础上编辑即可,如果此文件里有的参数没有出现,就会用默认值。

ngen= 50   // 期望生成的团簇构型的数目。默认为1
ngenmust= 0   // 若>0,则ngen选项无效,genmer将无限循环不断产生构型,直到成功产生的构型达到ngenmust。此选项默认为0,代表不做这个要求
step= 0.1   // 设定每一步移动单体的距离。对于单体数少的情况不用改这个,但如果要产生的团簇里单体数甚巨,想节省时间可以把这个步长值设大,但出现不合理接触的可能性会增加。默认为0.1,一般不需要改
ishuffle= 1    // 0: 不打乱加入单体的顺序  1: 打乱顺序。如果当前只有一种单体的话,此选项设0还是1都不影响结果。默认为0
imoveclust= 1    // 1: 将产生的团簇的中心挪到原点,即(0,0,0)的位置  0: 不挪动。默认为1
irandom= 1   // 0: 输入文件相同的情况下让每次运行结果都相同  1: 让每次运行结果都不同(因为此时每次运行时会使用不同的随机数种子)。默认为1
iradtype= 1    // 0: 使用CSD共价半径来检测原子间接触,适合生成原子团簇  1: 使用Bondi范德华半径检测原子间接触,适合生成分子团簇。默认为1
radscale= 1.0   // 原子半径乘的系数。显然设得越小,单体之间的原子间距离会越近;设得越大,单体之间的原子间距离会越大。默认为1.0,一般不需要改
maxemit= 300   // 生成每个团簇时尝试发射单体的最大次数。默认300,一般不用改。但当约束条件设得比较严、比较复杂而不容易成功产生团簇时,可以加大这个参数增加成功几率
----           这一行必须有且包含"----"字符,告诉程序代表接下来的内容是每类单体的数目(可以附带选项)和.xyz结构文件,单体类型可以有无数多个。
4    //第1类单体数目
examples/H2O.xyz   //第1类单体的结构文件。可以写绝对路径,也可以写这样的相对路径
3    //第2类单体数目
examples/ethanol.xyz   //第2类单体的结构文件
...略

在单体数目后头可以加入以下选项来设置加入此类单体时的特殊设定:
norot:单体不允许随机旋转。默认情况下每次加入单体时都是让单体随机旋转的,以让团簇构型间尽量有差异性
cenxyz [XYZ坐标]:定义单体的中心位置,坐标是相对于单体结构文件的坐标。默认情况下是以单体几何中心作为单体的中心位置
rmax [距离]:设定这个单体的中心与体系中心最远距离不能超过多少。默认不做限制
rmin [距离]:这个选项可以要求单体中心与体系中心距离不得小于多少,这实际上是发射单体时与体系中心的最初距离。默认为0
zmax [数值]:要求单体里的所有原子的z坐标不允许大于这个值。默认为不限制
zmin [数值]:要求单体里的所有原子的z坐标不允许小于这个值。默认为不限制
显然,这些选项对于精确控制团簇构型很有用,非常灵活。比如用cenxyz可以直接对某些单体进行定位;rmax可以令单体构成球形团簇;rmax和rmin联用可以实现构建空心团簇;zmax和zmin联用可以产生扁片状团簇;rmax、rmin、zmax、zmin一起联用可以生成环状或者桶装团簇。

以下是单体设置中加入上述选项时的格式示例
3 norot cenxyz 0 1.5 0.2
examples/CH3CH2CH2NO2.xyz
20 rmax 6
examples/H2O.xyz



4 用genmer产生团簇初始构型例子

下面给出几个genmer的例子展现其生成团簇的能力,用到的单体结构文件都是genmer自带的。为了简洁,基本上所有默认的参数都舍掉了。以下图像都是把genmer产生的traj.xyz文件直接拖进VMD程序(http://www.ks.uiuc.edu/Research/vmd/)显示的。分子单体的结构都是用GaussView建的,然后再按照http://sobereva.com/477里说的方式转换为xyz格式。

4.1 6个水团簇

用以下genmer.ini文件可以得到由6个水组成的20个团簇初始构型
ngen= 20
----
6
examples\H2O.xyz


结果如下:


6H2O.gif


4.2 3000个水团簇

用以下genmer.ini文件生成1个多达3000个水组成的团簇初始构型,并且水距离体系中心不能超过40埃。由于分子数较多,生成需要一定时间,大约十秒左右。注意rmax不能设得太小,否则会导致有的水在各个方向尝试加入均不成功,最后宣告失败(由于算法有一定随机性,可能这次没生成成功,但下次就成功了)。
ngen= 1
----
3000 rmax 40
examples\H2O.xyz


结果如下:
3000H2O.png


4.3 6个水和3个乙醇的团簇

用以下genmer.ini文件可生成6个水和3个乙醇构成的10个团簇,并且打乱乙醇和水的出现顺序
ngen= 10
ishuffle= 1
----
6
examples\H2O.xyz
3
examples\ethanol.xyz


结果如下:
6H2O 3ethanol.gif

4.4 顺式和反式1,3丁-二烯混合团簇

看过上面的例子,肯定知道怎么获得顺式和反式1,3-丁二烯的混合团簇了,只要提供一个顺式的.xyz文件,一个反式的.xyz文件,设定好各自的数目即可。但如前所述,也可以把多构象的分子的各个构象合一起作为单个.xyz文件提供。这里产生一个由30个1,3-丁二烯构成的团簇,并让顺式、反式出现概率为30%和70%。

1,3-丁二烯.xyz文件结构如下
10
0.3 cis   //顺式的出现比例为0.3(30%)
[顺式的坐标]
10
0.7 trans   //反式的出现比例为0.7(70%)
[反式的坐标]

genmer.ini内容为
ngen= 1
----
30
examples\butadiene_cis_trans.xyz


结果如下:
30butadiene_cis_trans.png

4.5 多巴胺二聚体

用以下genmer.ini文件可以生成15个多巴胺二聚体初始构型,这算是稍微大一些的分子了。
ngen= 15
----
2
examples\dopamine.xyz


结果如下:
dopamine_dimer.gif

4.6 水层包围多巴胺

用以下genmer.ini文件可以生成200个水分子紧密包裹在多巴胺附近的10种构型。由于ishuffle设为了0(其实这本身也是默认的),而且多巴胺的结构文件先出现,因此多巴胺被包围在内侧。水的rmax不能太大也不能太小,为了产生接近球形且较为致密的包围多巴胺的水球,可以先尝试比较大的rmax,然后逐渐缩小,直到程序无法成功产生团簇结构为止。
ngen= 10
ishuffle= 0
----
1
examples\dopamine.xyz
200 rmax 15
examples\H2O.xyz


结果如下:
dopamine_water.gif

4.7 水分子与多巴胺的氨基形成氢键二聚体

以下genmer.ini生成一个水与一个多巴胺形成的二聚体,共生成20个,要求水必须与多巴胺的氨基接触,目的是优化后可以得到水与氨基通过氢键形成的二聚体。
ngen= 20
imoveclust= 0
----
1 norot cenxyz 4.33709192    0.16551438    0.29362167
examples\dopamine.xyz
1 rmax 4
examples\H2O.xyz


注意这里对多巴胺用了norot,这样它在团簇里的朝向就和单体结构文件里相同了,而不会随机变化朝向。cenxyz后面的值正是dopamine.xyz里氨基的氮原子的坐标。由于多巴胺是第一个出现的分子,而且其中心被设定在了氮原子上,因此在产生的团簇里,多巴胺的氮原子恰好就是在团簇坐标的原点的位置。由于rmax 4让水的中心(默认是水的几何中心)不能与团簇原点距离超过4埃,所以保证了产生的团簇中水肯定在氨基附近。此例用了imoveclust= 0,使得生成的团簇的笛卡尔坐标原本是什么就在traj.xyz文件里记录什么,而自动不把整个团簇的几何中心在写入traj.xyz之前挪到(0,0,0)的位置,这使得traj.xyz中多巴胺的位置是固定的。如果你把这个例子弄懂了,就差不多算是完全理解genmer的算法了。

运行结果如下。为了看得清楚,这里把所有20个团簇结构叠加在了一起显示,可见水很好地包围了氨基,各个角度都有出现。如果你打开traj.xyz,会看到氮原子的坐标恰好是(0,0,0)。

1dopamine water_NH2.png

4.8 Li6团簇

用以下genmer.ini文件可以生成6个Li组成的团簇的15个初始构型。注意这里用了共价半径
ngen= 15
iradtype= 0
----
6
examples\Li.xyz


结果如下:
Li6.gif

4.9 环状的Li30团簇

用以下genmer.ini文件可以生成10个由30个Li组成的环状团簇的初始结构。这些Li只允许出现在径向距离r=5埃至9埃的范围内,并且要求z坐标必须在-1~1埃之间,因此得到的就是螺丝垫片形状的团簇了。由于当前设的约束条件比较多,Li只能出现在相对比较狭小的限定的空间中,为了增加产生的成功几率,这里将尝试发射单体次数maxemit设的比默认值300大得多得多(凡是像这种对单体出现的空间区域要求的比较严的情况,都建议改大maxemit。如果设得很大还是不能成功产生,那要么减少加入的单体数,要么把对空间范围的约束放宽点)
ngen= 10
iradtype= 0

maxemit= 1000
----
30 zmin -1 zmax 1 rmin 5 rmax 9
examples/Li.xyz


结果如下(为了看着更清楚,在VMD里显示的时候令间距小于特定距离的Li-Li显示出了键)
Li30_ring.gif

4.10 笼状的碳团簇包夹水分子

用以下genmer.ini文件可以生成100个碳组成的笼状团簇,其中包夹一个水分子。由于当前允许碳出现的空间范围很窄,即r=4~5埃,所以maxemit设得较大,但即便设为3000,还是有个别团簇产生失败,当前要求产生10个,实际成功产生了7个。如果把maxemit设得更大,成功产生的会更多。这里C.xyz是只含一个碳原子的xyz文件。
ngen= 10
iradtype= 0
maxemit= 3000
----
100 rmin 4 rmax 5
C.xyz
1
examples/H2O.xyz

结果如下:
C100+H2O.gif

4.11 CB6主体分子包夹客体分子乙醇

用以下genmer.ini文件可以生成环状主体分子CB6包夹乙醇分子的10个构型。CB6的中心被要求固定在0,0,0位置不动。由于CB6里的孔洞很窄,为了让生成容易成功,这里将原子半径设为了默认的0.8倍。由于乙醇如果不出现在主体分子孔洞内的话就没意义了,所以用了rmax 2让它的中心偏离CB6中心(即0,0,0)不得超过2埃,否则有可能会出现乙醇结合到CB6外部的结构。
  ngen= 10
  radscale= 0.8
----
1 norot cenxyz 0 0 0
CB6.xyz
1 rmax 2
examples/ethanol.xyz


结果如下
CB6-ethanol.gif

5 结合molclus做团簇构型搜索

上一节演示了用genmer得到了各种团簇的traj.xyz文件,这一节我们基于这些团簇初始构型,用molclus做构型的批量优化、获得能量,从而筛选出较能量较低,也因此较合理的结构。如果下文的例子没看懂,应当先把前述的《使用molclus程序做团簇构型搜索和分子构象搜索》(http://bbs.keinsci.com/thread-577-1-1.html)一文仔细看一遍了解molclus的基本用法。由于genmer每次产生的团簇构型可能不同(取决于irandom),因此你尝试重复下文的例子的话,结果可能和下文有所不同。

5.1 多巴胺二聚体

这里我们用molclus调用MOPAC在PM6-DH+级别下优化多巴胺二聚体来得到它们的团簇构型。PM6-DH+是一种适合计算弱相互作用的半经验方法,虽然定量结果算不上精确,但是对弱相互作用体系给出定性正确的结果是足够的。

把前面4.5节得到的含有15个多巴胺二聚初始构型的traj.xyz文件放到molclus目录下,把molclus的settings.ini里的MOPAC路径配置好,将iprog设为2,并确认MOPAC的模板输入文件template.mop里的关键词是PM6-DH+ precise。启动molclus,程序自动调用MOPAC对15个结构依次优化,很快就算完了,结果输出到当前目录下的isomer.xyz里。然后启动isostat,连敲回车三次,即看到按照由低到高排序的能量值,这里只列出前5个:
#    1 Count:    1  E=     -0.243275 a.u.  DGmin=   0.28  DE=    0.00 kcal/mol
#    2 Count:    1  E=     -0.242410 a.u.  DGmin=   0.26  DE=    0.54 kcal/mol
#    3 Count:    1  E=     -0.241462 a.u.  DGmin=   0.40  DE=    1.14 kcal/mol
#    4 Count:    1  E=     -0.240707 a.u.  DGmin=   0.27  DE=    1.61 kcal/mol
#    5 Count:    1  E=     -0.239548 a.u.  DGmin=   0.36  DE=    2.34 kcal/mol

然后可以按回车退出,也可以输入温度,比如300,得到各种构象的波尔兹曼分布比例(如果不知道原理,看此文http://sobereva.com/165。isostat是把单点能近似作为自由能代入公式计算的)
#    1   Count:   1   Ratio:  60.30%
#    2   Count:   1   Ratio:  24.28%
#    3   Count:   1   Ratio:   8.95%
#    4   Count:   1   Ratio:   4.04%
#    5   Count:   1   Ratio:   1.19%


按照能量排序后的结构都输出到了当前目录下的cluster.xyz当中。拖入VMD观看,前5个结构图如下,与上面列出的5个依次对应


1.png

为什么第1个结构能量很低是很显然的,因为两个氨基都各自和两个羟基同时形成了氢键。第2个结构是每个氨基和一个羟基形成氢键。第3个结构形成了明显的pi-pi堆积,而且两个氨基之间有氢键。第4个结构有两对氨基与羟基形成的氢键,而且还有一点pi氢键的特征。第5个结构只有一对氢键,但是pi氢键的特征更明显。

由于本例只考虑了15个多巴胺二聚体结构,不算很多,而二聚体极小点构型有很多,还不能确保二聚体最稳定构型就在其中,要想更稳妥应该考虑更多结构。另外,半经验方法精度终究比较有限,获得能量接近的结构的能量差定量不准却,甚至顺序也不对。如果要得到准确的二聚体结构和能量,建议在PM6-DH+优化完之后用molclus调用Gaussian在B3LYP-D3(BJ)/6-311G*级别下对各个结构再优化一下,然后用M06-2X-D3/6-311+G(2d,p)级别(要求更精确可以用B2PLYP-D3(BJ)/jun-cc-pVTZ)计算能量最低的数个结构的单点能。这种“半经验预优化->DFT结合一般基组优化->高精度方法算单点能"是我很推荐的搜索分子团簇的做法,这用molclus实现极其容易,可参考《使用molclus程序做团簇构型搜索和分子构象搜索》中计算水四聚体团簇构型的例子。


5.2 乙醇+水团簇

这一节用Gaussian里的UFF力场优化4.3节得到的6个水+3个乙醇的团簇的10个初始构型。力场计算极其便宜,可以比半经验优化更大体系、更大批量构型。(由于MOPAC和xtb速度很快,其实对于当前这样不大的体系,Gaussian做UFF计算比之快不了多少)。

把4.3节的traj.xyz拷到molclus目录下,把template.gjf里的关键词设为UFF=qeq opt,然后确认settings.ini文件里iprog设为了1。启动molclus,约一分钟就把所有构型都优化完了。用isostat统计结果如下
#    1 Count:    1  E=     -0.063488 a.u.  DGmin=   0.50  DE=    0.00 kcal/mol
#    2 Count:    1  E=     -0.063090 a.u.  DGmin=   0.38  DE=    0.25 kcal/mol
#    3 Count:    1  E=     -0.062488 a.u.  DGmin=   0.50  DE=    0.63 kcal/mol
#    4 Count:    1  E=     -0.062053 a.u.  DGmin=   0.38  DE=    0.90 kcal/mol
#    5 Count:    1  E=     -0.058900 a.u.  DGmin=   0.34  DE=    2.88 kcal/mol
#    6 Count:    1  E=     -0.058547 a.u.  DGmin=   0.34  DE=    3.10 kcal/mol
#    7 Count:    1  E=     -0.057159 a.u.  DGmin=   0.53  DE=    3.97 kcal/mol
#    8 Count:    1  E=     -0.056104 a.u.  DGmin=   1.38  DE=    4.63 kcal/mol
#    9 Count:    1  E=     -0.055407 a.u.  DGmin=   2.19  DE=    5.07 kcal/mol
#   10 Count:    1  E=     -0.049182 a.u.  DGmin=   0.39  DE=    8.98 kcal/mol


300K下的分布比例
#    1   Count:   1   Ratio:  44.57%
#    2   Count:   1   Ratio:  29.34%
#    3   Count:   1   Ratio:  15.56%
#    4   Count:   1   Ratio:   9.85%
#    5   Count:   1   Ratio:   0.36%
#    6   Count:   1   Ratio:   0.25%
#    7   Count:   1   Ratio:   0.06%
#    8   Count:   1   Ratio:   0.02%
#    9   Count:   1   Ratio:   0.01%
#   10   Count:   1   Ratio:   0.00%


按能量排序后的结构如下。和4.3节中genmer初始产生的结构对比,可看到经过力场优化后团簇变得明显紧凑了。之后我们可以再让molclus调用量子化学程序进一步优化、计算能量获得更准确的结果。

6H2O_3ethanol_UFF.gif


5.3 Li6团簇

这一节我们用molclus调用Gaussian在B3LYP/6-31G*级别下寻找Li6团簇稳定构型。

把4.7节得到的包含15个Li6团簇初猜构型的traj.xyz放到molclus目录下。将Gaussian的模板输入文件template.gjf里的关键词填入B3LYP/6-31G* opt(cartesian,maxcyc=50),备用模板文件template2.gjf里的关键词填入B3LYP/6-31G* opt(cartesian,maxstep=5,notrust,maxcyc=50,gdiis)。将iprog设为1,igaucontinue设为1。含义就是说,对每个结构先用第一个模板文件的设定优化,如果不收敛或报错,则取最后结构用第二个模板的设定继续优化(若不明白为什么这么用关键词看此文http://sobereva.com/164)。模板文件步数改为50是因为对当前体系Gaussian默认的步数上限太少,所以适当增大。优化过程中可能会由于构成键角、二面角的冗余内坐标原子间排成直线而报错,所以这里用笛卡尔坐标优化来避免。

启动molclus,程序对15个结构依次优化,结果输出到isomer.xyz里。然后启动isostat,连敲回车三次,看到有4种唯一的结构
#    1 Count:    2  E=    -45.110127 a.u.  DGmin=   0.02  DE=    0.00 kcal/mol
#    2 Count:    3  E=    -45.109176 a.u.  DGmin=   0.02  DE=    0.60 kcal/mol
#    3 Count:    7  E=    -45.105556 a.u.  DGmin=   1.55  DE=    2.87 kcal/mol
#    4 Count:    2  E=    -45.103271 a.u.  DGmin=   2.03  DE=    4.30 kcal/mol

相应的结构输出到了cluster.xyz中。从结构看会发现前两个结构其实是一样的,之所以能量有细微差别其实是收敛精度问题而已。所以实际上有三个唯一的团簇结构,如下所示,从左到右能量依次升高

2.png

这些团簇从图上看都有对称性,而当前优化过程由于收敛精度原因可能最后一步的构型不严格满足实际对称性,故建议自行在gview里把结构做对称化后再优化一下以得到更准确结果(把对称化后的坐标写成traj.xyz的格式,然后用molclus再跑一次就行了)。


5.4 环状Li30团簇

在4.9节中提到了怎么产生环状Li30的初始构型。xtb程序(基本使用介绍见http://sobereva.com/421)速度极快,耗时比MOPAC还低,而且xtb程序实现的GFN-xTB方法比MOPAC支持的那些半经验方法普适性好得多,因此明显更适合处理电子结构相对来说比较复杂的原子团簇体系。为了对Li30构型进行快速搜索,在这一节我们用molclus调用xtb对已产生好的那10个Li30初始结构进行优化。之后,如果你想要更准确、更可靠的结果,可以再调用Gaussian或ORCA在更准确的DFT下优化和算能量。

把molclus的settings.ini里的iprog设为4,并且在xtb_arg选项里加上--gfn 1,这代表让xtb使用GFN-xTB方法的1代。这点很重要,目前的(2019-Mar-18发布的)的xtb程序默认用的是GFN2-xTB,即GFN-xTB方法的2代,虽然比1代在原理上更好,但是碰见复杂的原子团簇很容易出现SCF不收敛,因此得强制要求用1代。

启动molclus,在普通4核机子下xtb仅仅几秒钟就可以优化完一个结构,xtb速度超级快!!!然后我们照常用isostat对结果排序,其中能量最低的一个结构如下所示

Li30_xtb.png

这个结构非常理想,特别有序,凭直觉就能知道应该是个很稳定、很有意义的结构。此例反映出molclus+xtb的组合对于原子数较多的原子团簇的搜索也很合用,可谓黄金组合!

上面是研究Li30的环状团簇,笔者还用molclus结合xtb程序对30个Li构成的笼状团簇做了构型搜索,genmer产生的团簇初始结构和优化后的结构如下所示,细节详见http://bbs.keinsci.com/thread-12509-1-1.html,如图可见得到的结构意想不到地好!genmer一开始给出的初始结构比较乱,但是经过神器xtb一优化,就变得非常有序、漂亮了。

3.png


6 总结

从以上介绍和例子可见,用genmer+molclus的方式做团簇构型搜索虽然思路简单,没有用到什么高深的技术,但是非常实用,结果理想,操作很简便,对使用者的知识要求很低,完全不用掌握什么复杂的理论知识,强烈推荐大家在实际研究中使用,并且向同行推广!

有一些用户问让genmer产生多少初猜团簇比较合适,这取决于体系大小、特征、研究目的、可接受的计算花费。如果你想得到能量最低的团簇结构,那么显然产生的初猜结构数目越多,最终能获得能量最低结构的几率就越大。对于结构比较简单的小分子构成的二聚体团簇,想获得能量最低的构型,一般产生20个初猜结构就够了。但是如果是复杂、较大的柔性分子构成的二聚体团簇,或者简单小分子构成的比如六聚体团簇,那要产生的初始构型就得很多,起码一百个。对于原子团簇,也是原子数越多,要产生的初猜数就应该越多。受限于用户的实际计算条件,考虑的初猜结构数目得量力决定。

评分

参与人数 21eV +106 收起 理由
Novice + 5 牛!
qwoop + 5 好物!
ZMINER + 4 好物!
qingqing1223 + 5 好物!
yihanxu + 2 谢谢
Jack + 5 好物!
kunkun + 5 谢谢
一颗赛艇 + 4 赞!
archer + 4 とてもいい!
tjuchan + 5 好物!
zidu113 + 5 牛!
zhou + 5 赞!
captain + 5 好物!
幻七熏 + 5 赞!
zsu007 + 5
aqhuangry + 5
chrischen1128 + 12 好物!
978142355 + 5 GJ!
hcxytpp@163.com + 5 谢谢
zhanfei + 5 赞!

查看全部评分

北京科音自然科学研究中心:http://www.keinsci.com  致力于计算化学的发展和传播,不定期开办各层次量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯
思想家公社QQ群1号:18616395,2号:466017436。超过4000人,用于交流理论、计算化学。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一概不批。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最流行的量子化学波函数分析程序)

此账号为诸Sobereva共用
Money and papers are rubbish, get a real life!

2422

帖子

10

威望

4514

eV
积分
7136

Level 6 (一方通行)

首席卖萌官

发表于 2015-12-15 15:22:29 | 显示全部楼层
沙发是娃娃的!!
She doesn't love me.
Even so,
my heart has been taken away by her.

164

帖子

0

威望

1497

eV
积分
1661

Level 5 (御坂)

Miss X

发表于 2015-12-15 19:36:57 | 显示全部楼层
谢谢老师介绍这么方便使用的软件!!!
真诚相待,淡淡相处,方天久地长^^

295

帖子

0

威望

3190

eV
积分
3485

Level 5 (御坂)

发表于 2015-12-15 22:52:54 | 显示全部楼层
让我想到了packmol !  Sob老师

76

帖子

0

威望

2184

eV
积分
2260

Level 5 (御坂)

发表于 2015-12-16 09:56:45 | 显示全部楼层
不用跑分子动力学了,减少了不少麻烦啊。

30

帖子

0

威望

878

eV
积分
908

科音成员

发表于 2015-12-16 14:21:44 | 显示全部楼层
大赞Sob!

35

帖子

0

威望

997

eV
积分
1032

Level 4 (黑子)

发表于 2015-12-22 13:10:38 | 显示全部楼层
这软件真好,赞!

33

帖子

0

威望

239

eV
积分
272

Level 3 能力者

发表于 2015-12-23 07:19:35 | 显示全部楼层
做团簇的福音啊!大赞!

76

帖子

0

威望

2184

eV
积分
2260

Level 5 (御坂)

发表于 2015-12-28 17:05:52 | 显示全部楼层
对于电中性体系非常不错,但是带电荷的体系就不太好,有什么办法来计算带电体系吗?

13

帖子

0

威望

105

eV
积分
118

Level 2 能力者

发表于 2015-12-30 17:04:40 | 显示全部楼层
      之前都是用Packmol来构建团簇的,那个软件的逻辑是你自己设定一个三维空间的尺寸,然后再设定该空间内的原子数目。所以用Packmol经常还要自己预先估算一下多大空间适合放多少分子,有时候数量填多了会报错,填少了又太疏了,总之是用起来不太方便。这个算法只要输入数量就可以自动生成团簇了,确实简单高效啊。
      不知道用该软件可否构建适用于分子动力学模拟的初始构型?如果可以的话,那么用该软件所构筑的团簇的Box尺寸有一起输出吗,此外除了正方体的Box,可否构建适用于其他更复杂形状box的初始构型,比如长方体或截角八面体?我觉得再多加一些限定条件应该就行了。

1万

帖子

25

威望

2万

eV
积分
41955

管理员

公社社长

 楼主| 发表于 2015-12-31 00:40:31 | 显示全部楼层
幻七熏 发表于 2015-12-28 17:05
对于电中性体系非常不错,但是带电荷的体系就不太好,有什么办法来计算带电体系吗?

如果是正、负电荷离子混杂,需要额外修改程序,增加一些规则,尽量让两类交替出现。
北京科音自然科学研究中心:http://www.keinsci.com  致力于计算化学的发展和传播,不定期开办各层次量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯
思想家公社QQ群1号:18616395,2号:466017436。超过4000人,用于交流理论、计算化学。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一概不批。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最流行的量子化学波函数分析程序)

此账号为诸Sobereva共用
Money and papers are rubbish, get a real life!

1万

帖子

25

威望

2万

eV
积分
41955

管理员

公社社长

 楼主| 发表于 2015-12-31 00:41:50 | 显示全部楼层
tzweir 发表于 2015-12-30 17:04
之前都是用Packmol来构建团簇的,那个软件的逻辑是你自己设定一个三维空间的尺寸,然后再设定该空间 ...


此程序目的不在于建立动力学模拟的初始构型,只是用于产生团簇构型搜索的初始构型,还没打算在建模的功能设置上做扩展。
北京科音自然科学研究中心:http://www.keinsci.com  致力于计算化学的发展和传播,不定期开办各层次量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯
思想家公社QQ群1号:18616395,2号:466017436。超过4000人,用于交流理论、计算化学。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一概不批。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最流行的量子化学波函数分析程序)

此账号为诸Sobereva共用
Money and papers are rubbish, get a real life!

10

帖子

0

威望

19

eV
积分
29

Level 2 能力者

发表于 2016-1-5 14:05:01 | 显示全部楼层
多谢分享,对于计算化学初学者相当有用。

24

帖子

0

威望

1004

eV
积分
1028

Level 4 (黑子)

发表于 2016-1-15 11:44:48 | 显示全部楼层
虽然genmer产生的团簇初始构型可能看起来略微松散,但这完全没关系,一做优化就自动收缩变得致密了。

这样是说结构优化可以起到分子动力学一样的效果吗?我用genmer产生了几个分子团簇,用高斯结构优化,gview看优化的过程觉得结构优化的过程与分子动力学不一样,这样的未经处理的随机初始构型用来结构优化真的可以做到构象搜索吗?
2015毕业季

1万

帖子

25

威望

2万

eV
积分
41955

管理员

公社社长

 楼主| 发表于 2016-1-15 21:03:40 | 显示全部楼层
zhou 发表于 2016-1-15 11:44
这样是说结构优化可以起到分子动力学一样的效果吗?我用genmer产生了几个分子团簇,用高斯结构优化,gvie ...

无论是genmer这样随意堆积,还是分子动力学,目的都是产生随机分布的构型(团簇要说构型而非构象),在势能面上随机采样,这样优化后就能找到最近的极小点。这种做法搜索分子团簇没有任何问题,又方便效果又好。
北京科音自然科学研究中心:http://www.keinsci.com  致力于计算化学的发展和传播,不定期开办各层次量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯
思想家公社QQ群1号:18616395,2号:466017436。超过4000人,用于交流理论、计算化学。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一概不批。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最流行的量子化学波函数分析程序)

此账号为诸Sobereva共用
Money and papers are rubbish, get a real life!
您需要登录后才可以回帖 登录 | 现在注册!

本版积分规则

手机版|北京科音自然科学研究中心|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949-1号 )

GMT+8, 2019-5-26 14:27 , Processed in 0.188699 second(s), 28 queries .

快速回复 返回顶部 返回列表