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

计算化学公社

 找回密码
 现在注册!
查看: 3908|回复: 10

[VASP] VASP 固定基矢优化结构方法

[复制链接]

4

帖子

0

威望

343

eV
积分
347

Level 3 能力者

发表于 2018-3-21 11:21:24 | 显示全部楼层 |阅读模式
本帖最后由 chempeng 于 2018-3-21 19:41 编辑

一直以来在本站学习颇多,感谢 sob 与大佬们的答疑解难。此贴为经验整理总结帖,希望对萌新们(我算一个~)有所帮助,如有疏漏还请各位批评斧正

一般对于三维材料的结构晶格参数优化有两种方法(教程参见 Learn VASP The Hard Way):
1. Birch-Murnaghan 状态方程拟合 http://www.bigbrosci.com/newsitem/278011633
2. VASP 计算中调节 ISIF 参数 http://www.bigbrosci.com/newsitem/278019758
而上述方法对于优化表面,二维或者一维体系不适用,优化此类材料时,须固定一个或两个晶胞方向不优化。一般可通过晶格参数扫描再拟合函数实现,引侯柱峰老师对该类问题解答回帖。
如果二维参数里有一个参数(也就是a,b相等),那通过固定 c,改 a 或 b 的值,对每一个 a 值下的结构进行 ISIF=2 的计算,得到总能,然后根据一系列的 总能-晶格常数 a 的值进行一元二次函数的拟合,得到最小值的位置。
如果是含两个常数(即 a,b 不等),通过地,固定 c,改变a和b地值,对每一组 a 和 b 值下的结构进行 ISIF=2 的计算,得到总能,然后根据一系列的 总能-晶格常数 a 和 b 的值,进行二元二次函数的拟合得到最小值的位置。

这样的方法普适性很高,但却是繁琐。于是有了在特定条件下,重新编译 vasp 使得该优化简化的方法,即当材料为正交晶系或单斜晶系时,可通过修改 constr_cell_relax.F 编译软件实现固定基矢进行结构优化。
对于 VASP 5.4.1 之后的版本,软件的安装编译已经方便了很多,本文以 intel 编译器 + mkl + impi 编译 VASP 5.4.1 为例进行介绍。VASP 的安装编译无需 root 权限,普通用户在家目录即可安装。

1. 查看编译器,mkl 等安装及环境变量设置正确
  1. export LD_LIBRARY_PATH=/share/intel/composer_xe_2015.6.233/mkl/lib/intel64/:$LD_LIBRARY_PATH
  2. export LD_LIBRARY_PATH=/share/intel/composer_xe_2015.6.233/compiler/lib/intel64:$LD_LIBRARY_PATH
  3. export PATH=/share/intel/composer_xe_2015.6.233/bin/intel64:$PATH
  4. export PATH=/share/intel/impi/5.0.3.049/intel64/bin:$PATH
复制代码
2. 修改 constr_cell_relax.F
以固定 c 方向为例,constr_cell_relax.F (src 目录下)修改后如下,显示部分修改内容
  1.       SUBROUTINE CONSTR_CELL_RELAX(FCELL)
  2.       USE prec
  3.       REAL(q) FCELL(3,3)
  4.       DO I=1,3
  5.         FCELL(3,I)=0
  6.         FCELL(I,3)=0
  7.       ENDDO
  8. !     just one simple example
  9. !     relaxation in x directions only
  10. !      SAVE=FCELL(1,1)
  11. !      FCELL=0   ! F90 style: set the whole array to zero
  12. !      FCELL(1,1)=SAVE
  13. !     relaxation in z direction only
  14. !      SAVE=FCELL(3,3)
  15. !      FCELL=0   ! F90 style: set the whole array to zero
  16. !      FCELL(3,3)=SAVE

  17.       RETURN
  18.       END SUBROUTINE
