Ex74 过渡态的计算（五）

发表于 更新于 分类于 阅读次数： Disqus：

前面2节我们讲了如何通过vtst的脚本进行插点,得到过渡态计算中特定数目的images，以及一个非常实用的小脚本来固定slab模型中底层的原子。本节我们就准备一下过渡态计算中：POSCAR、KPOINTS和POTCAR文件的准备工作。

关于Images

前面我们插了8个点，得到10个文件夹，分别为00， 01， 02，….09

细心的你会发现，每个文件夹中都会有一个POSCAR文件。00 中的POSCAR对应的是nebmake.pl IS FS 8命令中的IS，两个一模一样，同理，09 中的POSCAR是nebmake.pl IS FS 8命令中的FS。01到08位我们设想的一个反应路径。（通过p4vasp 0*/POSCAR 这个小技巧，来查看整个路径是不是合理）所有IMAGES中的元素行都是一样的，所以我们任取一个用来生成POTCAR文件。

结合之前VASP的使用经验，每个文件夹中都有一个POSCAR，那么算NEB的时候：是不是各个文件夹也中都要有INCAR， KPOINTS和 POTCAR以及提交任务的脚本呢？

如果你这么想，是很正常的。但我们算过渡态的时候，却不是这样操作的。

实际上，我们将00 到 09 这10个文件夹作为一个整体，看作一个POSCAR。和在这些文件夹相同的目录下， 有INCAR, KPOINTS 和 POTCAR文件。

如果你还不理解，看下面的两个图：首先我们放一个错误的准备工作：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
iciq-lq@ln3:/THFS/home/iciq-lq/LVASPTHW/ex74/wrong$ ls
00  01  02  03  04  05  06  07  08  09
iciq-lq@ln3:/THFS/home/iciq-lq/LVASPTHW/ex74/wrong$ ls *
00:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

01:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

02:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

03:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

04:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

05:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

06:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

07:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

08:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

09:
INCAR  job_sub  KPOINTS  POSCAR  POTCAR

这就是刚刚说的，把INCAR, KPOINTS, POTCAR 和脚本都复制到images对应的文件夹中。知道了错误的准备文件，那下面我们看一下正确的做法:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
iciq-lq@ln3:/THFS/home/iciq-lq/LVASPTHW/ex74/wrong$ cd ../right/
iciq-lq@ln3:/THFS/home/iciq-lq/LVASPTHW/ex74/right$ ls
00  01  02  03  04  05  06  07  08  09  INCAR  job_sub  KPOINTS  POTCAR
iciq-lq@ln3:/THFS/home/iciq-lq/LVASPTHW/ex74/right$ ls *
INCAR  job_sub  KPOINTS  POTCAR

00:
POSCAR

01:
POSCAR

02:
POSCAR

03:
POSCAR

04:
POSCAR

05:
POSCAR

06:
POSCAR

07:
POSCAR

08:
POSCAR

09:
POSCAR
iciq-lq@ln3:/THFS/home/iciq-lq/LVASPTHW/ex74/right$

把0{0..9}这些文件夹作为一个整体，看作POSCAR。这些文件夹 和INCAR，KPOINTS, POTCAR，以及提交任务的脚本放在同一个目录下。

备注：

1） 上面本人插了8个点，纯属自己的日常计算习惯，并不属于什么过渡态的技巧。大家根据自己服务器的运算能力，设置合理的数目。

2） POTCAR与这些images中的元素数目要对应，这个大家应该都会了。

3） KPOINTS: 划重点

如果你从网上搜索相关的资料，很多时候会发现大家经常说，先粗算一下（或者预优化一下），然后再逐渐提高精度。算过渡态的时候也是这样子的，我们不能100%确认自己预想的反应路径（Images文件对应的结构）就是合理的，或者刚刚插的这些点结构都很粗糙，直接用高精度的计算会造成资源的大量浪费，也就是费力不讨好。再严肃点，纳税人的钱不是让你这样随随便便烧的。所以，我们要掌握一些合理的技巧，争取用最少的计算资源，最短的时间，来获取最合理的过渡态结构。前面说到了粗算，到底什么是粗算？粗算需要改什么参数？

