计算化学公社

标题: VASP达到收敛标准，但计算仍在继续 [打印本页]
作者
Author: 种花家的蓝兔    时间: 2022-9-28 19:54
标题: VASP达到收敛标准，但计算仍在继续
各位老师前辈好，我在OUTCAR文件中看到TOTAL-FORCE (eV/Angst)下边的数值已经小于|EDIFFG|，以及dE<EDIFF，这应该是已经达到收敛标准了吧，为啥计算仍在计算没有停止呢？望各位老师不吝赐教
作者
Author: jing@lll    时间: 2022-9-28 20:17
grep RMS OUTCAR 看一下
作者
Author: 种花家的蓝兔    时间: 2022-9-28 20:32
jing@lll 发表于 2022-9-28 20:17
grep RMS OUTCAR 看一下

啥也没有
作者
Author: 含光君    时间: 2022-9-28 21:54
是不是这原子没收敛来着
(, 下载次数 Times of downloads: 31)
作者
Author: 种花家的蓝兔    时间: 2022-9-28 22:04
含光君 发表于 2022-9-28 21:54
是不是这原子没收敛来着

看力收敛，是看上边这些数值<|EDIFFG|吗，如果是的话，我的EDIFF设置的是-0.02，您画的那一行满足小于0.02吧
作者
Author: 含光君    时间: 2022-9-28 22:43
种花家的蓝兔 发表于 2022-9-28 22:04
看力收敛，是看上边这些数值

VASPwiki（https://www.vasp.at/wiki/index.php/EDIFFG）上没有说明是每个方向的力还是合力，我个人理解是每个离子所受合力达到收敛标准才行。以及，在VASPwiki上一个实例（https://www.vasp.at/wiki/index.php/H2O）上出现如下表述：
(, 下载次数 Times of downloads: 36)
从表达式看是对每一个离子判定收敛条件，所以我想大概是要达到合力收敛才行？以上内容仅为个人理解，如有不正确之处欢迎大佬指正。
作者
Author: 卡开发发    时间: 2022-9-29 07:58
种花家的蓝兔 发表于 2022-9-28 20:32
啥也没有

不要看RMS的Force，应该看Max的Force，但只有装了Henckelman他们的VTST插件这两个值才会打印到OUTCAR。

作者
Author: 卡开发发    时间: 2022-9-29 08:02
本帖最后由 卡开发发 于 2022-9-29 08:59 编辑
含光君 发表于 2022-9-28 22:43
VASPwiki（https://www.vasp.at/wiki/index.php/EDIFFG）上没有说明是每个方向的力还是合力，我个人理解 ...

应该是看每一个分量的，也就是N个forces（F=sqrt(Fx^2+Fy^2+Fz^2)）当中的最大的那一个分量的绝对值要≤-EDIFFG。LZ给出的力信息不完善其实没办法判断，另外有固定的原子也应该被排除。
作者
Author: 种花家的蓝兔    时间: 2022-9-29 08:24
卡开发发 发表于 2022-9-29 08:02
应该是看每一个分量的，也就是3N个forces（N个离子*3坐标）当中的最大的那一个分量的绝对值要≤-EDIFFG。 ...

那老师可以帮忙看看吗，我这个应该是达到了收敛标准，为啥计算仍在继续呢
作者
Author: 卡开发发    时间: 2022-9-29 08:37
种花家的蓝兔 发表于 2022-9-29 08:24
那老师可以帮忙看看吗，我这个应该是达到了收敛标准，为啥计算仍在继续呢

你不如发完整OUTCAR上来？
作者
Author: 卡开发发    时间: 2022-9-29 08:55
种花家的蓝兔 发表于 2022-9-29 08:24
那老师可以帮忙看看吗，我这个应该是达到了收敛标准，为啥计算仍在继续呢

我去翻了下源码，应该6楼说法是对的，应该是看N个原子的|F|=sqrt(Fx^2+Fy^2+Fz^2)最大的小于EDIFFG才行。
作者
Author: 种花家的蓝兔    时间: 2022-9-29 09:07
卡开发发 发表于 2022-9-29 08:55
我去翻了下源码，应该6楼说法是对的，应该是看N个原子的|F|=sqrt(Fx^2+Fy^2+Fz^2)最大的小于EDIFFG才行。

学到了学到了，谢谢老师的回复
作者
Author: 种花家的蓝兔    时间: 2022-9-29 09:07
含光君 发表于 2022-9-28 22:43
VASPwiki（https://www.vasp.at/wiki/index.php/EDIFFG）上没有说明是每个方向的力还是合力，我个人理解 ...

