VASP - Special functions - fix basic vector

A word spoken is past recalling.

一言既出, 驷马难追

Updated time: 01/30 2024.

Ref

  1. VASP固定基矢优化结构最新方法更新 (qq.com)

  2. http://blog.wangruixing.cn/2019/05/05/constr/

  3. 固定晶格部分无缺大师兄贡献:https://zhuanlan.zhihu.com/p/34750103

  4. https://github.com/Chengcheng-Xiao/VASP_OPT_AXIS

  5. VASP 固定基矢优化结构方法 - 第一性原理 (First Principle) - 计算化学公社 (keinsci.com)

  6. VASP固定基矢优化结构方法 | 世事如棋 (shishiruqi.com)

  7. Cd Si volume relaxation - VASP Wiki

  8. LATTICE_CONSTRAINTS - VASP Wiki

1. 固定基矢的作用

对于二维材料,薄膜材料,乃至块体材料的拉伸/压缩模拟中,施加指定应变的基矢方向需要固定,同时由于泊松效应,其他方向的cell长度也会发生相应改变。

除了固定基矢长度的作用外,某些情况下也可以固定基矢方向。

Ref 3 Fixing specific stress tensor element(s)

2. 数学描述以及固定基矢方向

/src/constr_cell_relax.F source code

!-----------------------------------------------------------------------
!
! At present, VASP does not allow to relax the cellshape selectively
! i.e. for instance only cell relaxation in x direction.
! To be more precisse, this behaviour can not be achived via the INCAR
! or POSCAR file.
! However, it is possible to set selected components of the stress tensor
! to zero.
! The most conveninent position to do this is the routines 
! CONSTR_CELL_RELAX  (constraint cell relaxation).
! FCELL contains the forces on the basis vectors.
! These forces are used to modify the basis vectors according
! to the following equations:
!
!      A_OLD(1:3,1:3)=A(1:3,1:3) ! F90 style 
!      DO J=1,3
!      DO I=1,3
!      DO K=1,3
!        A(I,J)=A(I,J) + FCELL(I,K)*A_OLD(K,J)*STEP_SIZE
!      ENDDO
!      ENDDO
!      ENDDO
! where A holds the basis vectors (in cartesian coordinates).
!
!-----------------------------------------------------------------------

      SUBROUTINE CONSTR_CELL_RELAX(FCELL)
      USE prec
      REAL(q) FCELL(3,3)

!     just one simple example
!     relaxation in x directions only
!      SAVE=FCELL(1,1)
!      FCELL=0   ! F90 style: set the whole array to zero
!      FCELL(1,1)=SAVE
!     relaxation in z direction only
!      SAVE=FCELL(3,3)
!      FCELL=0   ! F90 style: set the whole array to zero
!      FCELL(3,3)=SAVE

知乎评论区,找不到出处了,我觉得这里的矩阵表示有错误,这只是变化值

在源代码中交代的很清楚,cell的更新是通过下列循环完成的。FCELL和上图中的F是一样的,只是一个是左乘,一个右乘。

!      A_OLD(1:3,1:3)=A(1:3,1:3) ! F90 style 
!      DO J=1,3
!      DO I=1,3
!      DO K=1,3
!        A(I,J)=A(I,J) + FCELL(I,K)*A_OLD(K,J)*STEP_SIZE
!      ENDDO
!      ENDDO
!      ENDDO