这里我们需要了解芝麻和西瓜的故事。

通过前面的学习，我们知道KPOINTS以及slab（计算体系）的大小对于计算的影响最大，它们就是西瓜。

而其他的参数，比如收敛标准，阶段能的大小等，都是小芝麻（当然，你收敛精度设置地极其高那也是个别情况，对于粗算的时候我们不考虑）。

说白了， 就是要抓住影响机时的关键因素。所以粗算的时候，我们需要改的主要有2个地方: KPOINTS 和 slab的模型（集中力量啃西瓜）。

Slab的修改难度较大， 对于新手不太实用，老司机也经常翻车，故先不做讨论。这里我们主要介绍通过调节KPOINTS的一个技巧来实现快速粗算的目的。

先声明：

以下适用于slab表面上的反应计算。如果你算的是体相里面的，比如原子的扩散等，只可作为参考，不能照搬。

在粗算的时候，建议使用低密度的的K点。上来就用高密度的K点计算，即使是土豪，也是不建议的。那么K点密度怎么个低法呢？ 直接来个极端就好了，用gamma点。前面我们还讲过，如果用gamma点计算的时候，需要把slab的原子固定住。

为什么呢？

你可以使用4层的p(2x2)的Cu(111)的slab模型，放开上面2层，然后gamma点优化一个CO在上面的吸附,最后会发现表层的Cu原子驰豫得如同波浪一般起伏。也就是用gamma点会导致表层原子过分地驰豫，从而造成错误的结算结果。

但这也不是绝对的，如果你的slab在xy方向的尺寸很大，使用gamma点的时候满足我们之前讲过的 k*a 的经验规则，即使不固定原子，也是很安全的。（为保险，需要自己测试）如果slab在xy方向的尺寸很小，使用gamma点的就一定需要把原子固定住。固定slab的时候，这里面就有很多窍门了。

A）首先给大家介绍一个最简单的办法。也就是上一节的内容，使用塔马斯的POSCARtoolkit.py脚本，固定和放开任意slab里面的底层原子。比如，目前我们的例子，4层的Ru原子上面，H原子的扩散。我们就可以通过一个for循环来实现slab原子的快速固定。为了方便大家，写了一个基于POSCARtoolkit.py 的bash脚本，该脚本命名为fix.sh, 内容如下：

1
2
3
4
5
6
7
8
9
10
11
12
#!/usr/bin/env bash
## To use it
## fix.sh 4 POSCAR  will fix the bottom 4 layers of  POSCAR  in all image filders
## fix.sh 2 CONTCAR will fix the bottom 2 layers of  CONTCAR in all image folders
## The old POSCAR or CONTCAR will be replaced.

for i in 0* ; do
 cd $i;
 echo "$1" | POSCARtoolkit.py -i "$2" -f;
 mv  ${2}_C "$2"
 cd $OLDPWD
done

注意：

1) 这里大师兄将默认的阈值设置为0.5 Å（修改的脚本）。因为原来1.5 Å有点大， 使用的时候会将表面的H和表层的Ru看作一层。

运行到这里的时候， 我们将所有Images中POSCAR的slab对应的原子全部固定住了。也就是粗算的POSCAR设置好了。

2) 不要偷懒直接复制代码就运行，由于网页上的格式问题，可能直接复制过来运行的时候回出错。

B) 除了使用脚本外，还有另外一个方便的办法：

在课题刚刚开始的时候，就把slab模型中的原子按照层数（也就是z方向的大小）进行排列，前面我们讲过怎么sort坐标了，这里就不再啰嗦。

排列后的一个优势就是各个层的原子在POSCAR中的序号是连续的，可以很方便地进行选择然后通过使用sed 命令将POSCAR里面的T T T 批量转化为 F F F.

比如本例中：10-45为Ru原子，10-18， 19-27, 28-36， 37-45 从上往下数，分别对应第一、二、三、四层的Ru原子。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
Ru   H
   1.00000000000000
   8.1377999784000004   0.0000000000000000   0.0000000000000000
   4.0688999892000002   7.0475415119999996   0.0000000000000000
   0.0000000000000000   0.0000000000000000  21.5631999968999999