复制代码
对于其他方向晶格基矢的修改同理:对于a方向基矢,将 FCELL(3,I) 和 FCELL(I,3) 分别改为 FCELL(1,I) 和 FCELL(I,1) ;对于b方向基矢,则分别改为 FCELL(2,I) 和 FCELL(I,2) ;固定两个基矢则应该同时将两个方向对应的矩阵元设置为0。或参考 VASP并行可执行软件包,可对晶胞参数进行部分优化(自行搜索) 可实现添加参数文件实现不同基矢方向固定。(未做测试)
3. 编译
建议重新解压源码编译软件
  1. [yourname@localhost ~]$ tar -zxvf vasp.5.4.1.05Feb16.tar.gz
  2. [yourname@localhost ~]$ tar -zxvf vasp.5.lib_.tar_2.gz
  3. [yourname@localhost ~]$ cd vasp.5.4.1
  4. [yourname@localhost ~]$ cp arch/makefile.include.linux_intel makefile.include
  5. [yourname@localhost ~]$ make std
复制代码
普通用户在家目录安装编译前在 makefile.include 文件中添加 $MKLROOT,位置如下
  1. MKLROOT    =/share/intel/composer_xe_2015.6.233/mkl
  2. MKL_PATH   = $(MKLROOT)/lib/intel64
复制代码
编译结束后,进入 bin 目录,可将 vasp_std 重命名为 vasp_fixZ ,将软件移入系统 bin 目录下或将该目录添加至环境变量。
  1. [yourname@localhost ~]$ echo 'export PATH=/home/yourname/vasp.5.4.1/bin:$PATH' >> ~/.bashrc
  2. [yourname@localhost ~]$ source ~/.bashrc
复制代码
完成编译后,在 INCAR 文件中设置 ISIF=3,便可实现固定 Z 轴的结构优化。
修改 constr_cell_relax.F 资料来源科大 @用 在遇见大师兄群内的答疑,另附有 VASP 基矢矩阵的更新流程

评分

参与人数 2eV +13 收起 理由
sobereva + 8
obaica + 5 谢谢

查看全部评分

Blog: http://chempeng.coding.me/chempeng/ or https://chempeng.github.io/

4

帖子

0

威望

13

eV
积分
17

Level 1 能力者

发表于 2018-11-13 09:38:20 | 显示全部楼层
学长你好,请问你帖子里发的优化二维材料用侯老师说的方法对a=b改变晶格常数用ISLF=2优化是否要用到脚本,如果手动改是不是太麻烦了。

69

帖子

1

威望

667

eV
积分
756

Level 4 (黑子)

发表于 2018-11-13 22:11:58 | 显示全部楼层
楼主的方法只能限制Z方向,可以简单写一个3*3的矩阵文件OPTCELL限制任意的晶格矢量。
      INTEGER ICELL(3,3)
      INQUIRE(FILE='OPTCELL',EXIST=FILFLG)
      IF (FILFLG) THEN
         OPEN(67,FILE='OPTCELL',FORM='FORMATTED',STATUS='OLD')
         DO J=1,3
            READ(67,"(3I1)") (ICELL(I,J),I=1,3)
         ENDDO
         CLOSE(67)
         DO J=1,3
         DO I=1,3
            IF (ICELL(I,J)==0) FCELL(I,J)=0.0
         ENDDO
         ENDDO
      ENDIF

4

帖子

0

威望

343

eV
积分
347

Level 3 能力者

 楼主| 发表于 2018-11-14 11:37:34 | 显示全部楼层
Superconduct-BB 发表于 2018-11-13 09:38
学长你好,请问你帖子里发的优化二维材料用侯老师说的方法对a=b改变晶格常数用ISLF=2优化是否要用到脚本, ...

是的,一个 for 循环脚本就可以了。脚本可以参考候博的 VASP 软件包的使用入门指南
Blog: http://chempeng.coding.me/chempeng/ or https://chempeng.github.io/

4

帖子

0

威望

343

eV
积分
347

Level 3 能力者

 楼主| 发表于 2018-11-14 11:48:41 | 显示全部楼层
啦啦黑还黑 发表于 2018-11-13 22:11
楼主的方法只能限制Z方向,可以简单写一个3*3的矩阵文件OPTCELL限制任意的晶格矢量。
      INTEGER ICELL ...


补充提示下,计算前编写输入文件 OPTCELL ,对应格矢,1为优化,0为固定,如下
  1. $ cat OPTCELL
  2. 100
  3. 110
  4. 000
复制代码
Blog: http://chempeng.coding.me/chempeng/ or https://chempeng.github.io/

4

帖子

0

威望

13

eV
积分
17

Level 1 能力者