从一开始的只能对正交盒子固定基矢,扩展到单斜盒子固定基矢,再到对三斜盒子固定方向(Ref 5 Theta-Al2O3的优化

3. 操作中的误区

对于 Ref 3 中 Fixing specific stress tensor element(s) 方法重新对 VASP 编译

cd to_vasp_root_directory
cp path/stress_relax_finner.patch ./
patch -p0 < stress_relax_finner.patch

module load mpi/intel/17.0.7-thc
#module load intel/17.0.7-thc

cp arch/makefile.include.linux_intel  ./makefile.include
make std

对于VASP官网例子测试 Ref 6

POSCAR cell
5.5
 0.0    0.5     0.5
 0.5    0.0     0.5
 0.5    0.5     0.0

CONTCAR-000000000
   5.50000000000000     
     0.0000000000000000    0.5000000000000000    0.5000000000000000
     0.5000000000000000    0.0000000000000000    0.5000000000000000
     0.5000000000000000    0.5000000000000000    0.0000000000000000

CONTCAR-100010000
   5.50000000000000     
     0.0000000000000000    0.4958558152381574    0.5000000000000000
     0.4958558152381574    0.0000000000000000    0.5000000000000000
     0.4958558152381574    0.4958558152381574    0.0000000000000000

CONTCAR-110110000
   5.50000000000000     
     0.0000000000000000    0.4958559358684163    0.5000000000000000
     0.4958559358684163    0.0000000000000000    0.5000000000000000
     0.4958559358684163    0.4958559358684163    0.0000000000000000

CONTCAR - turn off IOPTCELL tag
   5.50000000000000     
     0.0000000000000000    0.4971253367745418    0.4971253367745418
     0.4971253367745418    0.0000000000000000    0.4971253367745418
     0.4971253367745418    0.4971253367745418   -0.0000000000000000

算到这里我觉得很奇怪,因为 manual 中写到

IOPTCELL = xx yx zx xy yy zy xz yz zz
for example:

IOPTCELL = 1 1 0 1 1 0 0 0 0
will relax the xx, xy, yx and yy components of the stress tensor, while keeping the other components fixed (to zero).

我以为对应某个分量的值设为0,那 cell 中这部分便被固定住了。但这是错误的!cell的更新是一个矩阵运算,不存在点乘。首先 POSCAR 给出的是正交晶系,但基矢并不是 conventional cell 的常规写法。接下来,即使我们设定

100010000 =》 FCELL = [[1 0 0] [0 1 0] [0 0 0]],经过 FCELL[i, k]·A_OLD[k, j] 的遍历几乎所有的元素都发生了改变。

误区的核心在于:没有理解了矩阵运算

当然,我还关闭对称性算了一下,我以为是对称性和IOPTCELL造成的约束冗余。

4. 固定基矢的作用对象

4.1 正交晶系

A = [[a 0 0] [0 b 0 ] [0 0 c] ]

可以实现固定任意基矢,也可以固定基矢的方向优化。

对于 A ,第一行对应 a1, 二 a2, 三 a3。以对 a1 操作为例, 固定该基矢,则 FCELL 的第一行为 [ 0 0 0] 。固定基矢的方向优化 [1 0 0]

4.2 单斜晶系

A = [[a 0 0] [b1 b2 0 ] [0 0 c] ]

可以实现固定任意基矢, 可以对 a, c 轴固定基矢方向优化,b 轴可以固定在 xy 平面内。

a c 轴同上。b 轴固定基矢,[0 0 0], 固定方向优化 [1 1 0]

4.3 三斜晶系

A = [[a 0 0] [b1 b2 0] [c1 c2 c3]]

a b 轴同上,c 轴固定基矢,[0 0 0],没法固定方向。

5. VASP 6 - LATTICE_CONSTRAINTS Command

Mind: LATTICE_CONSTRAINTS in combination with IBRION=1,2 is available from VASP 6.5.0.

LATTICE_CONSTRAINTS TAG 可用于MD 模拟和静态弛豫。尤其适用于液体以及二维材料的计算。

但是我没有 6.5.0 版本的VASP,所以没有进行测试。

结论

manual 中写的是对的,确实只要设为 1 这个地方就会被优化,设为 0 就会固定。关键在于,A 的形式是否和这里的一致以及,不是因为我们设为 0 而固定,而是得通过矩阵运算检验。例如对于单斜晶系 b 轴,设 [0, 1, 0], b1 仍会被优化

这是我的粗略见解,肯定有疏漏之处。


Please indicate the source when reprinting. Please verify the citation sources in the article and point out any errors or unclear expressions.