link管理
链接快照平台
- 输入网页链接，自动生成快照
- 标签化管理网页链接
相关文章推荐
|
留胡子的水龙头 · 黑龙江大学计算机与大数据学院（网络安全学院）· 2 月前 · |
|
|
腼腆的八宝粥 · 人工智能抢走饭碗？美国科技公司开年已裁员万人 ...· 2 年前 · |
|
|
霸气的骆驼 · 征服14座高峰：凡事皆可能 - 知乎· 2 年前 · |
|
|
乐观的火锅 · 何超琼直言自己爱惜许晋亨 ...· 2 年前 · |
|
|
逃课的风衣 · 如何评价《极限挑战》第三季第九期「暗花」? ...· 2 年前 · |
本节介绍一下ASE的另一个骚操作—-扩胞。说到扩胞，这个大家都不陌生，可以有很多种方式来实现，本教程中就介绍了p4vasp的一种操作，今天就看一下ASE的骚操作吧。废话不多说，直接上例子和脚本。
例子：
先展示一下效果吧：
1 |
qli@bigbro:~/Desktop/A19$ ls |
扩胞前：（CONTCAR）
扩胞后：（POSCAR）
脚本：
https://github.com/BigBroSci-LVTHW/LVTHW/tree/master/source/_posts/A19
1 |
#!/usr/bin/env python3 |
前面的扩胞命令：
expand.py CONTCAR 4 4 1
CONTCAR
对应的为
sys.argv[1]
4 4 1
分别为在
x,y,z
三个方向扩胞的倍数。
打赏
如果感觉本文对你的相关研究有帮助，欢迎打赏，支持作者的热心付出。如果你也有自己的骚操作，热烈欢迎无私分享，可以通过QQ（122103465）或者邮件（[email protected]）联系。
ASE 在 DFT 计算中的骚操作令人印象深刻。类似 ASE 的科研工具还有很多，比如 VASPKIT 和 QVASP。熟练使用这些科研工具肯定会助力大家快速完成科研任务，发出更好更高质量的文章。
此次，本文向大家隆重介绍另一款出色的科研工具：Pymatgen。
Pymatgen 是 python materials genomics 的缩写，它是一款基于 python 的、开源的、强大的材料分析软件（ https://pymatgen.org/）。
Pymatgen 包含一系列能够表示元素（Element）、位点（Site）、分子（Molecule）、和结构（Structure）的类（Class）。它具有为很多计算软件提供前处理和后处理的能力。这些计算软件包括VASP，ABINIT，exciting，FEFF，QCHEM，LAMMPS，ADF，AIIDA，ASE，Gaussian，Lobster，Phonopy，Shengbte，Pwscf，和Zeo++等等。它能实现科研狗的众多后处理需求，包括生成相图（Phase diagram）和布拜图（Pourbaix diagrams），分析态密度和能带等等。
Pymatgen 还提供了很多数据库（Materials Project REST API，Crystallography Open Database，and other external data sources）的接口，方便大家从数据库中查询结构和其他数据。
真是科研狗快乐科研之利器呀！
以下是Pymatgen官网提供的后处理的例子：
Top: (left) Phase and (right) Pourbaix diagram from the Materials API. Bottom left: Calculated bandstructure plot using pymatgen’s parsing and plotting utilities. Bottom right: Arrhenius plot using pymatgen’s DiffusionAnalyzer.
此次，本文就介绍一下如何使用 Pymatgen 的 DiffusionAnalyzer 类去计算锂离子固态电解质中锂离子电导率。
计算离子电导率的理论与公式
如何进行AIMD计算
AIMD计算通常非常耗时，所以，为了减少计算成本，我们可以适当放宽计算精度。如果用 VASP 进行计算，具体的，大家可以
后处理的基本公式
一旦AIMD计算完成，大家就可以着手计算离子电导率了。本文首先先介绍以下计算过程中使用的公式，方便有兴趣的同学自己开发脚本。
平均均方位移（average MSD）可以通过以下公式计算：
$\mathbf r_i(t)$ 是第 $i$ 个离子在 $t$ 时刻的位移。
自扩散系数（$D_s$）可以通过以下公式计算：
$d$ 是离子在材料中的扩散维度（一般地，$d=3$），$t$ 是离子扩散的时间。
最后，离子电导率（$\sigma$）可以这样计算：
$n$ 是材料中的离子密度，$e$ 是元电荷，$Z$ 是离子的价态，$k_B$ 是玻尔兹曼常数，$T$ 是温度。
电导率计算的例子
现在我们通过一个 Li_Sn_S 材料的例子来详细了解一下整个计算和处理的过程。该材料的结构显示如下：
本例中采用单胞做计算，INCAR 设置如下：
1 |
[test@ln0%tianhe2 li_sn_s]$ vi INCAR |
AIMD 计算结束之后会得到 XDATCAR 文件。很多时候，由于超算的时间限制，一个完整的AIMD计算需要提交两三次，从而产生两三个 XDATCAR 文件，这时，我们只要把它们按顺序通过 cat 命令合并在一起就行。例如我们有三个 XDATCAR 文件，分别命名成 XDATCAR01，XDATCAR02，和 XDATCAR03。
1 |
[test@ln0%tianhe2 li_sn_s]$ cat XDATCAR01 XDATCAR02 XDATCAR03 > XDATCAR |
新得到的XDATCAR文件，注意删掉重复的与晶格信息相关的行
，一般续算的次数也不多，在使用上面命令的时候，手动把
XDATCAR02
,
XDATCAR03
中的删除即可。
Pymatgen 大显身手
安装pymatgen
1 |
[test@ln0%tianhe2 li_sn_s]$ conda install --channel conda-forge pymatgen |
安装完成后，我们可以试着运行 python，导入 Pyamtgen 模块，如果像下面一样没有出错，就是安装成功了。
1 |
[test@ln0%tianhe2 li_sn_s]$ python |
查看 DiffusionAnalyzer 的类
大家可以通过官方文档（ https://pymatgen.org/pymatgen.analysis.diffusion_analyzer.html）查看接下来要使用的类，熟悉一下代码的用法。
1 |
class DiffusionAnalyzer(MSONable): |
这段代码显示，运行这个类需要一系列的输入信息，包括材料结构（structure），位移（displacements），要研究的离子（specie），温度（temperature）等等。
但是这个类提供了很多方法让大家可以通过读取 XDATCAR 或者 vasprun 文件的方式来实例化，例如
1 |
@classmethod |
好了，废话不多说，直接上代码，开始进行后处理。
代码示例
1 |
''' |
在终端运行该文件
1 |
[test@ln0%tianhe2 li_sn_s]$ python li_conductivity.py |
一段时间后就会得到MSD图像和离子电导率
1 |
[test@ln0%tianhe2 li_sn_s]$ vi result.dat |
可见，该材料在 900K 时的锂离子电导率为 884.05 mS/cm。
例子下载：
链接：
https://pan.baidu.com/s/1WGzOVJBoe6Ym8mvR1uWanA
提取码：jhc5
本文只是浅谈离子电导率的计算，欢迎大家指出计算过程中的不足之处。
打赏
如果感觉本文对你的相关研究有帮助，欢迎打赏，支持作者的热心付出。如果你也有自己的骚操作，热烈欢迎无私分享，可以通过QQ联系（122103465）
对科研狗们来说，ASE不是Automotive Service Excellence的缩写，而是Atomic Simulation Environment。
是一款基于Python语言的工具软件，可以方便地处理DFT计算中结构搭建，准备输入文件准备，读取输出文件，既可以可视化，又可以用来转换文件格式。本文主要介绍ASE在DFT计算中的骚操作。这个”骚”字让我想起来小学的时候，同桌读了屈原的《离骚》，然后转过头来用这个字来笑话我，说我真骚。我被憋的满脸通红却不知如何反驳。骚就骚吧，骚前面加个
离
可以千古流传，那么我相信，它后面加个
操作
也可以对大家日常的计算有所帮助。
首先声明：本教程具有自动过滤功能。每次过滤的时候，我都会提醒一下，以防你被滤纸卡住。但如果你读完教程却不能顺利进行操作，请不要找我，请分解自己让自己滤下去。
滤纸-1 ： ASE成功安装了没？ （ https://wiki.fysik.dtu.dk/ase/install.html）被卡住就不要继续喽！
今天我们讲的骚操作就是用ASE来进行格式转换。这也是本人日常用的最多的一个，最稳定也最方便的一个操作。既然是操作，那么就举个栗子吧：
例子1）：将mol文件转化为xyz。
首先在Chemspider数据库中，下载一个乙烷3D构型的mol文件。
滤纸-2 ：Chemspider 是啥玩意，怎么下载mol文件？ 之前教程有讲过，不会的找百度。被卡住就不要继续喽！
操作过程：
1 |
qli@bigbro:~/Downloads$ ls |
我们成功用ase gui 将mol转化为了xyz文件。
1） 每次运行都会出现那行警示信息，让你用ase convert转化。命令如下：
ase convert -i mol -o xyz 6084.mol C2H6.xyz
这个命令比较长，输入起来烦人，远远不如ase gui 干脆利落。警告信息就用你那高度近视镜自行滤过吧。
2)
-o
中的字母
o
是out的缩写，意思是输出。
3） ase gui 后面跟一个文件，可以直接可视化。
学会了
ase gui
的命令操作，下面我们需要进一步分解自己，以便透过更多的滤纸。
滤纸-3 ：
DFT计算中有很多的文件格式，VASP有POSCAR(CONTCAR)，Gaussian有gjf或者com，Material Studio导出来可以有cif啥的。这些格式之间转换着实让一些新手头疼。以至于经常在群里看到有人询问类似的问题。如果你没有用文本编辑器打开这些格式的文件，并认真分析下它们的数据结构，那么你被滤纸卡住了。请把自己变小一些再继续阅读。
分析一下上面的操作：
ase gui 6084.mol -o C2H6.xyz
6084.mol 是我们读取的对象，C2H6.xyz 是输出的文件。回到大家常用的VASP计算上来，cif文件文件都不陌生，那么我们就转化一下它吧。
例子-2： 将CIF转化为XYZ，POSCAR
首先从Materials Project 下载一个干冰的cif文件。
滤纸-4 ：如果你不知道啥是MP。打开这个网址： https://materialsproject.org，注册，随便下载个cif文件。
操作过程：
1 |
qli@bigbro:~/Downloads$ ls |
OK，大功告成。可以拖到服务器提交任务了。值得注意的是，ASE输出的POSCAR总是把元素行放在第一行的位置，类似于VASP4的POSCAR格式，如果不爽的话，
1） 自行修改为VASP5/6的POSCAR格式。
2） 也可以修改ASE的源代码，
2.1） 找到ASE的安装目录，打开
ase/io
目录下的
vasp5.py
文件；
2.2） 将
vasp.py
中
vasp5 = False
全部改为
vasp5 = True
；
2.3） 保存
vasp.py
即可。
1 |
qli@bigbro:~/anaconda3/lib/python3.8/site-packages/ase/io$ pwd |
如果有其他的需求，可以继续研究ASE或者等后面的骚操作。再继续分解自己：怎么才能把其他格式的文件转化成我需要的格式呢？
上文中的例子已经给你足够的启发了。
妹妹你大胆地往前走
往前走 莫回呀头
通天的大路 九千九百
九千九百九呀
……
本节主要给大家介绍一下中间体计算的具体流程。
计算中间体前的准备：
1）催化剂表面模型的结构以及性质自己多多了解；表面原子排列，谁带正电荷，谁带负电荷，原子间距等等，这些基本的东西自己心里要有数；没事可以整理到表格或者ppt里面，方便以后写文章用；
2） 表面上可能发生的反应路径，各个基元反应，自己列举出来，整理归类方便以后分析。
3）反应路径中各种可能分子的吸附结构先优化完成，通过这些基本分子的吸附能，吸附结构，对体系大体了解一下，一方面参考文献中别人提到的反应路径，另一方面添加自己认为可能的情况。
中间体的优化和前面我们讲到的分子或者原子在表面上的优化是一样的。意思就是计算分子吸附的时候，基本的操作同样适用于中间体的优化。可以分布算，先从粗糙的精度出发，然后慢慢提高精度。注意的一点：
对于计算速度的影响：
1） KPOINTS：先把slab固定住，用（1x1x1）的K点对刚刚搭建的中间体结构预优化一下（基本上跑个NSW = 60就差不多了）得到一个理想的初始结构，然后再增加K点密度，放开slab的表层原子，继续优化；
2）EDIFF，EDIFFG：这两个参数主要影响电子步和离子步的步数，当我们用低密度的K点进行计算的时候，电子步很快收敛，离子步收不收敛我们也不在乎，所以在这一步的时候，不用对这两个参数太较真；关键在于增加K点后，我们选取一个什么样的标准，既得到理想的结构和能量，也避免了服务器做无意义的计算；
3）ISPIN：体系需要考虑磁性的话，粗优化的时候，可以先关掉。后面提高精度后，记得加上就可以了。
4）NCORE等并行参数也不要忘记加上，详细参考：
熟练掌握不同精度下这几个参数的使用，可以极大地提高你的计算效率，节省机时。前面说的这些，基本上试用大部分的体系，但还是要靠大家自己多多测试，摸索一下最适合自己体系的方法。
本节主要取表面上：CH3OH —> CH3O —> CH2O 这两步来简单介绍，下图中圈出来的两步；其他的路径照着葫芦画瓢就行了。
第一步：CH3OH —> CH3O + H
1）反应产物为表面上的甲氧基（CH3O*）和氢原子（H*)， CH3O* 可以继续进行C—H键断裂，或者C—O键断裂，生成CH3*+O 或者CH2O\ + H*。
注意：这里你应该可以想象到反应网络的复杂性，从一个反应物出发可以有多个不同类型的基元反应，产物作为下一步的反应物，也会有不同类型的反应，导致整个反应网络变得越来越复杂，最后交叉合并，汇总到最基本的分解产物上来。这里我们用的是甲醇，当分子C链增长的时候，基元反应的数目就会急剧增加，增加到我们挨个算不能完成的情况，这也是我们为啥要研究BEP关系的原因。另一个原因就是，成千上万个反应，整天瞎算过渡态，有时候你并不会学到太多新鲜事物，只是重复搬砖的工作，对我们以后的成长不利。
2）H原子的吸附，我们就跳过吧，相信大家经过前面的学习，应该都会了。计算CH3O的结构，我们可以从CH3OH的POSCAR出发，进行下面几点思考：
A) 删除POSCAR中对应H原子的坐标，保存为新的POSCAR，即表面上只有CH3O，提交任务进行优化；
B) CH3OH的结构中，O是在一个top的位置上，但O—H 键断裂之后，O的成键能力增强了，还会继续呆在top位置上么？所以，我们要把OCH3挪到bridge, hcp, fcc 位点上，得到几个不同的结构，提交任务进行优化；
C) 前面我们的结构OCH3是斜着在表面上；会不会有直立的结构？ OCH3是通过O与表面结合，我们是不是可以参考O原子在表面的吸附结构；搭建OCH3直立的吸附结构？
D）前面说的这几点，我们都可以通过前面所说的关键点，用一个粗糙的模型扫描一遍。然后大致分析下结果，然后进行下一步的计算，最终得到稳定的OCH3在表面上的稳定结构。
第二步：CH3O —> CH2O + H
这一步，从甲氧基进行了一步C—H反应，得到了甲醛和H。H的话继续跳过；甲醛的话，这里我们应该也需要跳过；因为甲醛是一个分子，在计算中间体前的准备这一步中，我们应该已经完成。安全起见，我们也列举一下几点需要注意的结构搭建工作：
1） 甲醛分子以C，以O还是C=O双键与表面结合？
2） 以不同部分与表面结合的时候，会在哪个位点上？
3）垂直？斜着，还是平行吸附？
总之，尽可能多思考不同的表面结合方式，采用粗糙的模型快速进行筛选，也要实时思考，粗糙模型会不会对这样的快速筛选产生影响？会不会漏掉什么可能的结构？
后面的过程依次类推，直至甲醇分子分解为最基本的化学单元(H，C, O)。
下面介绍另外一种中间体模型，首先回顾一下上面搭建CH$_3$OH—> CH$_3$O + H这一个反应中间体的过程：（这里我们先定义这个中间体模型为：A）
1) 从CH$_3$OH的POSCAR出发，
i) 删除POSCAR中与O相连的H原子的坐标，此时表面上只有CH$_3$O的结构，且CH$_3$O位于top位上；
ii) 稍微修改一下OCH$_3$的结合位点，尝试一下附近的hcp，fcc，或者bridge位点上CH$_3$O的结构，并分别保存到不同文件夹里面；
iii) 提交任务计算，并对比不同OCH3结构的能量，取最低的那个。
iV) 第ii）步中我们暴力直接删掉H原子，基于的假设是一旦O—H键断裂后，产物之间（CH$_3$O , H ）互不影响。
另外一种中间体模型（设为：B）长什么样子呢？
1） OCH$_3$的结构与前面A的模型一样，有top，fcc，hcp，bri各种尝试的结构；
2） 除此之外，CH$_3$OH 中与O相连的H也在OCH$_3$ 附近。
也就是说，在A 模型中，我们直接删掉H原子，而在B模型中，我们只是断开O—H 键，H还留在表面上。也即是B模型保持了与反应物分子相一致的原子种类以及数目。
这两个办法有什么区别呢？哪个更好些?
我们先分析下B)的模型：它比A模型多了一个H，可以认为H 吸附在A模型上；因此A中CH$_3$O的存在会对H的吸附能产生影响，一方面来源CH$_3$O对催化剂电子结构的影响，另一方面来自与空间的影响，类似于前面我们学习过的覆盖度对O原子吸附的影响。因此我们可以计算并对比H在纯净slab和A上面的吸附能。
$\Delta E {ads}$ =$E {ads}^{slab}$ -$E {ads}^{A}$ = [$E {H }$ - $E {Slab}$ - $E {H 2}/2$] - [$E {CH_3O +H }$ - $E_{CH_3O }$ - $E_{H_2}/2$]
= ($E {H*}$ + $E {CH 3O*}$) - ($E {Slab}$ + $E_{CH_3O +H }$)
A和B模型的区别为： 将$CH_{3}O$ 和 H 分离到无穷远时候的体系能量的变化。
再看一下A和B两个模型中反应热的计算：
A)模型：$\Delta$E(A) = $E {CH_3O*}$ + $E {H }$ - $E_{CH_3OH }$ - $E_{slab}$
B)模型：$\Delta$E(B) = $E {CH_3O + H }$ - $E {CH_3OH*}$
反应热的差：
$\Delta$E(A) - $\Delta$E(B) = $E {CH_3O*}$ + $E {H }$ - ($E {Slab}$ + $E {CH_3O +H*}$) = $\Delta E_{ads}$
通过上面的分析，我们得到同样的结果，即两个体系模型的区别主要在于产物之间的相互作用。如果产物之间相互作用弱，可以认为彼此相差甚远，互不影响，也就是A模型，这也是大部分人在计算过程中常常采用的方法。例子有很多，就不一一介绍了，随便取一个： ACS Catal. 2014, 4, 11, 4178–4188。有兴趣的可以去看看。
对于B模型，可以从表面覆盖度对反应的影响这一角度来分析，在键断裂后，如果产物与催化剂结合较强，来不及扩散至互不影响，那么局部的覆盖度就会增加，从而对后面的反应产生影响。有兴趣的可以参考一下这两篇文章：
ACS Catal. 2015, 5, 1, 104–112； ACS Catal. 2014, 4, 6, 1991–2005。
覆盖度对表面反应的影响很复杂，一般人也不愿意碰，但是在计算过渡态相关的操作过程中，B模型却发挥着重要的作用。可以从下面2个角度来分析：1）过渡态插点的技术角度 以及2）获取合理过渡态的角度。这些将会在下一节分析。
上一节我们简单介绍了一下搬砖的原则性问题。这一节，我们主要分析下搬什么砖的问题。继续采用甲醇在金属表面上分解的这个例子。物化书中我们学过表面反应的基本步骤： 1）反应物分子的吸附；2）表面反应； 3） 产物脱附。
甲醇分解反应第一步是先从气相吸附到催化剂金属表面。 对于这一步， 我们自然需要从甲醇的气相分子出发，研究甲醇分子在表面上的吸附结构。其中涉及到一些基本的计算细节：
1） 如何优化气相分子的结构，计算频率呢？
2）如何优化体相结构，切面，优化，计算表面相关性质：比如表面能、功函数等。
这些我们前面也已经讲过了，所以直接跳过。
这里还有一个问题在于，你去用什么表面去模拟催化剂？ 比如说Cu，可以有(111), (110), (100), (211)….等各种不同的指数面，你选哪个来算？为什么选这个来算？这涉及到：我们的研究对象是什么？通过研究这一对象，我们要解释什么？期望获得什么样的结果？对于催化反应来说，我们的研究对象有2个：催化剂和它所催化的反应。
如何获得催化剂的合理模型就是一个至关重要的问题了。一般来说，有下面几个供大家参考的办法：
1）选稳定的表面，也就是表面能低的那个，对于纳米颗粒来说，尺寸如果够大，那么表面能低的表面原子在所有暴露的表面原子中的比例会高，那么可以近似用这个面来模拟催化剂的情况；那么尺寸多大才算足够大呢？
我们可以尝试着从两个极端来进行思考：我们的催化剂如果是金属的话，如果采用一个单金属原子作为催化剂模型，可以想象吸附能肯定很大。当原子慢慢堆积成团簇，团簇再慢慢长到无穷大的时候，吸附能也随之发生变化。我们取无穷大时（最稳定的表面）的吸附能（E0）作为参考，看看金属团簇多大尺寸的时候，吸附能跟E0接近。那么该尺寸以及更大的纳米颗粒就可以用无穷大的催化剂模型来模拟，也就是最稳定的表面。一般来说，纳米颗粒的直径大于2nm（胆儿肥的）的时候，可以用最稳定的表面来模拟催化剂。保守些，直径大约3nm的时候可以放心地去用，不会有什么问题。
下面推荐一本书： Fundamental Concepts in Heterogeneous Catalysis 。书中第21页， 作者用O在Pt(111)表面的吸附能作为参考，测试了不同纳米颗粒尺寸时的吸附能，发现直径在2nm的时候，已经非常接近了。前面说的那个2nm，也是从这里来的。如果你想问我这本书从哪里可以下载，可以关注公众号：BigBroScience，然后回复：ex84 。
2）对于前面的一番说法， 也会有人说，催化剂的活性位不一定在平坦的表面上，边边角角，台阶上可能更容易发生反应。如果你单纯把工作聚集在最稳定的表面上，肯定会遗漏些什么？这也是很多审稿人喜欢折腾我们这些小白鼠的地方。一个金属有很多不同的指数面，单纯考虑最稳定的是一个保险折中的办法。因为不管你的体系搞的再复杂，跟真实情况相比，模型终究是简单的。我们只是尽自己可能去尝试着理解，解释实验的现象，描述下体系的性质，以及进行一些比较简单直接的预测。但是，如果你要深究这些，兴趣也恰好落在在这一块，那么可以对体系的各个相关结构进行一个全面的计算，分析不同表面结构，性质，反应活性有什么不同。在这个全面的计算体系中，
除了各种台阶面，最稳定的表面你还是要算的，这样才能有个对比。原则上说，只要你力气足够，资源充足，可以选任何你想要的表面来计算；但问题是你总会有吐的时候，如果想多算几个面的时候，要留意下面的几点：
注意-1：开放的表面在计算表面能，slab优化或者吸附的时候，层数的影响会很大，这个前期的测试工作要做好；
注意-2：生存法则：
A）灌水第一条：同一个反应， 这个指数面算一下发一篇文章，换个面再发篇文章，再换个面继续发文章….
B）灌水第二条：同一个面： 算一个反应发一篇文章，换个反应算一下再发一篇文章，再换反应；
C）灌水第三条：反应，表面和反应交替，排列组合式发文章。
3）还有一种情况，就是我们手上有实验数据，想模拟一下实验的结果。这个时候，就可以参考下催化剂的表征结果，有针对性地选取催化剂的模型，使其尽可能地去描述催化剂的结构。
注意的有以下几点：
A) 如果做实验的人找你帮忙算个东西，能推掉就推掉，别在上头浪费自己的功夫，顶死也就是个共同一作，国内基本不认可，做无用功；如果刚开始推不掉，就看看他的实验表征结果，一般般的表征以模型搭不出来为由直接拒掉；好的模型是成功的重要因素。前面我们说了，模型再复杂也是简单的。再考虑到实验表征的结果不确定性，你有很大的几率是跟实验上对不上的，到头还是白忙活。很多人会拿工作要冲子刊，杰克斯，俺狗娃为由诱惑你上船，要掂掂自己的分量，别头脑发热用屁股做决定。
B) 对于做实验的人来说，99%会低估计算的复杂性，总以为鼠标点点就OK的事情，其中也不乏一些所谓的大牛们。话虽然难听，却也是大实话：不要瞎几把评论自己不懂的领域。真想找人算一下，要找个靠谱的合作伙伴，找到不靠谱的瞎算，投文章更会给整个工作拉分。
4）最后一个及其重要的就是，认真进行文献调研。看看前人都是怎么选取的，为什么这样选？根据什么依据来选？最终确定自己的计算模型。
关于模型的一些就先到这儿， 本篇的工作主要是重复一遍前人已经算烂的反应，然后在此基础上，分析BEP和TSS线性关系。就直接拿各个金属最稳定的表面来计算了。那么， 什么是BEP线性关系？ 这里的BEP的全称是：Bell–Evans–Polanyi principle. 指的是对于同一类型的N个反应，这些反应的能垒和反应热之间存在的一个线性关系: $E_a$= $\alpha \Delta E$ + $\beta$ 。更详细点的介绍，参考维基百科搜索：Bell–Evans–Polanyi principle。起初这个关系是基于均相反应来说的。在2000年的时候，Neurock将这个关系引入多相催化体系，发现表面反应也存在这一关系。随后，经过了Norskov等人的继续研究，相关的文章就越来越多了。有兴趣的可以通过下面的几篇经典文献开始接触，然后再逐渐深入：J. Catal. 191, 301 (2000)， J. Catal. 197, 229 (2001)， J. Chem. Phys. 114, 8244 (2001)
BEP关系有啥用呢？主要有两点:
关键的一点就是用来节约机时：过渡态不是一个稳定的结构，它的优化需要的计算时间比直接优化反应物种的要多得多。插8个点，也就是8倍的单个优化。再加上计算还经常容易失败。因此我们可以通过BEP关系，跳过TS的优化，直接预测反应的能垒，从而达到节约机时和生命的作用。这个方法简单粗暴，与直接进行反应路径扫描相比，在时间成本上具有压倒性的优势，但精度和准确度上就会稍逊一筹。
另外一点是来预测反应的能垒，这个其实和上面的有些重复，但角度不同。比如一个复杂的催化体系，基元反应数目成千上万，我们挨个算过渡态显然不可行。但是中间体的数目明显少的多，那么我们通过中间体获取反应热之后，通过BEP关系得到反应能垒，也就有了动力学的一些基本数据，从此出发，可以做一些机器学习以及微观反应动力学的模拟工作。从而跨越理论计算的原子分子尺度来到达一个更为宏观的尺寸和时间的尺度，也就是所谓的多尺度模拟。
继续前面的分析，如果中间体的数目也很多，算都算不完，那该怎么办？这个也可以通过一些线性组合，数据分析的方法来进行预测。比如通过Group Additivity的方法，可以将对分子中不同group的能量的求和来获得。具体的内容在后面其他章节再慢慢进行介绍。由于篇幅关系，本节（ex84的上半部分）主要是让大家对自己的研究对象有一个明确的认识，为什么研究它？为什么不研究别的？下半部分介绍中间体计算的一些注意事项。
在线弹性理论中，体系被施加一个无穷小的应变，从而体系应变后的能量可以以小应变为自变量泰勒展开，并忽略二阶导以上的高阶项。通过对不同特定的独立的应变模式求解能量应变方程或者应力应变方程，我们可以获得材料所有的二阶弹性常数（SOECs）。SOECs反映了体系的简谐弹性特征。通过SOECs，一方面我们可以获得体系的基于Vogit-Reuss-Hill平均的多晶力学性能，包括杨氏模量，剪切模量，体积模量，泊松比和维氏硬度等，以及在此基础之上的各向异性系数等等特征；另一方面，我们也可以通过SOECs获得单晶力学性能在空间中的具体分布，并将其用图的模式直观的表示出来。
图中我们展示了某种金属的杨氏模量和剪切模量（每个面的平均值）的空间分布图，感兴趣的读者可以猜猜这是哪种布拉维晶型的金属。
在实际负载下，材料承受的往往是有限形变，这时候材料的非简谐弹性将称为不可忽视的因素。很多实验数据有去测量体系的三阶弹性常数（TOECs），甚至四阶弹性常数（FOECs）。注意到TOECs和FOECs相对比SOECs在实验中更加受到实验精度的影响，具体的测量非常困难。如果我们能在理论上对其进行预测，将变得非常有意义。实际上已经有若干文献进行了相关研究，笔者也抱着好奇的态度尝试进行了TOECs和FOECs的计算。原理上来说相对简单，我们只需要将能量对应变进行泰勒展开，将阶次保留到四次，即可同时获得SOECs，TOECs和FOECs。
将所有的应变η设为同一个数值，我们就可以用一个包含四次，三次和二次项的多项式来拟合计算的数据，获得包含各个弹性常数的多元一次方程组。求解这个这个方程组即可获得所有的弹性常数。
不过实际操作起来并没有这么简单。即使是最简单的立方体系，其TOECs也达到了6个，而FOECs更是达到了11个。对于更低对称性的体系，其计算量更是一个天文数字（开玩笑的，天文倒不至于，很大就是了）。一种避让的方法就是，我们可以用应力对应变进行展开，这样相对需要的独立应变数少一些。笔者使用Au作为测试体系，根据文献对体系施加了11种独立应变，如下图。
感兴趣的童鞋可以直接阅读文末给的参考文献。笔者通过测试发现，首先，体系的应变量必须足够大。笔者一开始犯上了计算SOECs的习惯，使用的很小的应变，计算结果简直南辕北辙。后来将应变加到15%后，获得的TOECs相对理想。细想也是，应变很小时，体系还处于线弹性阶段，这时候的非谐弹性非常微弱，极容易被数值误差干扰。将应变加大到非谐区后，数值就非常稳定了。同时，尽管我们使用了大应变条件，但是实际上高阶弹性常数的具体数值依旧会对数值精度的波动非常敏感。
同文献一样，笔者使用了非常高精度的k点半径和截断能。随着精度的提高，结果趋于较好的收敛。笔者在此展示了笔者自己算的TOECs与实验的比较，两者相对来说还是比较接近的（相对实验的误差与文献差不多）。总的来说，如果有小伙伴需要计算材料的TOECs或者FOECs，笔者建议一定要对应变，k点，截断能三者都要做仔细地收敛测试。
Wang,Hao, and Mo Li. “Ab initio calculations of second-, third-, andfourth-order elastic constants for single crystals.” PhysicalReview B 79.22 (2009): 224102.
本文经ponychen授权发布，版权属于ponychen。
在上文中，我们介绍了用于解释Hall Petch关系的位错阻塞模型和晶界位错源模型。两种理论都是较为早期出现的模型，原理上相对简单，因此文章中仔细介绍了相对应的推导过程。本文将继续介绍剩下的相关模型，由于推理过程相对复杂，读者如果对具体推导过程有兴趣可以阅读对应参考文献。
上文介绍的两种模型都未在公式中直接体现材料变形中塑性应变对Hall Petch系数的影响。因此，位错阻塞模型和晶界位错源模型仅仅只能描述多晶材料在屈服极限之下的情况。我们知道，当多晶材料所受应变超过其屈服极限，材料由于位错的交截和阻碍出现了加工硬化。加工硬化势必会影响Hall Petch系数的数值。在Meakin和Petch（文献1）提出的work haraening（加工硬化）模型中，Hall Petch系数被看成是一个与应变有关的因变量。
由于Hall Petch系数中包含了应变，该模型同时成功预测了材料中的加工硬化行为。如果我们将晶粒尺寸d固定，可以发现该模型描述了一种应力与应变成抛物线形状的加工硬化关系，这与很多实验是相符合的。
目前为止介绍的所有模型都默认认为位错来自于位错源的开动。然而实际上，当晶体被弯曲，或者在晶粒之间的交叠以及孔洞处，由于非均匀变形，晶体为了降低变形带来的巨大应力，会自发形成一些几何必须位错来释放应力。下图展示:
当晶体被弯曲时，晶体中自发在易滑移面上出现的几何必须位错林。在Ashby（文献2）等提出的几何必须位错模型中，作者认为晶界作为晶体中的应力协调区，其上的位错类型主要为几何必须位错主导，并且晶体中绝大部分的几何必须位错都位于晶界上，因此材料的强度由晶界上的几何必须位错密度决定。而平衡位错密度又由晶粒尺寸决定，通过结合加工硬化模型，作者获得了如下的Hall Petch关系式。
位错阻塞模型将材料的强度归功于晶粒内部位错源开动的位错受阻，而目前介绍的其他模型则都将材料的强度归功于晶界处的位错密度受阻。在这两种极端情况的基础上，Meyersm（文献3）等人认为（复合模型），一个多晶材料应该被分为两部分：连续的拥有较高流变强度的晶界部分以及被晶界包围的离散的拥有较低流变强度的块体部分（见下图4，t为晶界厚度）。
当材料承受外力时，晶界部分首先发生变形，形成大量几何必须位错；随着进一步变形，晶界首先发生局部屈服形成加工硬化层，并形核大量位错向晶粒内运动；最后，随着晶粒内位错密度的上升，晶粒部分发生屈服。至此，材料整体完成了整个塑性变形过程。基于这一过程假设，作者得出了如下关系。
可以发现，该关系式并非符合正常的Hall Petch关系，式子中同时出现了一次方和二次方项。作者进一步假设随着晶粒尺寸的增大，其容纳几何必须位错的数量也在上升，因此其厚度t也在上升。基于这一假设，作者得到了下图的强度与晶粒尺寸的关系。
可以发现，复合模型在晶粒尺寸较大的时候可以较好的匹配Hall Petch关系，而当晶粒尺寸小到纳米级的时候，复合模型成功预测了超细晶材料的强度反常降低现象的存在。这是之前所有的理论模型都无法给出的结果。
实际上，在复合模型之后材料学家依旧提出了大量的理论模型，其目的也不仅仅是在于解释Hall Petch关系。由于篇幅所限，本文仅仅介绍了一些具有代表性的模型。实际上，从实验的角度来说，很多材料也并非严格遵循Hall Petch关系中-1/2指数，该指数存在一个很大的范围。对于我们来说，-1/2只是代表材料强度与晶粒尺寸联系的一个优美符号。
本文经ponychen授权发布，版权属于ponychen。
Hall Petch关系是材料学中描述多晶材料最为重要的一个相关关系，它反映了材料强度与晶粒的尺寸的开方成反比。同时材料在极小晶粒尺寸下表现出的反常Hall Petch 关系也一度成为材料学家的研究热点。关于反常Hall Petch关系的研究发现屡现于高水平期刊，本文将不做更多介绍。从Hall Petch关系的第一次报导开始，已经有大量的文章针对各种材料体系（纯金属，金属间化合物，多相合金等等）进行了相关研究。然而，从Hall Petch 关系提出至今，依旧没有出现一个完美简洁的理论可以解释这个现象。研究人员为了解释这个现象已经提出了大量的理论，在接下来的篇幅中，作者将从中提取了几个典型的理论做简要介绍
1951年，Hall（文献1）发表的论文中首次报导了钢中存在的Hall Petch关系:
$\sigma=\sigma 0 + k {HP}d^{-1/2} $
其中$\sigma$表示流变应力，d表示晶粒尺寸，k是Hall Petch 系数，其对不同材料，温度和应变条件都非常敏感。
与此同时，他也在文章中首次提出了试图解释Hall Petch关系的pile up（位错阻塞）模型。
该理论认为材料的强度来自于外力作用下导致的晶粒内部位错源开动，位错源启动后形成的位错林被晶粒之间的界面（晶界）阻碍后引起了材料的应力强化。简单起见，我们不妨假设两个近似六边形的晶粒以晶界AB接触，左边晶粒的中心O点存在一个位错源，在外加应力的作用下位错源沿着某一滑移面释放出一系列位错。随着第一个发射的位错冲撞到晶界AB上，整个位错列队被禁锢在位错源与晶界之间。由于晶界墙的阻力和新产生位错的弹性应变，整个位错列最后会形成一个稳定分布，假设晶粒的平均尺寸为d，则位错列中位错的数目为：
$n=\frac{(1-v)d\tau_e}{2\mu{b}}$
其中$\tau_e$表示位错的有效剪切应力，它可以表示为：
$\tau_e=\tau-\tau_i$
显然，它实际上是在加载应力中扣除了弹性应力部分。试想如果你要产生一根位错，首先晶体要发生塑性变形，因此塑性变形前的弹性部分并不是我们应该要考虑的部分。既然我们已经知道位错列中位错的数目，如果我们粗略地认为，所有的位错作用到晶界上的应力都是等效的，那么晶界上受到的应力$\tau_p$为
$\tau_p = n\tau_e$
Hall认为，材料屈服发生于位错列冲破晶界的时候，也就是说 $\tau_p$大于临界值$\tau_c$。如果我们取等号，再将前面的n代入，我们就可以得到：
$\tau = \tau_i +\surd{\frac{2\tau_c\mu{b}}{1-v}} \times d^{-1/2}$
该公式看起来非常合理简洁，-1/2次方的关系也出来了。但是如果材料的屈服依靠的是位错像裂纹一般冲破晶界，我们应该在实验中大量发现，但实际上，只有少数几篇文献报导过相关现象。更重要的是，依据该理论计算的屈服强度跟实际数据还是有较大出入。
在pile up模型的基础上，Cotrell（文献2）认为，材料的屈服并非来自于位错林冲破晶界，而是由于位错林在晶界的终端的应力场诱发了邻近晶粒中的位错源(如下图的R处)的开动。这样子需要的 临界应力要小于pile-up模型。
晶粒内部出现位错源固然可行，但是另一方面，晶界本身就是天生的位错源，在晶粒内部没有经受剧烈变形的情况下，晶界上的位错形核位置要远远多于晶粒内部。Li （文献3）于1963年提出了boundary source（晶界位错源）模型。为了方便起见，我们不妨假设晶粒是边长为d的立方体。假设晶粒的每一个面上的台阶密度为m，假设每一个台阶都有一定几率发射位错。进一步简化，我们不妨认为每个台阶都会发射位错，则一个面的总位错数目为$md^2$。每一个晶粒与毗邻晶粒公用一个晶面，则一个晶粒可以认为有三个晶界面属于这个晶粒。所以:
我们可以认为一个晶粒有3$md^2$个位错。这些位错在外加应力作用下继续向晶内扩散，形成位错林。位错林的密度为:
$\rho = 3md^2 / d^3 = 3m/d$
另一方面，晶体材料的位错密度和强度之间可以由著名的Taylor公式来描述:
$\sigma = \sigma_0 + M\alpha\mu{b}\surd\rho$
代入前者即得:
$\sigma = \sigma_0 + M\alpha\mu{b}3^{1/2}m^{1/2}d^{-1/2}$
我们再一次获得了Hall Petch关系。在原文中也有说明这个方式计算的理论屈服强度与pile up模型基本一致。目前介绍的这两个模型，前者认为位错完全起源于晶粒内部，后者则认为位错完全起源于晶界，两者分别是实际情况下的两个极端。但是殊途同归，两个极端最后得到了相同的结论和近似的屈服强度。不过，这两个模型的局限性也是明显的，实验中Hall Petch曲线的斜率对晶粒结构和化学成分非常敏感，但是上述模型中并未考虑这些实际因素。
值得注意的是，在2002年的一篇文章4中，该文作者成功在MD模拟中发现了晶界处位错的形核和发射过程。
（未完待续）
1 Hall, E. O. The deformation and ageing of mild steel: III discussion of results. Proceedings of the Physical Society. Section B 64.9 (1951): 747.
2 A.H. Cottrell,: The Mechanical Properties of Matter, Wiley, New York (1964)
3 Li, James CM. Petch relation and grain boundary sources. Transactions of the Metallurgical Society of AIME 227.1 (1963): 239.
4 Van Swygenhoven, H., P. M. Derlet, and A. Hasnaoui. Atomic mechanism for dislocation emission from nanosized grain boundaries. Physical Review B 66.2 (2002): 024101.
本文经ponychen授权发布，版权属于ponychen。
基于DFT+ U 研究应力调控CeO2(111)催化分解水制氢的活性
氢能清洁无污染，是最具潜力的能量载体。使用太阳能、风能和水能转换得到的清洁电能电解水制氢是人来未来最具发展潜力的、可持续的、绿色环保的制氢技术。传统低温碱性电解水制氢技术的效率受限于液体电解质较低的传质过程和氧电极较慢的反应动力。相比于液体传质技术，新兴发展的高温固体燃料电池 (SOEC)技术以固体作为传质媒介，充分利用工业废热提供反应所需的热能，具有极高的传质效率和水电解效率。然而人们对SOEC制氢反应机理尚缺乏深刻认识。计算模拟技术是揭示能源材料领域中多相催化反应过程及机理的重要手段，可以辅助了解SOEC电解水的关键反应步骤，以开发、低廉且稳定的催化剂用于增强氢电极的稳定性以及提高制氢反应在燃料电池中的效率。兼具优异的电子和离子传质能力的氧化铈（CeO2）是SOEC水分解重要的催化剂和稳定剂。CeO2优越的催化活性归因于Ce3+与Ce4+的氧化还原循环。基于密度泛函理论，该工作使用Hubbard- U （ U eff=4.5 eV）修正的方法来描述Ce 4 f 轨道电子的强关联作用，研究表面应力对CeO2(111)催化制氢的影响。
水在CeO2表面分解制氢包括三个基本反应步骤：（1）CeO2(111)表面生成一个氧空位；（2）H2O吸附于氧空位附近最终形成表面OH-；（3）表面OH-分解形成H2并留下一个晶格氧，如Figure 1所示水分解过程中形成的反应中间体以及最终生成H2的反应过程。该工作主要研究了三种可能的反应路径，即不同OH- 覆盖率下形成H2的过程如Figure 1所示。通过研究应力对三种不同反应路径中反应中间体形成的影响，以及应力对各反应步骤反应能垒的调控，来探究应力如何增强CeO2催化制氢的性能。
首先，研究发现拉伸应变有效增强了氧空位和OH-等中间体的稳定性，随之增加了各反应步骤的反应能垒。与之相反的是压缩应变可以通过减弱OH-与表面的结合能，从而降低生成H2的能垒。通过计算不同反应温度下中间体的自由能可以获得沿着2H，4H和5H的反应路径能量图。研究发现不管处于何种应力条件下生成氧空位和氧空位的扩散，以及OH-的形成不需要很高的能量势垒，而OH-分解形成H2是水分解制氢反应路径中能垒最高的一个步骤。然后，通过对比不同反应温度下，研究发现在温度低于1000 K时，沿着5H路径可以最快速地水分解生成H2 如图2所示，且最佳反应路径不受应力的影响。然而，当温度升高到1200 K时，应力有效地调控了制氢的反应路径，即拉伸应力下最优路径沿着5H，而在压缩应力下最优反应路径沿着2H或者4H。此外，研究还发现反应速度随着压缩应力的增大而增快，且当压缩应力大于3%时，反应速度获得显著提升。因此，计算模拟理论研究显示压缩应力（>3%）有利于增强CeO2(111)催化水分解制氢的反应活性。
图1. 2H，4H和5H覆盖的CeO2(111)表面生成H2的反应路径示意图。黄色, 红色，蓝色，白色，黑色和灰色原子球分别代表Ce，表面O，次表面O，H，次表面氧空位和表面氧空位。
图2不同反应温度与应力条件下最优反应路径生成H2的反应速度对比图。
本文第一作者，武甜甜为丹麦科技大学能源转化与储存系博士后研究员。
Google Scholar：
https://scholar.google.com/citations?user=rm84jfsAAAAJ&hl=en
指导老师为：Tejs Vegge教授。
课题组主页： https://www.dtu.dk/english/service/phonebook/person?id=5334&tab=1。
本文经作者授权发布，版权属于第一作者。
通讯作者：Philippe Sautet
胺是工业上生产药物、农用化学品、生物活性剂、聚合物和染料的重要中间体化合物，其中伯胺是衍生化反应的关键试剂。传统合成胺的方法包括氨与烷基卤的直接反应，或通过还原硝基络合物和腈，但是这些工艺的缺陷在于产生有害的副产物、废弃物并消耗大量氢气。而利用氨气通过借氢反应机理直接胺化脂肪类醇用于合成胺被认为是一种有潜力的绿色合成路线并且已经有数千吨规模的工业化应用。然而，合成胺工业仍然面临诸多挑战，尤其是选择性合成伯胺。在此背景下，探索催化剂的本征性质与合成胺活性和选择性的依赖关系变得尤为重要。近年来，理论模拟的飞速发展使得揭示催化过程中的构效关系成为可能。本文结合第一性原理微观动力学模拟和实验，揭示了决定醇的胺化活性和选择性的关键因素并构建了其与碳和氧原子吸附能的依赖关系，实现了该反应金属合金催化剂的理性高效筛选，为实验上醇胺化催化剂的设计和升级提供了重要理论参考。
内容和讨论
A. 醇的胺化反应机理: 首先甲醇的胺化被作为模型反应进行研究，该反应机理主要包括三部分(如图1)：醇的脱氢生成醛；醛与NH3迅速反应生成亚胺；亚胺的加氢生成胺。基于上述反应机理，DFT计算在九种过渡金属的密堆积表面系统展开。由于实际反应气氛中具有大量的NH3，因此催化剂的表面会覆盖一定浓度的NH3，早期的研究发现表面覆盖1/9 ML的NH3有助于准确描述胺化反应催化剂的活性顺序。其中，第二步醛与NH3反应生成亚胺被认为是气相里的快速平衡反应而且是理论模拟上的难点，该工作通过详细的反应器动力学模拟证明采用Eley–Rideal 反应机理能够准确描述该关键步骤。最终为微观动力学模拟提供了可靠的反应机理。
图1 : 甲醇与氨气的胺化反应机理
B. 线性降维: 基于Norskov课题组提出的d带中心理论和scaling关系，反应过程中表面物种的吸附能并非相互独立的而是受制于scaling关系，并且基元反应的能垒也可以通过BEP关系与反应热关联，因此通过详细的描述符的筛选以及线性拟合可以实现降维和最优描述符的确认。由于该反应涉及含C、N、O的三类中间体，scaling关系可以最终将该反应中涉及的中间体的能量与C、N和O原子的吸附能进行关联。然而，包含三种不同描述符导致无法利用典型的二维火山曲线进行活性描述。深入研究发现：N原子的吸附能同样依赖于C和O原子的吸附能并且具有完美的线性关系（如图2）。而该依赖关系的物理化学根源一方面为N(3.0)原子的电负性恰好介于C(2.5)和O(3.5)原子之间，另一方面N原子的最高占据轨道能量(HOMO)可以线性分解为C和O原子的相应能量。基于上述分析，甲醇胺化反应的能量维度被降低为二(C和O原子的吸附能)。
图2 : N原子的吸附能与C和O原子吸附能的线性关系
C. 醇的胺化活性和选择性与简单描述符的依赖关系: 基于图1的反应机理和图2的线性降维，通过构建微观动力学模型，最终建立了醇胺化的反应速率(TOF)与C和O原子吸附能的依赖关系。如图3a所示，由于对C和O原子吸附能力的区别，九种金属分布于火山图的不同位置而呈现出不同的反应活性。然而，针对醇的胺化反应，人们更加关心且关键的问题是选择性。因此，该工作进一步系统考察了导致伯胺选择性降低的不同可能性，最终发现CH3NH和CH2NH物种的耦合是导致选择性丢失的关键因素之一，并巧妙的将其与简单描述符进行了关联，最终构建了醇胺化生成伯胺的选择性与C和O原子吸附能的依赖关系如图3b。随着CH3NH/CH2NH耦合能垒( E a)的增加，副反应变得困难，进而催化剂的选择性增加。至此，通过计算催化剂表面C和O原子的吸附能，结合图3中的火山图便可定性预测相应催化剂的醇胺化反应的活性和选择性，从而大大提高催化剂的筛选效率。
图3: 醇的胺化活性(a)和选择性(b)与C和O原子吸附能的依赖关系
D. 实验室论证预测: 基于上述活性和选择性图，通过对元素周期表中不同过渡金属进行两两排列组合，最终计算了三百余种双金属合金的C和O原子吸附能，并预测出系列有潜力的醇选择性胺化催化剂。基于上述理论预测，索尔维公司研发中心对不同合金进行了详细的实验室测试以论证图3的可靠性。具体的：本文合成了系列Co基合金催化剂，一方面由于Co金属的相对廉价，另一方面Co具有较高的活性和伯胺选择性。根据图3a的理论预测，提高Co的活性需要设法降低其O原子的吸附能，由于Pt和Pd具有弱的O吸附能力，因此将Co与Pt和Pd形成合金将可以实现反应活性提高的，然而基于图3b的信息，这将同时伴随着伯胺选择性的降低。事实上，图4和图5的实验结果也证实了CoPt和CoPd合金相比Co金属的高活性和低伯胺选择性。此外，图3表明：将Co与Ag和Ru组合将会维持甚至提高活性且提高选择性。该预测结果同样与图4和图5中的实验结果一致，并且成功设计出活性选择性都大大提高的CoAg合金催化剂。本文预测的其他有潜力的合金催化剂也将进一步推动实验上对醇胺化反应研究。
图4: 实验上不同催化剂醇的胺化生成伯胺的TOF
图5: 实验上不同催化剂醇的胺化生成伯胺的选择性
总结和展望
通过精确的第一性原理计算和微观动力学方法，准确的模拟了醇的胺化这一复杂反应的详细机理，揭示了决定金属合金催化剂醇的胺化活性和选择性的关键因素，并与简单易得的C和O原子吸附能进行关联，系统的实验测试论证了理论预测的可靠性。最终大大提高了醇胺化反应金属合金催化剂的筛选效率，为工业上醇胺化催化剂的设计和升级提供了参考。此外，该工作采用的改进的动力学模拟方法能够引入反应器环境对复杂反应机理的影响，一方面实现对实验条件更加精准的描述，另一方面能够准确定性预测金属合金催化剂的活性和选择性的变化趋势，为后续高通量催化剂筛选提供了基础。
文章第一作者王涛博士目前为美国斯坦福大学SUNCAT中心博士后研究员，通讯作者Philippe Sautet为美国加州大学洛杉矶分校教授、ACS Catalysis副主编、法国科学院院士、美国能源局IMASC能源前沿研究中心副主任，长期致力于理论催化方向的基础研究。此外，法国里昂高等师范学院、法国里尔大学、比利时根特大学以及跨国公司索尔维研发中心E2P2实验室也为该工作提供了支持。本文经大师兄涛本人授权发布。
文章链接： https://www.nature.com/articles/s41929-019-0327-2
关注公众号：大师兄科研网，回复： 胺化 ，可以获取下载链接。
上海交通大学密西根学院朱虹老师课题组诚聘博士后 1~2 名，负责或协助团队进行金属合金体系的腐蚀性能的理论预测和腐蚀机理的研究，将与腐蚀实验团队以及世界各地的材料基因组计算研究团队开展紧密合作 (UC Berkeley, University of Michigan, UC San Diego, Georgia Tech, University of Maryland, etc)，并在现有的高通量计算平台上开发针对腐蚀相关基础参数的高通量计算工作流。诚挚邀请具有理论计算、电化学背景同学的加入。
朱虹 ( [email protected] ) 助理教授，上海交通大学材料基因组联合研究中心成员； 材料科学与工程学院，上海交通大学密西根联合学院双聘。
Google Scholar: https://scholar.google.com/citations?user=x1BGIfEAAAAJ
课题组主页: http://umji.sjtu.edu.cn/~hzhu 。
发邮件至 [email protected] ，并附上简历以及可以到站的时间。
Python现在越来越流行，连微信朋友圈都开始推广一些python培训课程广告了：俩美女一起上班，一个瞬间完成工作，另一个一脸懵逼地问道你咋这么快就干完活啦？另一个回答道：因为我学习了Python啊。是的，Python确实有这些神奇的功效，但在网络这么发达的时代，想学习python的途径有很多，参加什么割韭菜的培训班啦，找个计算机学院的男女朋友啦，向师兄师姐学习啦，自己学习琢磨啦等等。
学习Python的一个好处就是，很多懂python的都会github上分享一些自己的脚本，但这些脚本很多时候不是为你的课题量身定制的。稍微懂点可以加以修改就变成自己科研的一个利器了。更重要的好处就是广告里面所凸显的，由于计算中很多操作都是可以重复的，写个小脚本可以快速完成自己的工作任务，提高效率，省出时间该干嘛干嘛去。此外，很多复杂的任务，手动计算量很大，也得需要写点脚本来帮忙。那种瞬间完成任务的感觉，不会点程序语言的人是不会体验到的。而且，现在小学生，初中生都开始学习python了，我们也没有理由不跟上时代的潮流。
看到这里，应该能明白本文的目的。不要你交99元，也不给你优惠99元。大师兄推荐一本武林秘籍，专门针对科研人员的一本科学计算相关的Python书。里面大量例子简单易懂，解释详细，带你缓缓走入Python的程序世界，平日闲的没事加以练习，犹如武林高手给你灌输真气一般，达到计算的另一个境界，轻而易举超越那些培训班被割过的韭菜们。
获取方式，扫描右侧或者下方二维码关注公众号：大师兄科研网，回复： 武林秘籍 就可以获取下载链接了。为了避免别人学会了python跟你抢男女朋友，课题组争宠，大家下载后私下练习就OK啦，有对象的，一起练习讨论效果更佳。
（University of Massachusetts Lowell）
课题组简介
车芳琳老师课题组刚刚成立，主要研究催化以及材料科学方面的多尺度理论模拟，在多尺度范围内（从原子到反应器尺度）研究燃料电池，电催化及新型能源、光催化相关的催化剂性质，反应机理，动力学等。研究的相关课题有：轻质烷烃的活化，CO2电催化，多相催化的线性关系以及微观动力学模拟，以及多相催化反应器中的流体动力学模拟。课题组主页： http://sites.uml.edu/fanglin_che/ (或点击左下角阅读原文）
车老师本科毕业于大连理工大学（2008-2012），博士毕业于华盛顿州立大学（2012-2016），师从J.-S. McEwen教授。分别在加拿大多伦多E. H. Sargent课题组（ 2017-2018）和美国特拉华大学Dionisios G. Vlachos课题组（2018-2019）从事博士后研究。至今以一作（含共一）在Nat. Chem. (1), Nat. communs (2), Angew. Chem. Int., Ed.(1)， ACS Catal.(3), J. Catal.(1), ACS Energy Lett.(1), Adv. Mater(1), Appl. Catal. B (2)等高水平期刊发表文章若干篇。
Googl Scholar：
https://scholar.google.com/citations?hl=zh-CN&user=fiRHcIQAAAAJ
目前长期招聘硕士，博士，以及博士后。有兴趣的可以将个人简历发送至： [email protected]
我们计算的时候，经常会遇到扩胞的需求。比如，我们要优化（4x4）的Ag(111)的一个slab。直接切一个slab拿来优化，可能会比较耗时。另一个办法就是我们先优化一个（1x1)的Ag(111)的slab，然后在将优化完的结构扩成（4x4)的，最后再优化。这样能有效地减少工作量。
扩胞的话，有很多工具可以选择，MS, P4vasp,等等可视化的软件，用鼠标点点就OK了。也可以使用一些现成的脚本，小程序。本节就介绍一个通过ASE进行扩胞的小脚本，也是本人偶然在一个网站发现的。有兴趣的可以自己看下： https://www.nsc.liu.se/~pla/blog/2013/02/26/vaspsupercells/
废话不多说，直接上例子，例子完了是脚本的具体内容。
如果你有自己的脚本或者推荐的程序，也欢迎发送到本人邮箱：[email protected]。后面我们会逐渐扩展本节的内容。
1 |
qli@bigbro:~/Desktop/test_expand$ ls |
脚本内容如下：
1 |
#!/usr/bin/env python3 |
使用前提是你已经安装好了ASE。
再次啰嗦一下，欢迎大家分享自己扩胞的小脚本，推荐的小程序，一起完善本节的内容。请发送到本人邮箱：[email protected]。
汇集本人日常学习中的一些数据库，个人感觉不错的网站。也会不断地更新完善。如果你有好的数据库，链接，推荐的书籍，也欢迎发邮件至[email protected]留言帮助我们更新。
结构查找
The Materials Project
https://www.materialsproject.org/
使用方法：
1） 点击右上方：Login
2） 点击： Sign in with your email address
3） 输入自己的邮箱，然后点击：send login link
4） 打开邮箱，点击链接，即可查询数据库。
ioChem-BD
本人老板和其他几个老板一起做的数据库，涵盖VASP，高斯等主流程序，可以上传结果。数据库依托于巴塞罗那超算中心。本人博士，博后的所有工作都在这个数据库中。
ChemSpider
http://www.chemspider.com/ RSC旗下的一个数据库，非常方便查找下载结构。
Pubchem
https://pubchem.ncbi.nlm.nih.gov/
ICSD 数据库
http://www2.fiz-karlsruhe.de/icsd_home.html , 需要账号密码
CCDC: The Cambridge Crystallographic Data Centre
cccbdb 数据库
团簇数据库：
http://www-wales.ch.cam.ac.uk/~wales/CCD/CoCCD/cobalt.html
可能有些老了。
https://ww2.mathworks.cn/matlabcentral/fileexchange/33449-cluster-generator
Wulffman
https://www.ctcms.nist.gov/wulffman/
XPS
Cn Fullerenes
http://www.nanotube.msu.edu/fullerene/fullerene-isomers.html
热力学数据
NIST 数据库： https://physics.nist.gov/cuu/Constants/index.html
CRC Handbook of physical chemistry ： http://hbcponline.com/faces/contents/ContentsSearch.xhtml 查找热力学参数，晶格常数，entropy等等。
http://www2.ucdsb.on.ca/tiss/stretton/database/inorganic_thermo.htm 无机化合物的物理和热力学相关数据
http://www.knowledgedoor.com/2/elements_handbook/cohesive_energy.html 结合能
http://www.surface-tension.de/ 表面张力
http://www.wiredchemist.com/chemistry/data/thermodynamic-data
物理常数
https://physics.nist.gov/cuu/Constants/index.html
http://wild.life.nctu.edu.tw/class/common/energy-unit-conv-table.html
http://halas.rice.edu/conversions
书籍文献下载：
资源丰富，但不要贪多，更不要下一堆放到自己电脑里占地方。找几本好书，认真静下心来好好学习。
http://www.expaper.cn/ 科研速递论坛， 求助文献，书籍，成功率极高。
sci-hub 系列，现在网上各种域名，插件，桌面版的都有，我就不瞎凑热闹了。
5000人QQ群: 157099073 汇聚大家的力量，解决每一篇小文献的困难。
书籍和文献下载的原则，宁缺毋滥，认真阅读每一篇下载的文献，对自己的时间负责。别图多。
VASP 相关
http://www.vasp.at/ , VASP最权威，最专业的学习资料都在这里面了。英文刚开始看起来有点不适应，坚持下去，慢慢就好了。不建议看除官网以外的那些乱七八糟的说明（包括我自己写的教程）
VASP 教程视频
http://www.nersc.gov/users/training/events/3-day-vasp-workshop/ , A 3-day VASP Workshop at NERSC: VASP 开发者（长发大胡子那哥们）做的workshop, youtube 有视频。链接： https://www.youtube.com/playlist?list=PL20S5EeApOSumWZkzsaYxAvozjvFZ3ks4
VASP官方论坛
https://cms.mpi.univie.ac.at/vasp-forum/
VASP官方论坛，需要注册才可以提问，但一般来说你的问题已经有人提问过了，直接复制到google里面，基本都会显示该论坛的链接。可以不保存。
Henkelman 课题组
http://theory.cm.utexas.edu/henkelman/ , 大量VASP相关的脚本，VTST计算过渡态，bader电荷分析，CHGCAR处理脚本等，向Henkelman致以崇高的敬意。
ASE ： Atomic Simulation Environment： https://wiki.fysik.dtu.dk/ase/， 适合催化计算相关
Pymatgen ： http://pymatgen.org/ , 适合材料计算的
K-path ： https://www.materialscloud.org/work/tools/seekpath
能带Kpath选择的网站，上传自己的结构，便可以得到建议的路径。
Chemml : https://hachmannlab.github.io/chemml/
p4vasp ： http://www.p4vasp.at/
RDkit ： https://www.rdkit.org/
VASPkit ： https://sourceforge.net/projects/vaspkit/files/ 来自中国的良心软件, VASP计算前后处理，功能强大，使用简单，老少皆宜。
VESTA ： http://jp-minerals.org/vesta/en/ 来自日本的良心软件
Xcrysden ： http://www.xcrysden.org/， 用于批量做图，选择K-path做能带图，等。
Openbabel ： http://openbabel.org/wiki/Main_Page 强大的格式转换工具，它说第二，没人敢说第一。
RDkit ：化学信息学重要的开源包。
程序软件相关：
Bash ： http://tldp.org/HOWTO/Bash-Prog-Intro-HOWTO.html
https://learnpythonthehardway.org/
https://www.liaoxuefeng.com , 廖雪峰的网站，英文不想读的话，可以到这里学习。
https://matplotlib.org/ 用python画图
https://www.scipy.org/ python 科学运算相关的库
https://stackoverflow.com/ 一般谷歌搜索遇到出现的程序问题，都可以在这里找到答案。
https://stackoverflow.com/ 一般谷歌搜索遇到出现的程序问题，都可以在这里找到答案。
大牛们的Git-Hub 网站
很多懂程序, 喜欢分享的大牛们都有自己的github网址，下面会逐渐完善。
1） https://github.com/Ionizing （www.bigbrosci.com这个网站就是这牛人帮忙搭建的）
2） https://github.com/wangvei （VASPkit的作者，王伟老师）
3） https://github.com/QijingZheng （很屌的一个大牛，虽然没见过）
4） https://github.com/obaica （只要Follow了他，基本就可以找到其他大牛的github了）
5） https://github.com/tamaswells （脚本小王子，the king of sharing）
6） https://github.com/jensengroup （总有一些人，会默默地支持穷人们的科研）
7） https://github.com/LePingKYXK
8） http://home.ustc.edu.cn/~lipai/scripts/vasp.html
理论计算，推荐书籍：
1） Density Functional Theory: A Practical Introduction
http://onlinelibrary.wiley.com/book/10.1002/9780470447710
有权限的直接下载；没有权限的不要在网上瞎求，很多公式都不全，这个要注意。
2）Atkin’s Physical Chemisty （10th）
网上到处都是，最新的是第10版，如果下载不了，看第9版即可。
3）of Modern Catalysis and Kinetics
做多相催化的强烈建议！ http://onlinelibrary.wiley.com/book/10.1002/3527602658
以上书籍，如不提供下载链接，基本都可以在 http://gen.lib.rus.ec/ 这个网站获得。
溶剂相关的：
https://www.organicdivision.org/wp-content/uploads/2016/12/organic_solvents.html
http://murov.info/orgsolvents.htm
求职相关
Psik Network： http://psi-k.net/
注册后，一旦有人在网络里面发送海外职位信息，你便会收到邮件提醒。更为重要的是，网站里面有很多珍贵的学习资料。
Linkedin 领英 https://www.linkedin.com/，
有些老板也会在这上面发布招聘信息，本人当年的博士位置，就是在这里面发现的。
Research gate： https://www.researchgate.net/ ,
可以用来下文献，关注同行最新进展。
（七）理论计算，推荐书籍：
1 DFT: A Practical Introduction
http://onlinelibrary.wiley.com/book/10.1002/9780470447710
1）有权限的直接下载
2）没有权限的不要在网上瞎求，很多公式都不全，这个要注意。
2 Atkin’s Physical Chemisty （10th）
网上到处都是，最新的是第10版，如果下载不了，看第9版即可。
3 Concepts of Modern Catalysis and Kinetics 多相催化的强烈建议！
http://onlinelibrary.wiley.com/book/10.1002/3527602658
以上书籍，如不提供下载链接，基本都可以在 http://gen.lib.rus.ec/ 这个网站获得。
1 颜色数据库：
https://www.rapidtables.com/web/color/RGB_Color.html
2 元素周期表：
http://www.rsc.org/periodic-table
3 Latex相关：
本人的文章，除合作的外， 基本都是用latex写的，但是推荐网址的话，确实想不起有什么来。找一个latex入门pdf，练上几天，遇到问题google搜索，慢慢就好了。如果非要推荐一个网址的话，那么就选这个： https://tex.stackexchange.com/ 和 www.google.com
也欢迎大家关注微信公众号
前面我们介绍完了算表面吸附，以及过渡态的一些基本的操作，和注意事项。当我们面对一个新的课题时，就需要运用所学到的这些技能来完成所必需的计算，来验证我们的想法，思路等。后面几节，主要参考本人2014年发表的一篇关于甲醇分解反应计算的文章来介绍一下怎么运用所学习到的这些本领。（ https://pubs.acs.org/doi/abs/10.1021/cs501698w）
搬砖这个词用来形容我们完成一个课题的计算过程，简直不能再形象不过了。
首先，搬砖这一过程， 我们要明确几个要素：
1）老板：让你搬砖的人
2）搬砖工人：你自己
3）任务：搬砖
4）目的：盖房子
下面，你要明确一点：不会计算的搬砖工不是好搬砖工。
这里说的计算，不是我们的计算化学中的计算。而是把我们的工作如何一步一步分解，结合自己的实际情况，算算需要多少的时间完活。毕竟房子总有盖完的时候，搬砖工不能一直守着。所以，当老板告诉你要盖一个大房子的时候，需要多少砖这一点你要清楚。然而，很多人都是第一次搬砖，或者之前没搬过表面计算的砖，对上面说的这一点并没有一个很清晰的概念。所以，摸清自己的现实情况，这也就是本节的重点啦。
人有高低胖瘦之分，搬砖工也不例外，有的劲大，有的劲小。劲大搬的砖多，劲小搬的少。最后到底能搬多少砖，一方面取决于你的身体条件（自身因素，搬砖前的学习经历），另一方面取决于你的工资伙食情况（老板因素）。身体倍儿棒，老板钱儿多，这样的情况通常不会看本教程，有啥问题自己（组里）就能解决。身体倍儿棒，老板没钱，或者身体倍儿不棒，老板有钱的。就需要认真动动脑子，分析下当前的形势。没劲儿又没钱的，就应该更加注意了，需要更多地在脑力上下功夫。
比如，这篇文章中，老板让我们算一遍甲醇在四个金属上的分解反应（倒过来就是合成反应）。
第一步：分析下需要搬多少砖块
2）我们取的四个金属：Cu，Ru，Pt，Pd稳定的表面。所以前面的计算要乘以4。
3）前面说的是最理想的情况，而实际情况则是：
A) 对每一个表面吸附的物种来说，我们要尝试不同的位点，来获取一个稳定的构型；
B）对基元反应来说，这些过渡态，100%不可能都100%一次性找到，还要考虑不收敛，不同的可能性等；
这些都会使得计算量增加。所以，第一步，大体上有多少需要计算的东西，应该有个框架。
第二步：确定一个合理的搬砖方案。
A) 善于利用已经发表的文章的数据，比如，有些结构可以在支持信息里面找到，还有一些数据库里可以下载，也可以问作者要（可能性比较低）， 可以理解为找朋友一起帮忙搬砖；
B) 选择合适的slab模型：3x3还是4x4，slab取四层，还是取5层？这些是课题开始之前一些基本的测试工作，可以参考文献中别人的做法，也可以根据测试的结果自己合理选择。可以理解为：搬多大的砖块。
C) 选择合适的参数: 计算参数不对，很可能导致计算的结构或者能量有问题。 但这些都要具体分析，有些能量有问题，但结构还算OK的可以调整下参数，作为一个理想的初始结构继续用。可以理解为：半路翻车，捡起来那些没摔碎的转，继续搬。
D）善于使用前面我们介绍的快速获取理想初始结构的方法。可以理解为先用车把砖块搬到离工地最近的地方，省去往返来回跑的劲。
上面说的这些，尽可能在课题完全开展前做到位，因为它们不会花费很大的劲去做，但会节省后面很大的劲。而且，伙食好坏（计算资源给不给力），测试的过程一目了然。
第三步：结合自己的体力，伙食，认真搬砖。
这一步就简单啦，体力好同时进行，左手生擒中间体，右手活拿过渡态。体力不好，俩手搬一块砖，累了饿了（没资源）就一边凉快去。
小结
不过对于新手来说，课题刚刚开始就真正掌握上面的这几点，难度有些大，所以建议先做些前提的准备工作：
1） 熟悉自己的服务器，运算能力；
2）多查找文献，整理需要计算的框架，幸运的话，可以从支持信息直接找到结构；
3）多用小体系做测试，测试完了要对结果多思考总结。不要上来就狂交任务，最后把服务器累个半死，还不出什么好的结果。
伴随着Python2的落幕，Python3的兴起，很多基于2的软件也面临着着要么被遗忘，要么主动拥抱Python3的选择。其中就有我们最喜欢的p4vasp。 https://github.com/orest-d/p4vasp 或者 www.p4vasp.at 。本人不认识作者，也不知道啥时候出新的版本。使用p4vasp需要有python2，但其他软件又需要python3。为此我们不得不做出一些选择， 最残忍的莫过于弃用p4vasp啦。本节介绍两个办法：
1） 修改p4vasp运行默认的python版本
2）通过update-alternatives切换系统默认的python版本。
办法-1
第一步
： 修改
p4v
脚本：将最后一行中的
python
改成
python2
1 |
qli@lhz:~$ sudo gedit /usr/bin/p4v |
第二步
修改
p4v.py
， 将第一行中的
python
改成
python2
运行实例：
1 |
qli@lhz:~$ python --version |
办法-2
百度或者google里面搜一下：update-alternatives python 就会得到一系列的操作教程。本文不再瞎扯。
查看自己系统中python2和python3的版本，终端里面输入：python然后摁tab键，显示如下，本人安装了python2.7 和python3.7。
1 |
qli@lhz:~$ |
设置update-alternatives：在终端里面，依次输入下面两个命令行：
1 |
qli@lhz:~$ sudo update-alternatives --install /usr/bin/python python /usr/bin/python2.7 1 |
1 |
qli@lhz:~$ sudo update-alternatives --install /usr/bin/python python /usr/bin/python3.7 2 |
这样就设置OK啦。更加具体的，网上已经泛滥了
切换python版本，使用下面的命令：
1 |
qli@lhz:~$ sudo update-alternatives --config python |
下面是具体运行的例子，首先本人将上面的命令行保存到
.bashrc
文件中。
1 |
qli@lhz:~$ grep pversion .bashrc |
当前版本是python2.7，可以正常运行p4vasp。
1 |
qli@lhz:~$ pversion |
当前的python版本，会在最左侧那一列用
*
标出来。换成python3的话，输入 左侧所对应的数值就OK了。
1 |
qli@lhz:~$ python --version |
顺利将python切换到版本3，但p4vasp这时候就不能用了。
希望p4vasp赶紧出基于python3的版本，如果不出的话，大家就先这样凑活着用吧。如果牛牛们还有更好的解决办法，可以将教程稍微整理下，发送到邮箱[email protected]，本人会继续更新到本节当中。
前面一节，我们讲解了对于扩散这一类反应中，nebmake.pl这个命令生成IMAGES时候的一个坑。通过学习，要了解到一些脚本或者程序本身存在一些缺陷，我们在使用的时候，要避免盲目相信，直接把脚本的结构直接拿来用。总之，原则就是尽可能获取一些具有明确物理化学意义的，比较理想的初始构型。本节，讲解一下旋转过程中，nebmake.pl的坑。这里与其说是坑，不如说是nebmake.pl不适用的情况。因为涉及到旋转过程的时候，一般得到的IMAGES的结构都不咋地，需要自己认真检查，微调下结构或者重新搭结构。
这里我们讲一个极端的例子：乙烷分子的旋转
从一个交叉式的构象到另一个交叉式的构象，两个甲基绕C—C键旋转120°。如果每隔15°插一个点， 正常的结果应该如下图所示:
但是，当使用nebmake.pl这个方法产生的IMAGES结构如下图：
可以看到，
1）初始结构中，并没有旋转的效果，H原子走的是直线的路径。
2）IMAGE中，H和C原子的距离非常小，仅仅贴在一起了，显然这样的结构非常不理想。
如果使用这样的结构进行计算，第一个离子步结束后，计算出来原子间的作用力很强，会导致后面计算中分子直接散架，而这些散架的结构通常都不收敛，如果不及时检查结构，及时杀死，它会在服务器上一直就这样算着，而你还在傻傻地啃着西瓜，聊着QQ，等待结果。
解决办法：
知道有这个坑之后，怎么躲就好办多了。
1）对于类似的旋转结构，自己手动搭建；
2）使用其他的高级点的生成IMAGES的方法，例如IDPP。（ https://wiki.fysik.dtu.dk/ase/tutorials/neb/idpp.html）
3）自己写脚本实现旋转的过程。
这两节简单介绍了一下生成NEB计算IMAGES需要注意的地方，不管咋地，原则还是要再啰嗦一遍：尽可能得到具有理想物理化学意义的初始结构。毕竟好的开始是成功的一半。有兴趣的可以算一下这两节例子，加深一下自己的印象。有大佬公众号留言说贴自己的代码教程怒怼IDPP，希望看到的可以联系俺（[email protected]）。
最近整理材料，发现遗落在角落里的一个python脚本（
cssm.py
），可以批量切常见金属稳定表面的slab模型。
cssm.py
是：
cleave_stable_surfaces_from_metals
的缩写。废话不多说，下面是在天河II上运行的实例，大家照着敲一遍应该问题不大。
1 |
iciq-lq@ln3:/THFS/home/iciq-lq/cssm$ ls |
注：直接用
python cssm.py
会出错，因为默认的是
python2
版本，而ASE基于
python3
，换成
python3 cssm.py
就OK了。
脚本内容:
1 |
#!/usr/bin/env python3 |
简单介绍:
首先声明：最近很忙，没时间写脚本的详解， 如果已经安装ASE的或者有天河II号账号的，可以直接试试，如果不行，自己慢慢琢磨；如果不愿意琢磨，那么就按照前面练习中的笨方法用Material Studio切面。下面以
Ag
为例：
dict_metals
这个字典里面是本人计算的一些常见金属的bulk晶格常数，能量，bulk中原子数目等信息；仅供大家消遣，测试。如果你要通过这个脚本切用于发表文章的slab，运行前，先自己优化一遍bulk，然后修改成自己的结果再运行。
1）
POSCAR_Ag_1
，
POSCAR_Ag_2
，
POSCAR_Ag_3
分别是
1x1
,
2x2
,
3x3
的slab模型。
2）ASE默认把结构放在slab的中心，加真空层的时候也是在两侧加，因此脚本里面是7.5$\AA$，对应15$\AA$；
1 |
slab = fcc111(metal, a = lattice_a, size = (i, i, 4), vacuum = 7.5) |
是切面的关键行；也是脚本的核心部分；
4）每个
slab
有4层，如果想切厚点，修改
size = (i, i, 4)
中的4换成为你想要的层数。如果想切（5x5)，4层的
slab
，将前面一行的
for i in range(1,4)
改成
for i in range(1,6)
；
5) 通过
1 |
constraint_l = FixAtoms(indices=[atom.index for atom in slab if atom.index < i*i*2]) |
把
slab
的
top 2
层放开；
1 |
ase.io.write(name_out, slab, format='vasp') |
是通过ASE将结构输出到
vasp
的格式。
7）ASE的输出一般没有元素行，脚本里面自动加上了：
1 |
subprocess.call(['sed -i ' + '\'5a' + metal + '\' ' + name_out], shell = True) |
8）
bottom(name_out)
功能是把模型从中心拽到底部，输出的结果是
POSCAR_Ag_1_bottomed``，POSCAR_Ag_2_bottomed
，
POSCAR_Ag_3_bottomed
。（个人不喜欢ASE把结构放在中间）
9) 脚本感觉有些啰嗦，欢迎有兴趣的人改进，让它变得更加简洁易懂。
10）最后的bash命令，将所有的
bottomed
归类到各自对应的文件夹中,喜欢把模型放在中间的，可以不用运行
bottom
那个命令。
11）期待我们中国的大牛们写出比ASE更好的软件包出来，安利一波：有问题首先尝试
VASPKIT
。
被批评公众号长草了，简单介绍一下ASE转换文件的小技巧，顺便锄下草。顺利学会下面操作的前提:
1) 已经在自己电脑上安装好了ASE。
2) 电脑里面有POSCARtoolkit.py脚本。见（ https://www.bigbrosci.com/2018/11/11/ex73/）
终端操作
1 |
ase gui file_from -o POSCAR |
1
ase gui
是我们转化文件格式的命令。
2
file_from
是你想要转换的文件，可以为cif, xyz, xsd等。ASE支持很多格式文件的读取和输出，具体的大家可以参考下面网址：
a) https://wiki.fysik.dtu.dk/ase/ase/io/io.html#module-ase.io
b） https://wiki.fysik.dtu.dk/ase/_modules/ase/io/formats.html
3
-o
代表输出文件，这里我们写的是POSCAR，也可以写成CONTCAR。
1）想转换成xyz文件，改成 -o XXX.xyz;
2）想转换成cif文件，改成 -o XXX.cif;
想转什么文件就写什么文件；
想转啥格式的就写啥后缀。
如果失败的话，可能不支持，那就接着用自己的老方法吧。
个人使用经验：
本人的计算模型经常在
xyz
和
POSCAR
两个格式之间切换。下面分析下ASE转换的优点和缺点。
1）从POSCAR到xyz文件（优点）：
在xyz文件的第二行中，ASE会写入一些体系的周期性的信息，这一点非常棒。方便从xyz再转换回POSCAR或者CONTCAR。
1 |
qli@bigbro:~/Desktop/test$ ls |
Lattice=XXX
这一行可以帮助你从xyz转回POSCAR。另外，转换的时候，会出现一个warning信息，装看不见就行了。
2）从xyz到POSCAR
1 |
qli@bigbro:~/Desktop/test$ ase gui c2h4.xyz -o POSCAR |
从xyz转化成POSCAR的时候，ASE的优点：
1）按照Cartesian格式输出。
2）元素行在第一行。
1）输出文件还是老式的VASP格式（vasp 4.6）。第6行前面少了元素哪一行。
解决办法：ASE把这一行写到第一行了。我们直接在vim里面，把第一行复制到第6行前面就OK。
1 |
1 C H |
缺点2）少了固定元素的哪一行，以及固定坐标的’T’ 和 ‘F’ 。
解决办法：通过前面介绍的
POSCARtoolkit.py
脚本来解决，运行下面的命令，输入固定的层数即可。具体用法见（
https://www.bigbrosci.com/2018/11/11/ex73/）
1 |
qli@bigbro:~/Desktop/test$ ls |
小结：
推荐文章
|
|
留胡子的水龙头 · 黑龙江大学计算机与大数据学院（网络安全学院） 2 月前 |
|
|
腼腆的八宝粥 · 人工智能抢走饭碗？美国科技公司开年已裁员万人_中国网 2 年前 |
|
|
霸气的骆驼 · 征服14座高峰：凡事皆可能 - 知乎 2 年前 |
|
|
乐观的火锅 · 何超琼直言自己爱惜许晋亨 证实确未办离婚手续-搜狐娱乐 2 年前 |
|
|
逃课的风衣 · 如何评价《极限挑战》第三季第九期「暗花」? - 知乎 2 年前 |