发表于 2018-11-14 13:56:50 | 显示全部楼层
chempeng 发表于 2018-11-14 11:37
是的,一个 for 循环脚本就可以了。脚本可以参考候博的 VASP 软件包的使用入门指南

与侯老师的脚本有所不同。应该是直接把变量i×到a,b基失里吧,毕竟c不用变。

4

帖子

0

威望

13

eV
积分
17

Level 1 能力者

发表于 2018-11-14 21:16:14 | 显示全部楼层
chempeng 发表于 2018-11-14 11:37
是的,一个 for 循环脚本就可以了。脚本可以参考候博的 VASP 软件包的使用入门指南

学长能否发下你用来循环优化的脚本我刚才试着写了个,发现似乎把变量成到POSCAR的基失里并没有效果。

4

帖子

0

威望

343

eV
积分
347

Level 3 能力者

 楼主| 发表于 2018-11-16 16:53:08 | 显示全部楼层
Superconduct-BB 发表于 2018-11-14 21:16
学长能否发下你用来循环优化的脚本我刚才试着写了个,发现似乎把变量成到POSCAR的基失里并没有效果。

参考修改下
  1. #!/bin/sh
  2. rm WAVECAR
  3. for i in   1.00 1.02 1.04 1.06 1.08 1.10 1.12
  4. do
  5. rm POSCAR
  6. j=$(echo "$i*3" |bc)
  7. cat > POSCAR <<!
  8. 1-1
  9. 1.0
  10.         $j         0.0000000000         0.0000000000
  11.         0.0000000000        $j         0.0000000000
  12.         0.0000000000         0.0000000000         20
  13. xxx
  14. xxx
  15. xxx
  16. !
  17. echo "c = $i Ang"
  18. mpirun -np 8 vasp > out </dev/null
  19. cp CONTCAR CONTCAR-$i
  20. E=`grep "TOTEN" OUTCAR | tail -1 | awk '{printf "%12.6f \n", $5 }'`
  21. V=`grep "volume" OUTCAR | tail -1 | awk '{printf "%12.4f \n", $5 }'`
  22. cat OUTCAR | awk '$1=="in"&&$2=="kB" {print $3,$4,$5}' >> stress
  23. echo $i $V $E  >> energy-strain
  24. done
复制代码
Blog: http://chempeng.coding.me/chempeng/ or https://chempeng.github.io/

4

帖子

0

威望

13

eV
积分
17

Level 1 能力者

发表于 2018-11-17 19:27:59 | 显示全部楼层

多谢

2

帖子

0

威望

39

eV
积分
41

Level 2 能力者

发表于 2019-3-26 10:24:37 | 显示全部楼层
啦啦黑还黑 发表于 2018-11-13 22:11
楼主的方法只能限制Z方向,可以简单写一个3*3的矩阵文件OPTCELL限制任意的晶格矢量。
      INTEGER ICELL ...

请问一下,您这种矩阵写法的话,是不是不会限制晶胞类型啊?比我下面这种写法的优势在哪里呢?
     SUBROUTINE CONSTR_CELL_RELAX(FCELL)
      USE PREC
      REAL(Q) FCELL(3,3), SAVE(3)
      LOGICAL FILFLG
      INTEGER ICELL(3)
      INQUIRE(FILE='OPTCELL',EXIST=FILFLG)
      IF (FILFLG) THEN
        OPEN(67,FILE='OPTCELL',FORM='FORMATTED',STATUS='OLD')
          READ(67,"(3I1)") (ICELL(I),I=1,3)
        CLOSE(67)
        DO I=1,3
          SAVE(I)=FCELL(I,I)
        ENDDO
        FCELL=0.0d0
        DO I=1,3
          IF (ICELL(I)==1) FCELL(I,I)=SAVE(I)
        ENDDO
      ENDIF

2

帖子

0

威望

39

eV
积分
41

Level 2 能力者

发表于 2019-3-26 12:14:53 | 显示全部楼层
请问一下,您测试过三斜体系吗?其中gamma角度和c的大小会受影响发生改变吗?我测试出来c轴大小不变,但是gamma角度会发生变化...您有什么办法固定gamma角度不变吗?其中OPTCELL设置为110,110,000
您需要登录后才可以回帖 登录 | 现在注册!

本版积分规则

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

GMT+8, 2019-12-6 17:20 , Processed in 0.153695 second(s), 24 queries .

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