Ru   H
36  1
Selective dynamics
Direct
  0.2212968850718230  0.2239409561475650  0.2936169649913790 F F F
  0.2210036114118170  0.5561663540453879  0.2936884876706930 F F F
  0.2212954084343650  0.8880992570872690  0.2936171801115680 F F F
  0.5547059120221930  0.2208073369214740  0.2945957656788490 F F F
  0.5547051050620690  0.5578229108467880  0.2945953069057910 F F F
  0.5572530259134471  0.8880427972121010  0.2935860240494040 F F F
  0.8910551675972660  0.2211397177458350  0.2944546614649750 F F F
  0.8893412673473050  0.5559102687354450  0.2935156586054220 F F F
  0.8893404937786400  0.8880862995751930  0.2935168781920150 F F F
  0.1119432861914480  0.1106972088633780  0.1977748273385700 F F F
  0.1112675774095350  0.4441060896203480  0.1977506775579700 F F F
  0.1112675246777360  0.7779649556723670  0.1977514668340350 F F F
  0.4440066349954530  0.1106943439331120  0.1977106006833940 F F F
  0.4446549772136020  0.4443417041956370  0.1981365664070010 F F F
  0.4440061461204080  0.7786377563630820  0.1977110173547200 F F F
  0.7777009262273890  0.1113226147079440  0.1980782640541890 F F F
  0.7777006933574599  0.4443146361556100  0.1980780083860340 F F F
  0.7774653980933121  0.7779366397494840  0.1977260420610920 F F F
  0.2222233413377670  0.2222233413358840  0.0992873228095945 F F F
  0.2222233413377670  0.5555566746692190  0.0992873228095945 F F F
  0.2222233413377670  0.8888900080025550  0.0992873228095945 F F F
  0.5555566746710950  0.2222233413358840  0.0992873228095945 F F F
  0.5555566746710950  0.5555566746692190  0.0992873228095945 F F F
  0.5555566746710950  0.8888900080025550  0.0992873228095945 F F F
  0.8888900080044310  0.2222233413358840  0.0992873228095945 F F F
  0.8888900080044310  0.5555566746692190  0.0992873228095945 F F F
  0.8888900080044310  0.8888900080025550  0.0992873228095945 F F F
  0.1111100016643290  0.1111100016603930  0.0000000000000000 F F F
  0.1111100016643290  0.4444433349937280  0.0000000000000000 F F F
  0.1111100016643290  0.7777766683270571  0.0000000000000000 F F F
  0.4444433349976650  0.1111100016603930  0.0000000000000000 F F F
  0.4444433349976650  0.4444433349937280  0.0000000000000000 F F F
  0.4444433349976650  0.7777766683270571  0.0000000000000000 F F F
  0.7777766683310010  0.1111100016603930  0.0000000000000000 F F F
  0.7777766683310010  0.4444433349937280  0.0000000000000000 F F F
  0.7777766683310010  0.7777766683270571  0.0000000000000000 F F F
  0.6419725686006360  0.3456803556224560  0.3432239557346220 T T T

在上面POSCAR基础上，放开上面的2层：

1
sed –i  '10,27s/F/T/g’ 0*/POSCAR

固定所有的Ru原子，可以：

1
sed –i  ’10,45s/T/F/g’ 0*/POSCAR

大家随意练习折腾，本人也相信，如果从本书的Ex-0一直认真学下来，如果sed 命令还不会的话，说明你真的没有认真在学习。本节的例子以及脚本下载链接：

https://pan.baidu.com/s/1qqySBpqc0zMu1lvpLVZuWA 提取码: 3ryx

总结

本节，我们介绍了准备文件中POSCAR， KPOINTS，和POTCAR，以及他们在目录里面的排列方式。为减轻计算压力，过渡态可以先粗算一下， 通过使用gamma点结合固定slab的方式来减少计算量。通过POSCARtooklit.py或者sed命令实现批量固定slab的方法。