谢谢您的回复
作者
Author: 含光君    时间: 2022-9-29 10:06
卡开发发 发表于 2022-9-29 08:55
我去翻了下源码，应该6楼说法是对的，应该是看N个原子的|F|=sqrt(Fx^2+Fy^2+Fz^2)最大的小于EDIFFG才行。

谢谢发发老师~
作者
Author: wsz    时间: 2022-9-29 13:00
本帖最后由 wsz 于 2022-9-29 13:04 编辑

vasp的收敛标准就是看所有非固定原子中受力最大的原子的受力，三个方向的模长
之前写过一个脚本，对于没有添加VTST的可以用这个查看受力，但是写的比较粗糙，需要手动修改一下体系原子数（第24行），没有去除固定的原子（需要的话可以自己修改一下）python2和3都可以，只用了内置的math库
作者
Author: 卡开发发    时间: 2022-9-29 13:03
含光君 发表于 2022-9-29 10:06
谢谢发发老师~

没有，我一开始也记错了，后来看你那层提示去翻了一下源码确定了一下。
作者
Author: 种花家的蓝兔    时间: 2022-9-29 19:18
wsz 发表于 2022-9-29 13:00
vasp的收敛标准就是看所有非固定原子中受力最大的原子的受力，三个方向的模长
之前写过一个脚本，对于没有 ...

好的好的，非常感谢您，感谢您的回复
作者
Author: 种花家的蓝兔    时间: 2022-9-29 19:22
本帖最后由 种花家的蓝兔 于 2022-9-29 19:25 编辑
wsz 发表于 2022-9-29 13:00
vasp的收敛标准就是看所有非固定原子中受力最大的原子的受力，三个方向的模长
之前写过一个脚本，对于没有 ...
作者
Author: 卡开发发    时间: 2022-9-29 22:22
本帖最后由 卡开发发 于 2022-9-29 22:26 编辑

固定的信息可能需要通过CONCAR或者通过vasprun.xml获得，但是vasprun.xml只能在计算结束时才会写入。固定的原子在固定的分量置零就行了。
做法也有两种：
1、完全使用原生python，我就用了math和re两个模块。
  1. import re
  2. import math

  4. class Data:
  5.         def __init__(self):
  6.                 self.natms=0
  7.                 self.constr=[]
  8.                 self.forces=[]
  9.         def get_poscar_info(self):
  10.                 """
  11.                 vasp5 POSCAR format
  12.                 """
  13.                 with open('POSCAR') as f:
  14.                         data=f.read().lower()
  15.                         data=data.replace('t','1')
  16.                         data=data.replace('f','0')
  17.                         data=data.split('\n')
  18.                         self.natms=sum(int(d) for d in data[6].split())
  19.                         if data[7].strip()=='s':
  20.                                 self.constr=[[int(d) for d in line.split()[3:6]]
  21.                                                 for line in data[9:9+self.natms]]
  22.                         else:
  23.                                 self.constr=[[1 for i in range(3)] for j in range(self.natms)]
  24.         def get_outcar_info(self):
  25.                 with open('OUTCAR') as f:
  26.                         data=f.read()
  27.                         parttern='TOTAL-FORCE \(eV/Angst\)\n'
  28.                         parttern+='\s+--+\n(.*?)\n'
  29.                         parttern+='\s+--+\n'
  30.                         data=re.findall(parttern,data,re.S)
  31.                         self.forces=[[[float(d) for d in line.split()[3:6]]
  32.                                         for line in lines.split('\n')] for lines in data]
  33.         def get_forces(self):
  34.                 self.get_poscar_info()
  35.                 self.get_outcar_info()
  36.                 norm=lambda x:math.sqrt(x[0]**2+x[1]**2+x[2]**2)
  37.                 for s,force in enumerate(self.forces):
  38.                         max_force=max([norm([force[i][j]*self.constr[i][j] for j in range(3)])
  39.                                         for i in range(self.natms)])
  40.                         print('% 4d %12.6f'%(s+1,max_force))
复制代码

2、完全使用生态，如ase和numpy
  1. from  ase import io
  2. import numpy as np

  4. poscar=io.read('POSCAR')
  5. outcar=list(io.iread('OUTCAR'))
  6. for s,image in enumerate(outcar):
  7.         image.set_constraint(poscar.constraints)
  8.         forces=image.get_forces()
  9.         max_force=np.max(np.linalg.norm(forces,axis=1))
  10.         print('% 4d %12.6f'%(s+1,max_force))
复制代码







欢迎光临 计算化学公社 (http://bbs.keinsci.com/) Powered by Discuz! X3.3