一、电子迭代算法 ALGO =Normal selects IALGO =38 (blocked Davidson iteration scheme). ALGO =VeryFast selects IALGO =48 (RMM-DIIS). ALGO =Fast selects a faily robust mixture of the Davidson and RMM-DIIS algorithms. In this case, Davidson ( IALGO =38) is used for the initial phase, and then VASP switches to RMM-DIIS ( IALGO =48). Subsequencly, for each ionic update, one IALGO =38 sweep is performed for each ionic step (except the first one). 二、 LELF = .TRUE. | .FALSE. Default: LELF = .FALSE. Description: LELF determines whether to create an ELFCAR file or not. This file contains the so-called ELF (electron localization function) If LELF is set, NPAR =1 has to be set explicitely in the INCAR file in addition 三、静电势 计算静电势只需加上 LVHAR = .TRUE. LVHAR = Default: LVHAR = .FALSE. Description: This tag determines whether the total local potential (saved in the file LOCPOT ) contains the entire local potential (ionic + Hartree + exchange correlation) or the electrostatic contributions only (ionic + Hartree). This tag is available in VASP.5.2.12 and newer version. Note that in VASP.5.2.12, the default is to write the entire local potential, including the exchange correlation potential. 四、偶极修正 对于slab需要加上偶极修正，尤其氧化物这种极性强的材料。 LDIPOL = .TRUE. IDIPOL = 3
转载自 http://blog.sciencenet.cn/blog-567091-777871.html 1) 自旋极化、旋轨耦合、加U等哪几项需要同时设置 2）只有先自旋极化了，才可能考虑加U或旋轨耦合？ 手册说明： LDAUTYPE=2 (Default): The simplified (rotationally invariant) approach to the LSDA+U, introduced by Dudarev et al. . This flavour of LSDA+U is of the following form: LSDA 局域自旋密度近似，当然要设置ISPIN=2了！！ LDAUTYPE=4: Same as LDAUTYPE=1, but LDA+U instead of LSDA+U (i.e. no LSDA exchange splitting). In the LDA+U case the double counting energy is given by, LDA局域谜底近似，是否意味着可以不设置ISPIN=2，但U和J需单独设置，参考上面的公式 It is important to be aware of the fact that when using the L(S)DA+U, in general the total energy will depend on the parameters and . It is therefore not meaningful to compare the total energies resulting from calculations with different and/or (c.q. in case of Dudarev's approach). Furthermore, since LDA+U usually results in aspherical charge densities at and atoms we recommend to set LASPH = .TRUE . in the INCAR file for gradient corrected functionals (see Sec. 6.44 ). For Ce O for instance, identical results to the FLAPW methods can be only obtained setting LASPH = .TRUE . Note on bandstructure calculation: The CHGCAR file also contains only information up to LMAXMIX (defaulted to 2) for the on-site PAW occupancy matrices. When the CHGCAR file is read and kept fixed in the course of the calculations (ICHARG=11), the results will be necessarily not identical to a selfconsistent run. The deviations can be (or actually are) large for L(S)DA+U calculations. For the calculation of band structures within the L(S)DA+U approach, it is hence strictly required to increase LMAXMIX to 4 (d elements) and 6 (f elements). (see Sec. 6.63 ). 网络摘录： http://emuch.net/html/201004/2007448.html 作者: nkleof vasp计算DOS有一种方法，就是以较低k网格达到收敛，再以很高的k网格通过ICHARG=11读取自洽CHGCAR进行DOS计算。 我昨天进行了这样的计算，但是在计算中进行了+U，如下是INCAR PREC = Accurate LREAL = A IALGO = 48 NPAR = 8 NSIM = 8 ISTART = 1 ICHARG = 11 NELMIN = 4 LORBIT = 12 EDIFF = 1E-5 ISMEAR = -5 EMIN = -20 EMAX = 20 NEDOS = 11000 IS PIN = 2 LDAU = .TRUE. LDAUTYPE = 2 LDAUL = -1 2 -1 LDAUU = 0 2 0 LDAUJ = 0 0 0 LDAUPRINT = 2 其中自洽收敛是以11 11 11网格，非自洽DOS是以21 21 21网格，另外也进行了21 21 21网格的自洽计算进行对比。 如下是三次的DOS 图，分别为k11网格自洽的DOS，在k11网格自洽基础上k21的DOS结果，k21网格自洽的DOS 从图上可以看到，在k11自洽基础上进行k21的ICHARG=11的DOS结果和其他两个有点不同，费米面以上的部分有所变化。 所以，我猜测可能是因为是INCAR里面的+U设置在后来的非自洽计算中重复起了作用， 使得实际作用的U值成了4eV。即原来以U=2进行了自洽计算，得到的CHGCAR已经包含2eV 的Coulomb repulsion，后来进行非自洽计算时INCAR中的+U设置又起了作用，再加上了2eV的Coulomb repulsion，所以出来的DOS结果和原来U=2时的不同。【见上面的手册说明？】 希望有高手来证实一下，或者解释一下k21的非自洽结果与原自洽结果的DOS不同的原因。 /////////////////////////// U=4的计算结果出来了，证明我前面提到的猜测不正确，U=4的DOS和U=2，ICHARG=11的DOS并不相同。其实从上个帖子的DOS图可以看出来， 当U值增大到4时，费米面以上的那部分峰肯定会向着高能量方向移动一段距离，而不是几乎在原地不动。 锐利的碎片 (站内联系TA) Originally posted by nkleof at 2010-04-29 23:00:53: 意思是说在自洽CHGCAR的基础上进行ICHARG=11的计算也需要原来自洽计算时的+U设置？ 对啊，U值最好先算下带隙看下，不同关联势都不一样。一般要求整个过程都应该+U。 网络摘录： http://emuch.net/html/201007/2225276.html J的值的意义与LDAUTYPE有关的 。当LDAUTYPE为2时，U-J的差值才有有意义，即有效的U参数。 至于有效U参数的值怎么选取，通过线性响应的方式来估测，类似pwscf里面的。要么将U值从小测到大，并计算结果同实验结果(比如电子结构等)进行比较，然后选取一个合适的值。再要不然，参考文献中的值。 我看说明书里面说的是：in Dudarev’s approach the parametersU and J do not enter seperately, only the difference (U J) is meaningfull. 并且对于LDAUTYPE， 1 Rotationally invariant LSDA+U according to Liechtenstein et al. 4 Idem 1., but LDA+U instead of LSDA+U (i.e. no LSDA exchange splitting) 2 Dudarev’s approach to LSDA+U (Default) 那是不是说LDAUTYPE=2时只有U-J的值有意义呢？ 我还想请问一下：如果我不打开自旋（用GGA+U算）的话，那么可以设置ISPIN的值为1然后将LDAUTYPE能选择为2吗？还是说非自旋计算时只能将LDAUTYPE设置为4呢？ valenhou001 (站内联系TA) 有关U,J的意义，前一个回帖中说的有点小误，我已经更正过来了。 既然LSDA+U里面U，J参数在目前来说，用的时候，究竟什么样的值才有意义，难以说的清楚，所以尽量减少不明确性，推荐采用默认的方式，即LDAUTYPE=2，要给定的参数就是U-J的值了。 LD AUTYPE=2的效果和LDAUTYPE=1或4的都是差不多的。 valenhou001 (站内联系TA) Originally posted by xiaojie7783 at 2010-07-16 23:56:15: 谢谢，我还想请问：如果我不打开自旋（用GGA+U算）的话，那么可以设置ISPIN的值为1然后将LDAUTYPE能选择为2吗？还是说非自旋计算时只能将LDAUTYPE设置为4呢？ 谢谢 对非自旋选用LDAUTYPE=4，而自旋极化是选用LDAUTYPE=2这样不行 。为了比较自旋极化和非自旋极化的稳定性，最好要么都是用LDAUTYPE=2， 要么LDAUTYPE=1对自旋极化，LDAUTYPE=4对非自旋极化。 百度百科摘录： 在HUBBARD模型一级近似下，U就是考虑了 同一个原子上自旋相反的局域电子之间的库仑排斥，从而导致了能带的“重正化”。 强关联电子体系是指电子间的交互作用不可忽略的系统。在简单的固体理论中，固体中电子之间的静电相互作用被忽略了，不会出现在哈密顿算符里。故各个电子被看成是独立的，不会相互影响。然而，在许多物质中，静电能不能被忽略。当把这一部分能量写入哈密尔顿量，就得到强关联模型（或赫巴德模型(Hubbard model))。在强关联电子体系 ，由于电子之间的强相互作用，导致了许多新奇的物理现象。如高温超导体、二维电子气中的分数量子霍尔效应、锰氧化物材料中的巨磁阻效应、重费米子系统、二维高迁移率材料中的金属－绝缘体相变、量子相变和量子临界现象、一维导体中的电荷密度波等等 。 1937年，科学家就发现NiO，MnO，CoO 等氧化物并不是能带理论所预言的金属，而是能隙很大的绝缘体。Mott 引进了关联能来解释这一物理问题，认为 d 电子间库仑相互作用抑制了极化涨落，产生了关联能隙，后来这一类绝缘体即被称为莫特绝缘体。Mott 进一步讨论了VO2，V2O3等材料因温度或压力改变所引起的绝缘体到金属的相变，认定它们也是电子关联导致的相变，后来被称为 Mott 转变。莫特绝缘体几乎占了3d 过渡金属二元氧化物中的一半，还包括很多的多元复杂氧化物和 4f 稀土化合物及5f 锕系化合物。 钙钛矿结构的锰氧化物是强关联电子体系的一个例子。这类材料的显示出庞磁电阻效应，以及电荷有序、轨道有序、超导序和磁有序．在LaCaMnO系的材料中，加上磁场后的电阻变化率可达到103~106。这种材料的铁磁性的根源是双交换相互作用，而且磁性转变与绝缘体－金属转变相邻近。 重费米子体系是强关联电子体系的另一个例子。在重费米子金属中，存在RKKY相互作用与Kondo相互作用的竞争。RKKY相互作用是局域磁矩之间通过极化的传导电子云而发生的间接交换相互作用。Kondo相互作用是局域磁矩与周围传导电子的直接交换相互作用。在低温下，两种相互作用竞争的结果，使重费米子金属有多种基态：磁有序态、超导态、费米液体态和非费米液体基态等。另一些过渡金属氧化物（如LiV2O5）同样具有典型的重费米子特性。 铜基以及铁基高温超导体同样是强关联电子体系。以BiSrCaCuO 为例，在掺杂浓度x为零的材料是反铁磁序的绝缘体，随着掺杂的增加会发生绝缘体到金属的转变。而在低温就具有超导电性，随着掺杂的增加，Tc 达到一峰值之后，又逐渐下降，高温超导体的正常态的电子性质都十分异常。 参考资料 徐慧．强关联电子体系材料．中南大学出版社，2009-07 LDA+U 中的U最好参考实验值,如果找不到实验值就要自己去试,通常认为带隙或晶格常数不再发生变化时的U值作为近似参考值 我来说一下。 LDA+U或是GGA+U中的U指的是HUBBARD模型中的自旋相反电子的强关联排斥能。通常的DFT计算是单电子近似的，或者是平均场近似，在这种近似下，的确可以解释很多问题。但是，真正的凝聚态体系，分子或是固体，本质上是多体理论。很多问题，比如MOTT绝缘体，在单电子近似下的结论就是错的。因此必须超越单电子近似，考虑体系的量子多体效应。U就是在这种情况下出现的。 简单的说，或者在HUBBARD模型一级近似下，U就是考虑了同一个原子上自旋相反的局域电子之间的库仑排斥，从而导致了能带的“重正化”。 比如所，在普通能带理论中，如果原胞中的电子数是奇数个，则体系一定是金属。但是在考虑过在位格点（ON-SITE）的排斥能后，就出现了绝缘体。 因此，考虑U的情况一般是体系中存在局域化的d,f电子（比如Ce）。这些电子的存在是”标准“的DFT理论无法合理处理的，必须要考虑U的影响。 而通常U是一个试验值，不同的化学环境，同一个原子的U就可能不同，必须通过试验来确定。 呵呵，先说这么多，累死了。 最近在做DFT+U计算，有些收获分享一下。本人数学功底先天不足，难免有理解错误的地方，欢迎指正。 首先是对+U必要性的感悟。DFT（LDA和GGA）对一般体系的结构、能带计算还是很令人满意的。这些一般体系是指前3周期的体系和纯金属体系。但是，对于包含d电子和f电子的体系，特别是过渡金属氧化物或氮化物，在“金属/绝缘体”的定性判断上常常出错。LDA和GGA往往把一些绝缘体/半导体的gap算的太小，甚至最高占据轨道/能带(HOMO)在Fermi面之上，即成金属了。这些最高占据轨道往往是来自这些金属原子的d电子和f电子。那么为什么DFT处理d/f电子会出现较大误差呢？我认为这些d/f电子的特殊性在于它们不仅离原子核较远，轨道空间伸展形状奇特（比如穿透效应，即穿透更加外层的s，p轨道）而且多数会发生自旋极化，即原应自旋相反的电子对发生自旋跃迁，导致金属原子神奇地出现了净的磁矩。如果单胞内所有金属原子的净磁矩方向相同，则单胞总体产生净磁矩，材料呈现铁磁性；如果单胞内相邻金属原子的净磁矩相反，则单胞没有总体净磁矩，但存在一上一下的磁矩排列，材料呈现为反铁磁性。正是d/f电子这些神奇的特性，使得DFT计算的过渡金属氧化物或氮化物能带结构频频失误。为什么呢？原因在于DFT计算电子-电子之间的排斥作用所用的交换-相关泛函，是基于单粒子近似发展起来的。单粒子近似是把自旋配对的电子对看成单个粒子，这样Kohn-Sham轨道数量为总电子数的一半。但是当处理d/f电子时，不得不把电子对分开为自旋向上和自旋向下的两个粒子，即开壳层计算open shell。对开壳层体系使用单粒子近似下的交换-相关泛函时，其中的相关泛函就不能完全考虑d/f电子的相关效应。d/f电子的相关效应，就是自旋电子与自旋电子之间的排斥，和自旋电子与内层电子之间的排斥。这些排斥作用使得轨道/能带之间分割较远，轨道/能带较窄。但是DFT的相关泛函对这些排斥考虑的不够，结果轨道/能带过宽，轨道与轨道相互接近甚至重叠。这就是为什么DFT把绝缘体计算成了金属。一般把这些过度金属半导体/绝缘体体。
此前就想整理这套VASP的教程，一直没有时间。最近看到大师兄的QQ群上有人上传了，感谢匿名网友的整理工作。该教程整理了“ex0: 学习前的准备工作” ~ “Ex72- 过渡态的计算（三）”的内容，值得初学者好好研究。 Learn VASP The Hard Way.part1.rar Learn VASP The Hard Way.part3.rar Learn VASP The Hard Way.part2.rar Learn VASP The Hard Way.part5.rar Learn VASP The Hard Way.part4.rar Learn VASP The Hard Way.part7.rar Learn VASP The Hard Way.part6.rar 由于附件只能上传不超过8M的文件，因此原始PDF文件被压缩为了7个文件单独上传。 把上述7个文件全部下载放在一个文件夹里，再解压其中的一个文件即可得到PDF文件。 转载说明：《Learn VASP The Hard Way》转载自大师兄科研网，版权属于李强博士。 另说明：本文件下载自大师兄科研QQ群（217821116），感谢网友的整理和上传。
参考原文：1. http://www.52souji.net/how-to-import-the-model-built-using-ms-into-lammps.html 2. http://www.52souji.net/vasp-poscar2lammps-format-convert.html 3. https://sites.google.com/site/pmitev/academic/vasp-poscar2lammps-awk 需要注意的是： Dependencies. The script needs gawk in order to run. It expects that that awk is located under /bin. If this is different in your Linux distribution, just change the first line accordingly. 脚本需要gawk,需要安装此软件，指令如下： sudo apt-get install gawk 之后将VASP-poscar2lammps.awk和获得vasp格式文件放到同一文件夹下，因为awk文件在/usr/bin文件下，(我也不知道怎么到这里了)，需要将VASP-poscar2lammps.awk，第一行进行修改： 修改如下： 在终端输入： VASP-poscar2lammps.awk Cu2.vaspdata.lammps2 root@daishoujun-Virtual-Machine:/usr/bin# VASP-poscar2lammps.awk Cu2.vaspdata.lammps2 下行是对上行的解释， VASP-poscar2lammps.awkf.POSCARf.lammps 其中，VASP-poscar2lammps.awk为转换脚本名称，f.POSCAR为具有VASP POSCAR结构的文件的文件名，f.lammps为转换的LAMMPS结构文件的文件名。 其它：把脚本VASP-poscar2lammps.awk下载下来，拷贝到你的当前目录（或者添加到系统路径中），修改权限为可执行（chmod +x），+前面有个空格。
转自： http://www.v-suan.com/archives/690.html DOS，就是态密度，也就是每个轨道的电子云分布比例，通过态密度可以了解电子结构。能带跟DOS是对应的，通过能带可以得到一些重要的电子结构信息。 一、 用 VASP 计算态密度 使用的软件：VASP 辅助分析计算的小程序：read_dos.f90 操作流程及主要步骤： 主要分两步：一、静态自洽计算；二、非自洽计算（以fcc结构的Al为例） 第 一 步 静态自洽计算 （得到自洽的电荷密度CHG、CHGCAR和E-fermi，提供给下一步非自洽计算用） 准备INCAR：即定义PREC，EDIFF，ENCUT，ISTART= 0（默认），IBRION = 2，ISMEAR(对于半导体和绝缘体：ISMEAR= 0，SIGMA = 0.05；ISMEAR= 0 表示采用Gaussian smearing 方法； 对于金属：ISMEAR= 1，SIGMA = 0.2;ISMEAR= 1 表示采用1阶smearing方法) SYSTEM=Al-fcc PREC=Accurate EDIFF= 1e-5 ENCUT= 250.0 ISTART= 0 ；IBRION = 2 #【默认NSW=0，此处定义IBRION=2无意义】 ISMEAR= 1 ; SIGMA = 0.2 准备好POSCAR文件，或以优化的晶格参数作为基础，把结构优化产生的CONTCAR拷贝成POSCAR Al 3.975 0.0 0.5 0.5 0.5 0.0 0.5 0.5 0.5 0.0 1 Direct 0.0 0.0 0.0 POTCAR直接从赝势库中得到 KPOINTS Automatic generation 0 MohkorstPack 9 9 9 0.0 0.0 0.0 提交运行 记下OUTCAR中E-fermi的数值，下一步处理结果时会用到。 grep “ E-fermi ” OUTCAR 第 二 步 非自洽计算 准备INCAR：即定义ENCUT，PREC，ISTART= 1（存在WAVECAR文件,所以取1），ICHARG = 11 （ 表示从CHGCAR中读入电荷分布，并且在计算中保持不变 ） ，ISMEAR= -5 ( -5表示带Bloch修正的四面体方法；采用 tetrahedron method with Bl¨ochl corrections计算DOS更准确一些，在PWscf中也是如此) ，RWIGS (或LORBIT = 10，这时可不设RWIGS) SYSTEM = Al-fcc PREC = Accurate ENCUT = 250 #【截断能不能大于上一步非自洽计算？最好与上一步计算保持一致】 ISTART = 1；ICHARG = 11 ISMEAR = -5 LORBIT = 10 # 【对于PAW势,可设置LORBIT = 10，此时可不用设置RWIGS参数；具体解释见附】 #【一般 上述参数用来做非自洽计算、产生DOS已经足够；如NSW表示离子运动步数，默认为0(此时，IBRION默认为-1，表示 原子位置不移动)，即做静态计算。自洽还是非自洽，如何判断？ 】 准备好KPOINTS文件， 增加 k点网格(增加k点可加大积分计算精度，因为是非自洽计算，k点数虽增大但仍可以很快) Automatic generation 0 MohkorstPack 19 19 19 0.0 0.0 0.0 POSCAR、POTCAR同上一步自洽计算。 将上一步自洽计算得到的CHG、CHGCAR拷贝至同一目录下，提交运行。【 ICHARG=11，ISMEAR=-5 】 计算完后，得到包含了态密度值的DOSCAR文件，采用read_dos.f90小程序对态密度文件DOSCAR进行处理，得到总态密度dos_total.dat， 及分波态密度 *_dos.dat 。 运行read_dos.f90时会提示：please input the Ef，此时需要把第二步自洽计算得到的E-fermi的值给出。接着会提示：please input the number of atoms in the POSCAR file，此时需要把POSCAR中的原子数给出。【操作./ read_dos.x DOSCAR】 输出文件的能量值是以费米能级作为能量参考零点。dos_total.dat的第一列数据是能量值，单位为eV；第二列数据是总态密度的值，单位State/eV.unit cell；第三列数据是总态密度的积分值，也就是电子数，单位为electrons (用origin等软件绘图时把dos_total.dat的最后一列删除即可) 。 *_dos.dat是相应原子【unit cell中含有的每个原子都有各自对应的s、p、d分波态密度值】的分波态密度值，其中的第一列数据是能量值，单位为eV；第二、三、四列 数据分别对应于 s 、 p 、 d 态的分波态密度值 ，单位为State/eV.atom。 #【如果将s、p、d态密度按每种元素求和，则态密度和的单位是State/eV.atom?如何将单位转换为State/eV.el.】 #【为什么自洽完了，要做非自洽计算的英文解释见后面附中的英文解释。提问：非自洽 产生的DOSCAR 文件与上一步自洽计算产生的DOSCAR之间有什么区别？ 】 二、 用 VASP 计算 能带 使用的软件：VASP 辅助分析计算的小程序 ： read_band.f 、 get_kpoints.sh 操作流程及主要步骤： 计算材料的能带结构即色散曲线E（k），主要分两步：一 、静态自洽计算(fcc)；二、非自洽计算（以fcc结构的Al为例） 第 一 步 静态自洽计算 （同上） 【因此，可单独建一个scf的文件夹，用以产生并存储电荷密度CHG、CHGCAR和E-fermi】 第 二 步 非自洽计算 （ 在固定电子密度的情况下，得到选取 K 点的能量本征值。 ） 准备INCAR：即定义ENCUT，PREC，ISTART= 1，ICHARG = 11 （ 表示从CHGCAR中读入电荷分布，并且在计算中保持不变 ） ， 并增加NBANDS （默认值为NELECT/2+NIONS/2，NELECT和NIONS分别为电子数和离子数，可以上一步自洽计算产生的OUTCAR文件中找到这两个参数，如grep “NIONS” scf/OUTCAR，和 grep “NELECT” scf/OUTCAR。)【NBANDS的设置，详见后面的英文注释】 SYSTEM = Al-fcc PREC = Accurate EDIFF = 1e-5 ENCUT = 250 IBRION= 2 ISTART = 1；ICHARG = 11 ISMEAR= 1；SIGMA =0.2 NBANDS = 12 准备KPOINTS文件， 使用 Line-mode ，给出高对称性 k 点之间 的分割点数。【分割越密，则路径积分越准确，VASP 5.x格式有变吗？】 k-points along high symmetry lines !注释行，无特别的意义 20 ! intersections，沿G-X特殊点之间产生10个k点 Line-mode ! 程序自动产生特殊k点间的k点 Rec ! 各k点相对于倒格子基失来写的 Al的元胞的高对称性点 0.5 0.0 0.5 ! X 0.0 0.0 0.0 ! G 0.0 0.0 0.0 ! G 0.5 0.5 0.5 ! L 0.5 0.5 0.5 ! L 0.5 0.25 0.75 ! W 0.5 0.25 0.75 ! W 0.375 0.375 0.75 ! K 0.375 0.375 0.75 ! K 0.0 0.0 0.0 ! G #【对高对称性点的理解及注释见后附注部分。 】 说明：通过指定Line-mode，VASP会自动在起点和终点之间插入指定的K点数，比如上面的文件就是指定VASP计算沿着X到Gamma点再到L、W、K回到Gamma点，每个方向上各取20个K点。 POSCAR、POTCAR同上一步自洽计算。 将上一步自洽计算得到的CHG、CHGCAR拷贝至同一目录下，提交运行，计算完后可以从OUTCAR文件或者EIGENVAL文件里得到需要的每个K点的能级信息 【与运用 WAVECAR 得到的能带结构一样， CHGCAR 是 WAVECAR 的积分。是的，波函数代表一切】 运行小程序get_kpoints.sh，得到kpoints.dat（即高对称点的位置；）。 编译程 序read_band.f，运行，系统提示： please input thefermi energy:E-fermi from scf OUTCAR 将上一步得到的E-fermi值输入到屏幕，则输出bands.dat，并产生spoint.dat文件 (该文件给出了作图时高对称性点的横坐标位置)，因此，可用origin等绘图软件作图了。【read_band.f程序的解释见附】 lpf的经验和心得： 对于第一步自洽计算： 由于在能带计算时 k 点是一些在倒空间高对称线上的点 ，不能进行自洽计算【解释得好啊】，因此在进行能带计算前必须加一次自洽计算以得到精确的电荷密度值。 自洽计算得到的电荷密度文件CHGCAR是能带和DOS计算需要的输入文件。 对于非自洽计算：计算能带时ICHARG= 11，金属用ISMEAR=1；半导体或绝缘体，用ISMEAR=0 。 计算态密度：ICHARG = 11，ISMEAR=-5 【WAVECAR 文件在计算 DOS 和 BAND 过程中有用么？是万能的，但要开计算band和dos的源程序调用谁了】 用origin绘图技巧，用Vi打开spoint.dat文件，得到作图时高对称点的横坐标位置，将横坐标端点分别设为第一个、最后一个高对称点对应的坐标；画竖线，双击，修改Objectproperties—Cordinates:Units选为Scale，设置开始x和终点x均为spoint.dat文件中的高对称点数值；y值即设为纵坐标范围。【实际上图形中的横坐标范围已经定死了，无法修改；其实也恰好与spoint.dat中点的范围一致。如何断定spoint.dat文件中数值具体对应哪个高对称点符号？】 1：赝势 赝势文件夹各个目录对应起来分别是pot == PP, LDA ; pot_GGA == PP, GGA ; potpaw == PAW, LDA ; potpaw_GGA == PAW, GGA, PW91 ; potpaw_PBE ==PAW , GGA, PBE。 选择某个目录进去，我们还会发现对应每种元素往往还会有多种赝势存在。这是因为根据对截断能量的选取不同还可以分为Ga,Ga_s,Ga_h，或者根据半芯态的不同还可以分为Ga,Ga_sv,Ga_pv的不同。 2：non-scf ISTART=1;ICHARG=11时，PREC必须与做静态计算的设置一致（尤其是截断能），否则会报错 WARNING:dimensions on CHGCAR file are different ERROR: chargedensity could not be read from file CHGCAR for ICHARG10 3：POSCAR POSCAR中的那个缩放系数对于正格矢起作用如果原子坐标是用Car格式，它也起作用。如果原子坐标是Direct格式，则它不起做作用。【yexq评：缩放系数即对应PWscf中的alat，只影响晶格基矢量的大小，个人感觉不会影响原子的crystal位置】。 4：INCAR NSW是优化原子位置和晶胞个开启项，默认=0。如果没有设置NSW，即使ISIF=3，依然只做静态计算。 5：KPOINTS 采用Line-mode时，”0 0 0 ！gama“ 的最后一个0和！之间需要有空格 。【yexq评：真是经验总结，很仔细】 附录-yexq总结： 1. 为什么自洽计算完了，还要做非自洽计算，请见下面的英文解释 Calculating a DOS can be done in two ways: The simple one is to perform a static (NSW=0, IBRION=-1) selfconsistent calculation and to take the DOSCAR file from this calculation. However, the simple approach discussed above is not applicable in all cases: A high quality DOS requires usually very fine kmeshes. You should think at least in orders of 16x16x16 meshes for small cells and even for large cells you might need something like 6x6x6- or 8x8x8-meshes. The usual way, to do DOS or band-structure calculations in this case is the following: the charge density and the effective potential converge rapidly with increasing number of k-points. So, as a first step one generates a high quality charge density using a few k-points in a static selfconsistent run. The next step is to perform a non-selfconsistent calculation using the CHGCAR file from this selfconsistent run (i.e. ICHARG is set to 11, see section 6.14) . Mind, this is the only way to calculate the band structure, because for band-structure calculations the supplied k-points form usually no regular three-dimensional grid and therefore a selfconsistent calculation gives pure nonsense ! 【但是还是需要前面的自洽计算】 For ICHARG=11, all k-points can be treated independently, there is no coupling between them, because the charge density and the potential are kept fixed. Therefore there is also no need to treat the k-points within one single run simultaneously. 【可新建1个pdos文件夹，非自洽计算可安排在该文件夹中进行 】 2. NBANDS的设置问题 可见NELECT是所采用赝势中的价电子数吗？例如 对 TiH2 ：取 14/2 + 3/2 = 8.5，PdH： 取17/2+2/2=9.5 3. 对高对称性点选取的理解 采用与自洽计算所用的相同大小的胞(元胞或单胞)，来选取高对称性点； 高对称性点并不一定形成闭合的曲线； 应对结构施加正确的symmetry后，选高对称性点。这也说明，能带计算不是在uniform的格子上进行，而是沿高对称性点划定的路径积分。所以，如果结构没施加对称性，则选取的高对称性点将与施加了对称性操作的“相同”结构不一样；同理按元胞找的对称性点与按单胞找的对称性点将不一样。 对于具有相同对称性的结构，其高对称性点都有一致的约定的取法，即与具体物质无关。例如同为fcc结构的不同物质，其高对称性点将是一样的。因此可以通过下述网站输入结构的空间群号寻找高对称性点： http://www.cryst.ehu.es/cgi-bin/cryst/programs/nph-kv-list，或采用Material studio 得到高对称性点。 4.能带计算完成后，绘制能带图时，如何确定图中高对称性点的位置。 spoint.dat文件中数值按大小排列，其对应的高对称性点的顺序可能与计算前KPOINTS中输入的高对称点的顺序一致， 实际上高对称性点的crystal坐标(即KPOINTS中的三个数值)可在EIGENVAL文件中找到 (或OUTCAR中找到)。 其实，read_band.f程序很简单，即将路径上的各k点crystal坐标求均方根，然后逐点累加，并记路每一k点对应的engenvalue；当遇到特殊k点(高对称性点)时，同时将该点坐标输出到spoint.dat文件中。 该程序跟据KPOINTS文件中特殊k点之间的取样间隔来判断某点是否是特殊k点，如是则同时输出其坐标到到spoint.dat文件中。因此，当特殊k点之间插入的一般k点数变化时，则应修改源文件中的if(mod((i-1),20).eq.0 .or.i.eq.nks语句，例如将20修改为所插入的一般k点数或实际分割的k点数。因此，该程序缺乏通用性，然而程序虽简单，但却高度体现了作者对能带结构的理解。 实际高对称性点的位置可以直接求均方根确定吗？(不能，没有加和，还应该加上路径上其他点坐标产生的加和) 2.部分参数的英文注释 （1）ISMEAR = -5 | -4 | -3 | -2 | 0 | N SIGMA = width of the smearing in eV Default ISMEAR = 1 SIGMA = 0.2 ISMEAR determines how the partial occupancies fnk are set for each wavefunction. For the finite temperature LDA SIGMA determines the width of the smearing in eV. ISMEAR: -5 tetrahedron method with Bl¨ochl corrections For the calculation of the total energy in bulk materials we recommend the tetrahedron method with Bl¨ochl corrections (ISMEAR=-5). This method also gives a good account for the electronic density of states (DOS). The only drawback is that the methods is not variational with respect to the partial occupancies. Therefore the calculated forces and the stress tensor can be wrong by up to 5 to 10 % for metals. For the calculation of phonon frequencies based on forces we recommend the method of Methfessel-Paxton (ISMEAR0). For semiconductors and insulators the forces are correct, because partial occupancies do not vary and are zero or one. For the calculations of the DOS and very accurate total energy calculations (no relaxation in metals) use the tetrahedron method (ISMEAR=-5). For semiconductors or insulators use the tetrahedron method (ISMEAR=-5), if the cell is too large (or if you use only a single or two k-points) use ISMEAR=0 in combination with a small SIGMA=0.05. For relaxations in metals always use ISMEAR=1 or ISMEAR=2 and an appropriated SIGMA value (the entropy term should be less than 1 meV per atom). Mind: Avoid to use ISMEAR0 for semiconductors and insulators, since it might cause problems. For metals a sensible value is usually SIGMA= 0.2 (which is the default). （2）RWIGS or LORBIT If RWIGS or LORBIT (Wigner Seitz radii, see section 6.326.33) is set in the INCAR file, a lm- and site-projected DOS is calculated and also written to the file DOSCAR . One set of data is written for each ion, each set of data holds NDOS lines with the following data energy s-DOS p-DOS d-DOS or energy s-DOS(up) s-DOS(down) p-DOS(up) p-DOS(dwn) d-DOS(up) d-DOS(dwn) for the non spin-polarized and spin polarized case respectively. As before the written densities are understood as the difference of the integrated DOS between two pins. (3) DOSCAR file Mind: For relaxations, the DOSCAR is usually useless. If you want to get an accurate DOS for the final configuration, first copy CONTCAR to POSCAR and continue with one static (ISTART=1; NSW=0) calculation.【静态自洽计算？】 （4）LORBIT The default for LORBIT is .FALSE. (respectively 0). This flag determines, together with an appropriate RWIGS (see section 6.32), whether the PROCAR or PROOUT files (see section 5.21) are written. The file PROCAR contains the spd- and site projected wave function character of each band. The wave function character is calculated, either by projecting the wavefunctions onto spherical harmonics that are non zero within spheres of a radius RWIGS around each ion (LORIT=1, 2), or using a quick projection scheme relying that works only for the PAW method (LORBIT=10,11,12, see below). If the LORBIT flag is not equal zero, the site and l-projected density of states is also calculated. If the projector augmented wave method is used, LORBIT can also be set to 10, 11 or 12. This alternative setting selects a quick method for the determination of the spd- and site projected wave function character and does not require the specification of a Wigner-Seitz radius in the INCAR file (the RWIGS line is neglected in this case). The method works only for PAW POTCAR files and not for ultrasoft or norm conserving pseudopotentials. 即对每个能带的波函数进行s pd和site分解或投影；LORBIT = 10，不设RWIGS，针对PAW势的计算。 LORBIT = 10表示对波函数按原子site进行s 、p、d轨道分解或投影；并输出投影波函数到PROCAR or PROOUT 文件中，详见上面英文注释部分：This alternative setting selects a quick method for the determination of the spd- and site projected wave function character。什么样的quick method?
1, Restart the calculation with a pre-converged PBE WAVECAR (ISTART = 1, the default value when WAVECAR exists) 2, Decrease the mixing parameter AMIX from the default 0.4 to 0.2 and to 0.1, increase NELM to 12 0 or more. 3, Try the Anderson mixing by setting IMIX = 1, AMIX = 0.1, BMIX = 0.0001, AMIX_MAG = 2*AMIX, BMIX_MAG = 0.0001 4, The above methods are usually enough, and I would pray for a while if they are not ... 来自： https://sites.google.com/site/haoweipeng/softwares https://sites.google.com/site/haoweipeng/softwares
现有版本vasp只支持wannier1.2，但实际上vasp是可以和wannier2.0兼容的 具体操作如下： The problem seems to be due to the new version of wannier90-2.0.0 (compiled with the current Intel compiler) The variables list of wannier_setup (src/wannier_lib.F90) has two additions compared to wannier90-1.2: proj_s_loc proj_s_qaxis_loc which are optional, intent(out). VASP does not (yet) use these quantities, and as such, calls wannier_setup without them. This should not be a problem, since these variables are declared as optionals, but apparently the ifort compiler messes up (we have not tried other compilers). Solution: 1) delete these 2 parameters from the call list of wannier_setup in wannier_lib.F90 2) comment out all lines in wannier_setup referring to these variables,
《Learn VASP The Hard Way》由李强博士撰写，目前基础篇（ex0 ~ ex27）已经完结，感谢作者付出！ 有问题可以直接 到大师兄的QQ群（ 遇见大师兄 217821116 ）讨论， 或关注大师兄的 微信公众号 （ BigBroScience ）。 大师兄这套教程的编写方式受《learn python the hard way》（笨方法学Python）这本书的影响, (参考: https://www.douban.com/group/topic/73878230/)， learn python the hard way.pdf 这种教学模式确实非常好，特别对于没有基础或者基础不好的人来说，用该方法来学一门新技术有非常显著的效果，可以快速入门。
Ex24 ：乙醇分子的频率计算（二） 频率计算的输出与 POSCAR 原子的固定 按照前面一节介绍的方法，结构优化过程完毕后，准备频率计算的输入文件，提交任务等待结束。 1 查看结果 1.1 查看 OSZICAR ，你会发现一共计算了 55 步。 下面我们分析一下这 55 步是怎么回事： A）乙醇分子CH3CH2OH 含有9 个原子，每个原子在 xyz三个方向上均有一个自由度，共9*3 = 27个 B）我们设置的NFREE = 2 ，也就是在每个方向上 +POTIM 和 –POTIM都移动并算一下，这样就有了27*2 = 54 步：官网原文如下：大家自己去查阅IBRION和NFREE的相关内容。 The parameter NFREEdetermines how many displacements are used for each direction and ion, andPOTIM determines the step size. The step size is defaulted to 0.015 Å (startingfrom VASP.5.1), if too large values are supplied in the input file. Expertiseshows that this is a very reasonable compromise. NFREE=2 usescentral differences, i.e., each ion is displaced by a small positive andnegative displacement, ±POTIM, along each of the cartesian directions. C）还有一步：55-54 = 1, 这一步指的是第一个离子步，为频率计算前的单点计算。 所以当你设置了NFREE=2的时候，频率计算需要1+N*6步。N为体系中的振动的原子数。 1.2 我们这么算，默认的是所有的原子在 xyz 三个方向上均可移动。但很多时候，我们只需要振动感兴趣的原子或者某一特定的方向就可以了，也就是说选择性计算频率。比如我们只想算乙醇的羟基振动峰，那么其他的原子就可以固定了。 VASP 中怎么实现这个功能呢？首先我们看一下当前的 POSCAR ： 左下角的 : set nu 显示文本的行数，取消行数可以通过 :set nonu 如果想要选择性地固定某些原子，我们需要以下几个步骤： 1） 在第七行和第8行之间插入一行，内容为 Selective Dynamics ，前面讲过VASP只认第一个字母，也就是S是必须的，Selective Dynamics 和 See, Sea….等其他S开头的都是一个效果的。 2） 加入Selective之后，我们需要在每一行的坐标后面加上 T或者F表示允许和禁止移动。 这里我们需要加三个T或者F，表示在x y z 三个方向上选择性固定原子的移动， 三个方向都允许：T T T 三个方向上都不允许为 F F F x移动，y和z方向上固定为： T F F x和y方向上固定，z方向振动： F F T 以此类推，其他的大家根据自己的情况固定。 ------------------------- 师兄，那么我们要做的就是在POSCAR中逐行写上 T T T 或者 F F F 就行了吧 ? 是的，下面大师兄教给你的是通过 vim 的命令或者p4vasp实现这个功能。 2 通过 Vim 实现原子的固定和选择 2.1 加入 Select 的关键字母 S ，并在坐标后面全部加上 T T T 图中1：插入一行，告诉VASP我们要选择性的固定某些原子或者在某些方向上； 图中2： 10,18s 中 s 代表选择（select）的意思，这里表示我们选中了第10到18行，10和18之间有个逗号表示连续；10,18s后面用一个 / 分开，紧跟着你我们要替换的内容； $ 在这里是末尾的意思，$/T T T我们要把每一行的最后替换成 T T T 后面再用一个/分开，加上g 表示 global 全部替换的意思。 输入完毕后，回车，效果如下： 每一行的末尾都加入了 T T T 箭头指的地方告诉我们： 9 行中的 9 个地方发生了替换。 （自己复习下sed的用法，比较两者的区别） 2.2 下面我们要把 OH 之外的原子全部固定住： 通过 p4vasp 我们可以知道， OH 的两个原子为第一个（ O ）和最后一个（ H ），因此我们把第 11 行到 17 行中的所有 T 替换成 F 就可以了。我们可以使用 : 11,17s/T/F/g 来实现。最后效果如下： C ）扩展： 如果所有的原子后面为： T T T ，这与我们之前的计算效果是一样的（第 8 行没有 S ，坐标后面没有 T T T ）。此外， 固定之前要先找到哪些原子我们希望固定的，以及它们在坐标中的顺序。 3 使用 p4vasp 实现上述功能 基于 p4vasp 的可视化特性，直接用鼠标操作是很多 linux 小白最喜欢看到的，下面我们主要讲解一下p4vasp的操作，请务必将vim和p4vasp的操作关联起来，这样你会就会发现可视化和命令之间的微妙关系了。 3.1 打开 p4vasp ，导入 POSCAR ； System 显示 ethanol说明数据已导入。 3.2 鼠标系列操作 A ）点一下 Show 按钮，显示乙醇分子结构； B） 点击 Build 按钮：显示分子的坐标信息：这里的坐标顺序和 POSCAR 完全一致： C） 选择乙醇的羟基（使用空格键选择原子），下图中可以知道 O 和 H 在坐标顺序中为第一个和最后一个； D ）选中右上角的 Selective Dynamics ，效果如图中箭头 2 指出来的部分，这个动作和我们之前在第八行插入 Selective Dynamics 以及在坐标每一行加入 T T T 是等效的 E）选中第二行坐标，摁住 shift 键，然后再点一下倒数第二行，选中整个区域后，点击右上方的 Unselect ， 这个动作等效于我们将 11, 到 17 行中的 T 全部替换为 F 。 F ）操作完成后，坐上角： File à Save System As 保存新的 POSCAR 即可。 4 扩展练习： 4.1 掌握本节的两个 POSCAR 处理方法； 4.2 学会计算频率分析所需的步数； 4.3 浏览 OUTCAR ，查找频率输出结果 5 总结： 通过对比手动敲命令修改 POSCAR 和使用 p4vasp 进行鼠标操作，这里大师兄希望大家能掌握以下 3 点： 1） 学会 Vim 的使用技巧，当然了， Vim 及其强大，完全掌握基本不可能，但最基本的操作要了解； 2） P4vasp 的操作要熟悉，不知道怎么操作的，导入一个计算文件，随便点点，找找感觉； 3） 在使用 p4vasp 的操作中，你要学会思考：怎么将鼠标操作转化为命令语言来实现，为将来写脚本做好准备，毕竟很多时候，我们用不了可视化界面，就只能手动修改格式了,结合命令，比可视化操作更快，还可以批量进行。 ex24-乙醇分子的频率计算（二）.pdf 转载说明：《Learn VASP The Hard Way》转载自大师兄科研网(http://www.bigbrosci.cn)，版权属于李强 博士 （QQ: 122103465）
Ex-DOS1 DOS计算（一） 使用VASP计算，很多时候都逃不掉DOS，能带计算的相关问题，尤其是对于计算材料的童鞋们，更是家常便饭一般。群里很多人，很多新手们都时常在讨论DOS的计算。这里我们通过VASP官网的说明，解释一下算DOS的具体步骤。注意这一节没有指定练习的章节，因为这不属于菜鸟篇的内容，大师兄临时拿过来先打消大家的疑虑，等写到这一节的时候，再给分配章节数目。 1 KPOINTS： 1.1 K点数目 与结构优化相比，算DOS的时候，需要用到更多的K点数目，这是因为K点越多，画出来的DOS图质量越高。 引用官网的话： A high quality DOS requires usually very fine k-meshes. 1.2 K点数目的选取 K点数目越多越好，我们该如何设置K点数目呢？ 还记的前面我们讲到的K点选择的经验规则吗？那一个规则可以认为是我们平时计算时K点选择的标配。对于DOS计算，我们就需要把配置提高一个档次了。 一般来说，K * a = 45 左右之间完全可以满足你的要求，大伙可以根据这个经验来选择K点。 2 ISMEAR的选择： 2.1 看官网的话： For the calculation of the total energy in bulk materials we recommend thetetrahedron method with Blöchl corrections (ISMEAR=-5). This method also givesa smooth nice electronic density of states (DOS). 也就是说 ISMEAR = -5 的时候（Blöchl修正的四面体方法），我们可以得到一个非常平滑的DOS图。很多人问怎么才能获得平滑的DOS图，这就是答案！！ISMEAR = -5 并且使用较多的K点。 2.2 注意部分 1）不能使用ISMEAR= -5 的情况： 如果你的模型很大，只用了一个gamma点，或者K点的数目小于等于3 的时候，使用ISMEAR = -5 会导致计算出错。（the tetrahedron method is not applicable, if less than three k-points are used.）很多人只用了gamma点，然后用修正的四面体方法计算能带，最后到群里咨询出错的原因！！！ 2）对于金属体系来说，算能带，DOS的时候，结构保持不动（ 这是关键 ），可以放心使用ISMEAR= -5。但是，结构优化的时候不能使用ISMEAR= -5（注意：是 优化结构的时候不能用 ！），这是因为四面体方法不能很好地处理费米能级处的电子占据情况，导致算出来的力会有一定百分比的误差。参考： https://cms.mpi.univie.ac.at/vasp/vasp/Partial_occupancies_different_methods.html 3）半导体和绝缘体的体系可以尽情使用ISMEAR= -5，但绝对不能在这两个体系中使用 ISMEAR0。等于0 则可以。 4）使用ISMEAR= -5 的时候，SIGMA的取值没有影响。 5）总结：算DOS，只要K点不少于3，其他情况都用ISMEAR=-5。 ２.3 不能使用ISMEAR = -5 的情况怎么办？ 两个解决办法： 1） 既然Ｋ点不够，那么我就增加Ｋ点，然后再使用ISMEAR= -5 （ 简单粗暴，强烈推荐使用 ）； 2）服务器不给力，不能增加Ｋ点的时候，怎么办？又分两种情况：（少许有些绕弯） 2.1）半导体和绝缘体，使用ISMEAR= 0 （Gaussian Smearing，高斯展宽，后面统一简称GS方法） ；（再次注意：此时绝对不能大于0！！！） 2.2）金属，可以使用ISMEAR = 0，也可以使用Methfessel-Paxton (MP)方法：后面统一简称：MP方法。 ISMEAR = 1, 2….Ｎ，一般来说，ISMEAR =０和 1 基本就可以了。（注意：金属可以等于0，也可大于0） 2.4 SIGMA 既然选择了ISMEAR，就逃不开SIGMA的取值。 MP方法（ISMEAR=1..N）： SIGMA取值太大，计算出来的能量可能不正确；SIGMA取值越小，计算越精确，需要的时间也就越多。http://cms.mpi.univie.ac.at/vasp/guide/node124.html SIGMA的取值和KPOINTS密切相关，Kpoints确定之后，使用多大的SIGMA值，大家最好测试一下。原则如下： SIGMA取值在 保证 OUTCAR中' entropy T*S' 这项的能量平均到每个原子上小于 1 meV的 前提下 ： 尽可能地大 。这样做可以保证准确度的同时加快收敛速度。 记住： VASP学习最快的途径就是不停地看官网，然后亲自上手去测试，测试，测试！并观察分析结果。 高斯展宽,Gaussian Smearing (ISMEAR = 0)： 1) 对于大部分的体系都能得到理想的结果； 2) SIGMA取值比较大的时候会得到与MP方法相近的误差；但是误差多大，GS方法不可以得到，而MP方法可以。从这一点上来说，MP要比GS好些； 3) 从经验上来说：对于金属体系，使用MP方法（ISMEAR=1..N）时，SIGMA= 0.10 足够了，官网给的参考值是0.20。 4） 使用GS方法的时候（ISMEAR=0），SIGMA的数值要测试下，保证'entropy T*S'这一项平均到每个原子上小于0.001 eV也就是1meV。5）不想测试，对于金属体系：ISMEAR=0.05是一个很安全的选择。 对于半导体和绝缘体，SIGMA取值要小，SIGMA = 0.01 – 0.05 之间也是很安全的。 3 测试： 3.1 体系选取 以下是本人的一个例子： 84个Ni原子，17个C，20个H，6个O，共127个原子。（只是一个例子，不具有代表性，大家对于自己的体系，要亲自测试一番。） KPOINTS的选择如下： GS （ISMEAR=0）和 MP（ISMEAR=1）方法均选择了SIGMA=0.03 和0.10 两个参数。四个测试文件，命名为： gs-0.03，gs-0.10, mp-0.03,和 mp-0.10 3.2 结果分析： 图中列出来的是单点计算的能量。 3.2.1 GS方法分析： GS中 0.03 和 0.10 比较： 1） Without entropy 后面的能量差别很大， sigma -0 后面的能量差别很小；这也是大师兄为什么推荐使用sigma 后面能量的原因。 2）使用ISMEAR = 0.03，最后电子步的结果： entropy T*S 后面的能量为：-0.078, 平均到每个原子上为： -0.078/127 = 0.000614 eV/atom = 0.614 meV/atom 1 meV/atom。 说明取值没啥问题。 使用ISMEAR= 0.10 的结果： -0.9234/127 = 0.00727 eV/atom =7.27 meV/atom 1meV/atom。 此时，我们的取值就有问题了。说明是不可靠的，需要用一个更小的sigma值。 如果使用ISMEAR=0，会有两个能量，without entropy 和 sigma – 0。 SIGMA越小，二者之间的的差别越小！一般来说，不同的SIGMA值，Sigma – 0 的结果差别不大，而withoutentropy后面的能量则差别有时候会很大。在你的计算中，只能用一个能量，不可以混用。建议使用Sigma – 0后面的能量，且体系中所有的计算都用这个。如果使用了without entropy 后面的能量，那么所有的计算都用这个能量。 对于A和B之间的能量差：ΔE = EA-EB，使用without entropy 后面的能量计算出来的ΔE 和使用 sigma – 0 的出来的ΔE，结果是一样的，差别很小。 3.2.2 MP方法的结果分析： SIGMA = 0.03 的结果 SIGMA = 0.10 的结果 分析方式如GS的测试。 使用小的SIGMA值，同样可以得到更小的 entropy T*S的数值。也就是计算的越精确。 without entropy和sigma-0 后的能量差也越小。注意，当ISMEAR = 0.10 的时候，entropy T*S 值为0.0168 eV，比之前高斯展宽时ISMEAR=0.10的数值要小很多。这表明对于当前测试的体系，SIGMA= 0.10 的时候，使用MP比GS方法的效果要好。当SIGMA数值很小的时候（ SIGMA= 0.03 ），两者不相上下。 3.2.3小结一下： SIMGA小的时候，ISMEAR = 0 和 ISMEAR = 1 效果基本一样； SIMGA数值大的时候，ISMEAR= 1 效果更好。 所以，如果你首先设置一个很小的SIMGA数值，不用太担心ISMEAR的选取。但反过来，如果你首先设置了ISMEAR这一项，就需要考虑SIGMA的大小了。懒得考虑就设置一个比较小的SIGMA值，比如0.05，即可。 3.3 ISMEAR 和SIGMA对计算时间的影响。 官网里面说了 SIGMA的选择对计算的收敛会产生影响。我们分析下测试所花费的时间。 前面均是单点计算，计算所需的总时间用grep LOOP+ OUTCAR查询 如果你想查询一个电子步需要的时间用： grep LOOP+ OUTCAR （去掉+号）。 Grep RMM OSZICAR | tail –n 1 获取计算了多少电子步。 计算时间不是看Elapsed后面的那项吗？ 是的，Elapsed的时间是总的计算时间。对于我们这个单点计算，LOOP+后面的能量和 Elapsed 的能量差别不大。我们这里只考虑SCF收敛过程中所花费的时间。 从前面的计算分析看，平均到每一个电子步上面，使用同一个 ISMEAR 的时候，小的 SIGMA 值对应的单个电子步时间略大一些。使用 ISMEAR=0 时，两个 SIGMA 值对应的单电子步时间差别更大。 但是官网说， if σ istoo small the convergence speed with the number of k-points will deteriorate. 这里我们并没有明显观测到，反而在使用 ISMEAR=0 + SIGMA = 0.03 的时候，收敛步数减少了一些，从而略微节省了时间。 但是官网说， if σ is too small the convergence speed with the number of k-points will deteriorate. 这里我们并没有明显观测到，反而在使用 ISMEAR=0 和 SIGMA = 0.03 的时候，收敛步数减少了一些，从而略微节省了时间。原因在于 too small 是多么地小，什么程度才算 too small ，官网没有明确说清楚。为此，专门设置了一下 SIGMA = 0.001 来测试一番。 从这里可以看出来，SIGMA= 0.001的时候，计算时间明显增加了不少，为6112.84 s。而SIGMA = 0.1 的时候，同样计算了 93步，但花费的时间为：5913.63 s。平均到每一个电子步上，也有所增加。 但是，对于大伙来说，SIGMA=0.01就已经很小了，too small 的时候并不常见。所以，不想测试的话， ISMEAR的选择可以在0.1-0.01之间。 注意：通过SIGMA的测试，我们可以发现，它对计算时间的影响完全没有KPOINTS，ENCUT等的影响大。大家在选择的时候，不要太过焦虑和纠结。大师兄建议大家 最好还是要测试一下，选择一个合适的SIGMA数值。原因如下： ISMEAR = 0 的时候， SIGMA = 0.03 和 0.10 ，单个电子步的计算时间为： 64.337 和 64.222s。 假设我们优化一个结构花了200个离子步，每个离子步有20个电子步，总共需要 200*20 = 40000个电子步。与0.03相比，使用0.10节省的时间为： 40000*（64.337-64.222）/3600= 1.3 h。 也就是一个任务省出来 1.3小时的时间。当你有10个任务的时候，就能省出来13个小时的时间..... 4 扩展练习： 4.1 阅读VASP官网关于ISMEAR和SIGMA的所有说明： 4.2 下载VASP的pdf说明书，搜索书中所有的ISMEAR和SIGMA关键词，阅读所有相关的内容； 4.3 思考SMEAR方法的意义？SIGMA的意义？ 4.4 查看VASP说明书，查阅相关文献，了解MP和GS方法 4.5 分析下为什么算DOS的时候，要算两步： selfconsistent 和 none- selfconsistent calculations? 5 总结： 看完本节：你应该知道计算DOS的时候，KPOINTS设置的要大一些。 ISMEAR要用-5。 Kpoints因计算硬件限制不能设置的很大，数目小于3的时候，对于金属，非金属体系均可以使用ISMEAR=0，SIGMA的数值需要测试一下，一般来说在0.01-0.05之间足够了。 此外，金属体系还可以用ISMEAR=1..N，官网建议SIGMA为0.20，太小的SIGMA值对收敛会产生影响。使用0.01-0.10的数值都是很安全的选择。具体对收敛时间的影响大家要测试一下。本节的测试例子中，只有 SIGMA=0.001的时候少许增加了计算时间，且SIGMA的影响不如KPOINTS和ENCUT等的影响作用大。 非DOS计算的时候，对于金属来说ISMEAR不能等于 -5，优先使用ISMEAR= 1。非金属来说（半导体和绝缘体），不能 0 。对于所有的体系， ISMEAR= 0 则是一个很安全的选择，但SIGMA的数值要测试一下。说了这么多废话，还是官网简单明了： For further considerations on the choice for the smearing method see sections 9.4,10.6. To summarize, use the following guidelines: For semiconductors or insulators use always tetrahedron method (ISMEAR=-5), if the cell is too large (or if you use only 1 or two k-points) use ISMEAR=0. For relaxations in metals always use ISMEAR=1 and an appropriated SIGMA value (the entropy term should less than 1 meV per atom ). Mind: Avoid to use ISMEAR0 for semiconductors and insulators, it might result in problems. For metals a sensible value is usually SIGMA= 0.2 (that's the value we use for most transition metal surfaces). For the DOS and very accurate total energy calculations ( no relaxation in metals ) use the tetrahedron method (ISMEAR=-5). ex-DOS1.pdf 转载说明：《Learn VASP The Hard Way》转载自大师兄科研网(http://www.bigbrosci.cn)，版权属于李强 博士 （QQ: 122103465）
Ex20 谁偷走的我的机时？（四） 模型对计算时间的影响 上一节介绍的 KPOINTS 对计算的影响，相信大家已经认真阅读参考书的第三章部分了。本季我们讨论一下模型的大小对计算的影响。主要体现在晶胞的尺寸，对称性以及对K点的影响上。 1 测试工作： 为了方便处理，我们把 O2 计算的格子设置为长宽高均为 8.0 Å 。 重复之前KPOINTS的批处理操作，我们可以获得一系列不同大小格子的文件夹。如下图： 命令： for i in $(seq 10 2 20); do cp 8 $i; sed –i 3,5s/8.0/$i/g $i/POSCAR ;done 2 测试结果分析 2.1 模型大小对计算时间的影响 注意：在后面加入 sort –n 后输出的变化。 使用数据作图： 从图中可以看出来，计算时间随着格子的大小，需要的计算时间增加的很快。 注意：在测试中， Kpoint 一直保持不变（因为只有一个Gamma点）。而在我们实际的计算操作中，使用1x1x1 Kpoints的机会并不多。如果格子在某个方向增加了 2 倍，那么对应的改方向的 K 点就需要除以 2 。重复一下上节的经验指导。也就是在计算过程中，保持 k*a 保持不变。当然， k*a 是我们提前测试好的。 举例： 一个 10x10x10 Å 3 的体相材料，我们计算的时候 K 点设置为： 6x6x6 。 当我们将材料在 x 方向增加 1 倍，变为 20x10x10 Å 3 。为保持一致的精确度，那么我们的 K 点需要设置为： 3x6x6 。 这是因为倒易晶格矢量和实际的晶格矢量之间存在着倒数的关系： 注： 类似的图，不加说明，均出自我们的参考书 也就是说，我们选取的晶格越大，倒易晶格矢量越小。用同等数目的 K 点分布到倒易晶格中，网格的密度也会越大，从而造成计算量的增加。 2.2 体系的对称性对计算速度的影响： 2.2.1 K 点保持不变： 这一点前面关于氧原子的计算就已经介绍到了，降低体系的对称性会增加额外的计算时间。如图：将 12x12x12 Å 3 （计算需要 186.86 s ）的格子修改如下： 计算结束后，查看时间，为194.9 s, 计算时间增加了 8 秒。 2.2.2 对称性对 K 点的影响： 体系的对称性不仅仅提现在前面的计算中，更可以在计算中极大地减少 K 点的数目，从而加快计算，节省时间。这一点我们引用参考书中的一段话： 2.3.3 模型对称性与K点对称性的关系 在这里，体系的对称性与 K 点对称性的匹配问题，尤其是对于 hexagonal 的结构来说，必须要使用 gamma centered points. 也就是第三行的第一个字母必须为 G 或者 g 。我们看一下官网的原话： We strongly recommend to use only Gamma centered grids forhexagonal lattices. Many tests we have performed indicate that the energy convergessignificantly faster with centered grids than with standard Monkhorst Pack grids. Grids generated with the M setting in the third line, in fact do not have full hexagonal symmerty. 如果你不确定自己的体系，直接用 G 就可以了。 For reasons of safety it might be a good choice to use only meshes with theirorigin at (switch G or g on third line or odd divisions) if the tetrahedron method is used. 3 扩展练习： 1 认真阅读： Density Functional Theory: A Practical Introduction: 第三章的前两节； 2 VASP 官网查找 K 点相关的说明。 4 总结： 学习完本节，大家应该掌握的内容有： 4.1 晶格大小对计算时间的影响； 4.2 体系的对称性对计算时间的影响； 4.3 掌握 K 点和晶格大小的经验规则； 4.4 晶格对称性和 K 点对称性的一致性。 ex20-谁偷走的我的机时？（四）.pdf 转载说明：《Learn VASP The Hard Way》转载自大师兄科研网(http://www.bigbrosci.cn)，版权属于李强 博士 （QQ: 122103465）
Ex6 OUTCAR 的基本内容 OUTCAR是VASP的主要输出文件，但官网并没有给出太多详细具体的说明。大师兄对这一点也存在疑虑，思考了下其中的原因：VASP可以计算的性质很多，每一个特定的任务都有其对应的输出内容，这就导致了OUTCAR的复杂多样性。工作人员偷懒，也就把这个问题留给了用户自己去解决，另外一个原因就是里面的条目太多了，每一项单独解释都需要花费大量的时间，这一点大师兄在认真查看OUTCAR的时候发现的。 本节，我们通过前面的计算例子，简单介绍一下OUTCAR的基本结构以及输出内容，给大家一个大体的印象，不至于看到这么多行的结果眼花缭乱。在介绍之前要告诉大家的是，VASP各个输出部分之间用横杠分割（---------------------），当你看到横杠的时候，就知道要进入结果的下一个部分内容了。 VASP会列出来其版本，时间，以及服务器的相关情况，如下图： 输出我们的INCAR和POTCAR的部分信息： 每次运行的时候，（即使你的输入是正确的），VASP都会输出一个大大的WARNING来吓唬你，大家可以忽略它。但如果你的计算失败了，这个警告信息对你排查错误可能会有所帮助。 下面是POTCAR的基本信息，如果你想通过OUTCAR查看POTCAR中的元素时，可以使用下面的命令： grep POTCAR OUTCAR grep TIT OUTCAR grep ENMAX OUTCAR ….. ZVAL是该POTCAR中对应元素的价电子，这里氧原子含有6个外层价电子。 POSCAR的一些基本信息： 坐标格式，原子位置，以及晶胞的形状大小。 体系的对称性以及点群操作相关的信息，在这里我们的体系是立方体，为O_h的点群，有48个对称操作。群论的知识大师兄早已原封不动地还给老师了，在这里就不再详细介绍了… K点信息：想查看K点个数： grep reciprocal OUTCAR 计算中参数详情（默认的参数值以及代表的意义也列出来了，这个地方大家仔细看下, 一些不常见的参数，采用默认值就好，不用在INCAR里面额外再写一遍！） 本计算的文字描述，任务类型： 体系大小，K点数目等信息: 计算所需的内存信息: 开始电子步的的迭代过程：注意能量在不同迭代步数中的变化： 每一电子步完成后，输出结果同时在OSZICAR中更新一行。 其中各项的具体含义：alpha Z and the Ewald energy define the electrostatic interaction of the ions in a compensatingelectron gas. The alpha Z componentdeals with the divergent parts (G=0). The followingparts are the Hartree and exchange correlation energy as defined in theKohn-Sham Hamiltonian. The entropy part stems from the smearing(using the free energy as variational parameter, electronic entropy), EBANDS fromKohn-Sham eigenvalues, and EATOM is the reference energy for the potential(which is defined in the POTCAR file). 出自扩展阅读3 迭代结束，输出主要的结果：费米能级以及能带信息。Band 1 对应的是2个 2s 电子， band 2 34 对应的是4个2P电子。固体物理中，费米能级对应的是最高电子占据轨道的能量，也就是HOMO，大家可以对比下band 2 3 4 和费米能级的能量。再思考下，这个结果是否正确？ 体系的坐标，各个方向力的大小，以及总的能量，不包含Entropy的能量，以及sigma趋于0 时的能量。这里 without 和 energy 之前有两个空格。 计算的内存和时间等信息，看到下面，说明计算正常结束了。 分析完前面的内容，大家会发现：具体到里面各项的含义以及各个细节上，还有很多值得讨论的地方，比如群论，薛定谔方程求解过程，POTCAR的相关信息等。对于新手来说，看完本节，能大体浏览下来，知道各个部分包含什么内容就很不错了。后面如果有什么疑问，可以直接提出来，这一部分后续慢慢完善。 扩展训练： 1 结合氧原子的外层电子排布，思考下，本例算的结果是否正确？ 2 POTCAR 的相关信息，浏览官网，尽可能多的掌握信息： https://www.vasp.at/vasp-workshop/slides/pseudoppdatabase.pdf https://cms.mpi.univie.ac.at/vasp/vasp/POTCAR_file.html http://cms.mpi.univie.ac.at/vasp/vasp/Pseudopotentials_supplied_with_VASP_package.html https://cms.mpi.univie.ac.at/vasp/vasp/PAW_potentials.html https://cms.mpi.univie.ac.at/vasp/vasp/Recommended_PAW_potentials_DFT_calculations_using_vasp_5_2.html https://cms.mpi.univie.ac.at/vasp/vasp/Recommended_GW_PAW_potentials_vasp_5_2.html 3 https://cms.mpi.univie.ac.at/vasp-forum/viewtopic.php?t=273 总结： 本练习中，带领大家粗略浏览了一遍 OUTCAR 中各部分的信息。大家在浏览的时候时刻要思考，这部分包含的什么内容，具有什么物理或者化学意义？ 怎么用 grep 关键词获取有用的信息等。 由于计算内容的多样性，对于OUTCAR的详细解释，目前来说还需要很多的时间和精力去完成补充，当然，在后续的计算过程中，我们还会结合具体的例子进行讲解。 ex6-OUTCAR的基本内容.pdf ex6.zip 转载说明：《Learn VASP The Hard Way》转载自大师兄科研网(http://www.bigbrosci.cn)，版权属于李强 博士 （QQ: 122103465）
转载说明：《Learn VASP The Hard Way》转载自大师兄科研网(http://www.bigbrosci.cn)，版权属于李强 博士 （QQ: 122103465） 此书主要针对于0-6个月的VASP初学者，VASP小白，或者刚刚转换计算方向（从VASP计算一个性质到另外一个性质时），以及某一部分计算细节生疏需要复习的科研工作者。列出来一堆求解薛定谔方程中的各个公式定理等，肯定会对初学者造成一定的误导，因此本书不讨论过多量子力学的基本原理。旨在为初学者提供一个快速进入计算而又避免过多新手错误的方法。 由于对计算细节的不了解，且无人指导（导师啪啪啪打脸），很多人在计算了一个多月后发现自己的参数设置错了，并且需要重算，浪费了大把的时间和机时。大师兄遇到过很多种类似的情况，周围的朋友遇到这种情况的也不在少数。更有甚者，课题做完了才发现是错的。 出现错误并不可怕，可怕的是我们不知道长教训，后面继续犯错。这也是本书的一个出发点，首先保证大家提交的任务准确无误，可以尽最大可能避免遇到前面类似的问题，进而起到间接节约时间的作用。在正确计算的同时，大家可以从头学习密度泛函理论，阅读相关课题的参考文献等。下面链接列举了一些相关的书籍介绍： https://www.zhihu.com/question/20329413 http://muchong.com/bbs/viewthread.php?tid=8061842fpage=2 本书的另一个出发点就是，大师兄加了很多关于计算的QQ群，但是群里面很多问题都非常低级，令人费解，或者说是匪夷所思的，从最基本的建模都做不到，到计算结果不会分析等等。这些人简单而又低级的问题充斥在各个QQ群里。暂且不说这些人的导师有多么地不负责任。很多热心的人却在群里整天忙着应付这些问题，而对于自身，除了得到个活雷锋的标签外，对理论功底的提高，帮助甚微。可以说是花自己的时间替别人指导学生。 不论群主给自己的群定位有多高，高级群，中级群，精英群等等，都避免不了这样的问题出现。本书主要通过实例引导大家主动思考去解决这些最基本的常见问题，进而避免因自己的低级问题浪费他人的时间。 大师兄希望的是，对于求助或者应助的人，大家尽可能地讨论一些更高级，更深层次的科研问题，而不是浪费在这些低级的问题上。即使在新的计算中遇到了之前没有碰到过的小细节，自己也知道怎么去动脑子，主动解决。 引用一个牛人的话：量化入门一定要按照正确、合理的顺序，循序渐进。要从简单到复杂，从构建整体的知识框架并会用最常用的量化程序算最基本的任务开始，再到逐渐了解更多理论、深化对理论的理解以及掌握更多计算技巧和程序的使用。具体可以浏览下网站：http://sobereva.com/355。 学习计算化学的人，对解决科研问题都有着一种执着的态度，通过构建模型来阐明已知或者预测位置的结果。相对于做计算的科研工作者们，虽然我们没有实验技巧的提高，但我们可以通过训练自己的大脑来弥补。懂得思考的人永远站在社会发展的最前端。 如何学习本书，大师兄在学习程序时，受到learn_python_the_hard_way这本书的启发，务实是这本书的一大特色，开始学习语言，乱七八糟的先统统给我闭嘴，照着我的代码练习一番，然后再自己思考琢磨，出现问题拿自己的代码和作者代码比较找出原因。通过系统地学习，随着水平的提高，再逐步解释前面未讲解的内容。这一种学习方法非常适合零基础的菜鸟，因为一开始太多的概念根本不可能一股脑儿全部接受。从简单入手，指导着循序渐进，最后达到精通。 对于量化计算，本书也采用这样的思路，手把手先教会大家如何计算，如何避免错误。 从最基本的计算开始，通过示例讲解，结合一些脚本的使用，引导大家思考解决自己的问题。因此，在这本书的学习过程里，每一章节会对应一个例子，大家务必手动搭建模型，输入文件（切忌复制粘贴），然后进行计算，得到和大师兄一致的结果。为了引导大家主动浏览官网解决问题，很多都会采用官网的例子，大师兄会重新计算后放到章节里面（官网的结果版本比较陈旧）。 此外，论坛或者QQ群里，有很多无知或者stupid的回复，处在迷糊之际的菜鸟由于对自己的不自信，会一股脑儿去相信别人错误的观点，进而一路错下去，这是最可怕的。所以，如果有疑问，可以先酝酿一两天，然后再大胆提出来，改正就是进步。 由于时间有限，且保证书的质量，更新可能会比较慢（计划1年），跟不上同学的科研速度，在这段时间，大家有什么问题，对之前的章节有什么意见，可以尽管提出来，这对充实完善本书有着莫大的推动作用。你的一个问题，可以解决无数个师弟师妹们的疑惑，他们在后续学习中会受益匪浅。 There are kinds of questions you will find yoursel asking and not knowing where to get quick answers from. That's what BigBro(a)s are trying to fix. 你会发现自己在问各种各样的问题，但不知道从哪里可以得到快速解答。 这正是大师兄正在尝试解决的问题：结合最基本的化学常识和软件计算细节，做最好的快速基本入门书籍。 对于刚刚入门的小白，下面是大师兄推荐的两本书和一个网站： 参考书： 1 Density functional theory：A practical introduction, by David Sholl 2 Quantum Chemistry by Levine 3 Vasp Manual: http://www.vasp.at ex-0-学习前的准备工作.pdf 说明：本系列文章已经得到作者的转载许可。
VERY BAD NEWS! internal error in subroutine SGRCON 转自： http://blog.sina.com.cn/s/blog_650463380102vgxu.html VERY BAD NEWS! internal error in subroutine INVGRP: inverse of rotation matrix was not found (increase SYMPREC)2 VERY BAD NEWS! internal error in subroutine INVGRP: inverse of rotation matrix was not found (increase SYMPREC)2 或者 Very bad news! internal error in subroutine SGRCON: Found some non-integer element in rotation matrix: 3 等错误或警告，且不能运行。 SYMPREC=1E-04 就行了，不加的话就老报错，且不能运行。 vasp-2012manual.pdf中： 6.27 ISYM-tag and SYMPREC-tag The SYMPREC-tag (VASP.4.4.4 and newer versions only) determines how accurate the positions in the POSCAR file must be. The default is 10−5, which is usually suffiently large even if the POSCAR file has been generated with a single precision program. Increasing the SYMPREC tag means, that the positions in the POSCAR file can be less accurate.
The original method is introduced by PRL 82,2713(1999) 本文为本人工作记录。程序下载及计算方法请到 https://sourceforge.net/projects/ideal-strength-vasp/?source=navbar 感谢本程序作者： Dr. Hanyu Liu Email: hanyuliu801@gmail.com 准备工作 1.结构优化。将CONTCAR到出，在MS中导入对称性， 再转为POSCAR（for further calculations） 2.编译VASP。把VASP安装包cp到自己的目录。 for ideal tensile strength, add 'FCELL(1,1)= 0.0' to constr_cell_relax.F of vasp code. #the Stress at x axis is fixed and vasp not relax the lattice at x axis. for ideal shear strength, add 'FCELL(1,3)= 0.0' and 'FCELL(3,1) = 0.0' constr_cell_relax.F of vasp code. Here, you need to recompile vasp. 【make clean 】【make】 我在make的过程除了error： 没有lib。Lib文件放在mu'lu中，改下路径ming ，OK！ 输入文件 0. pbs.sh (集群上提交任务的脚本 ./strenth4.py Strength.log) 1.POSCAR 2. POTCAR 3. KPOINTS (不用太大： ./writekp.py 0.04) 4. strength4.py (核心程序，下载到Hanyu主页，见文章开篇) 5. input.dat (strength4.py 要读的输入文件) as followed: POSCAR #the name of POSCAR 0.02 #strain of distortion 100 #total step of distortion (虽说一般30steps即可，建议设置大一些，因为数据够了的话可以手动kill job) -45.0 35.264390 0.0 # rotate Z, Y and X. 1 # 1 tensile, 2 shear /home/zzu002/lhy/test/vasp/vasp.5.3-1/vasp # 提任务的命令。注意tensile和shear 的vasp路径不同！ 关于计算的方向： We set x axis as tensile strength direction (tensile (x)). Also we set x and z axis as shear direction (shear (x) ). For example, if you want to calculating ideal strength of Diamond along 100, you just set 0.0 0.0 0.0. If you want to calculate 110 orientation, you need set -45.0 0.0 0.0 If you want to calculate 111 orientation, you need set -45.0 35.264390 0.0. It is for rotating the x, y and z axis as a Right-handed helical rule. 我计算的tensile plane： TENSILE dire input#rotate Z, Y, X. GPa jobID jobState 100 0.0 0.0 0.0 103.115345 53421 killed 110 -45.0 0.0 0.0 100.23066 53422 killed 111 -45.0 35.264390 0.0 120.147138 53463 killed 001 0.0 90.0 0.0 75.300806 53464 killed 010 -90.0 0.0 0.0 122.560777 53465 killed 101 0.0 45.0 0.0 75.7 54584 011 -90.0 45.0 0.0 64 54585 +101 0.0 -45.0 0.0 82 54586 +110 45.0 0.0 0.0 54587 =110 killed 11+1 -45.0 -35.26439 0.0 10.961293 54953 SHEAR-001 1 0.0 90.0 0.0 55035 2 0.0 90.0 90.0 55036 3 0.0 90.0 180.0 55161 Research highlights: 1, Miao Zhang, Mingchun Lu, Yonghui Du, Lili Gao, Cheng Lu and Hanyu Liu, Hardness of FeB4?: Density functional theory investigation J. Chem. Phys., 140, 174505 (2014) 2, Yinwei Li, Jian Hao, Hanyu Liu, Siyu Lu and John S. Tse High energy density and superhard nitrogen-rich B-N compounds, Phys. Rev. Lett., 115, 105502 (2015) 3, Miao Zhang, Hanyu Liu, Quan Li, Bo Gao, Yanchao Wang, Hongdong Li, Changfeng Chen and Yanming Ma, “Superhard BC3 in cubic diamond structure” Phys. Rev. Lett., 114, 015502 (2015) 4, Quan Li, Hanyu Liu, Dan Zhou, Weitao Zheng, Zhijian Wu and Yanming Ma, novel low compressible and superhard carbon nitride: Body-centered tetragonal CN2 Phys. Chem. Chem. Phys., 14, 13081–13087 (2012)
BoltzTraP软件简介 BoltzTraP软件包是基于半经典玻尔兹曼理论，通过对能带的傅里叶积分，计算电子的群速度，进而求解体系的电子输运性能的可执行代码。BoltzTraP软件包由丹麦的Madsen教授和美国的Singh教授于2006年共同发布。该软件包的输入参数包括晶体结构参数、空间群对称操作、能量本征值以及控制参数。 因为BoltzTraP软件包里面有将不可约布里渊区能量本征值通过对称操作转换为整个布里渊区的能带色散关系 ，因此， 只需要通过计算 加密k点的不可约布里渊区的能量本征值 EIGENVAL作为BoltzTraP的输入 。由于第一性原理计算中通常会低估体系的带隙，而半导体热电材料带隙过小时计算输运性质会出现双极化现象。为了避免这个问题的产生，在输运计算中常 采用剪刀操作 ，即保持其带隙的形状不变而将其带隙设置为实验值，这个设置在BoltzTraP软件包的输入控制文件intrans中进行。控制文件intrans中还可以对计算的温度、浓度等参数进行设置。利用BoltzTraP软件包，可以计算半导体热电材料的 Seebeck系数、电导率弛豫时间的比值、电子热导率、电子热容以及霍尔系数 等。 https://www.researchgate.net/post/How_could_I_make_BoltzTrap_input_files_from_VASP_outputs/1 Dear Colleagues, Please try to change line 289 in vasp2boltz.py to if ('spacegroup' in ao.info) and (ao.info!={}): . And if you have factorization error (you will), delete duplicated k-points in case.energy. I (and my colleague who fixed it) cannot be sure it's correct. So I want you to follow this solution and please feed me back :) http://chuansong.me/n/1747114
VASP发布一种新泛函SCAN MetaGGA 转载 VASP即将发布一种新泛函SCAN MetaGGA。 功能介绍： 新泛函在诸多方面超越了LDA和PBE-GGA，在一些方面达到甚至超越了杂化泛函的可靠性，对弱相互作用体系优于LDA和PBE-GGA，计算量能与LDA和GGA相媲美，然而，SCAN的收敛性可能比GGA差。 发布时间：2017.4.18 参考文献: PRL115，036402（2015）及Nature Chemistry (2016) doi: 10.1038/nchem.2535 https://journals.aps.org/prx/pdf/10.1103/PhysRevX.6.041005 1. 介绍： The SCAN meta GGA has been shown15 to be superior to the PBE GGA for some standard molecular test sets and a small collection of solids. The mean absolute errors for SCAN15 are smaller than those for PBE by a factor of about 4 for the atomization energies of the 223 G3 molecules, a factor of 3 for the binding energies of the S22 set of weakly bound dimers of small molecules, and a factor of 4 for the LC20 set of lattice constants for solids. SCAN is also more accurate, by about 30%, in predicting the BH76 energy barriers to chemical reactions. Future studies will also show that the mean absolute errors of SCAN for the heats of formation of 94 binary solids are smaller than those of PBE by about 30%, or a factor of 3, for compounds with or without transition-metal elements, respectively. However, this Article shows that SCAN has an unexpected and striking performance for diversely bonded systems, many of which were believed to be out of reach of semilocal functionals, and is comparable to or even more accurate than a computationally more expensive hybrid GGA. 2. 如何在vasp5.4.4中实现scan+rVV10看下面链接： https://sites.google.com/site/haoweipeng/softwares
VASP手册里关于VDW的说明： 建议用 optPBE-vdW, or the optB88-vdW functional，比较精确（Peng Feng） http://cms.mpi.univie.ac.at/vasp/vasp/vdW_DF_functional_Langreth_Lundqvist_et_al.html vdW-DF functional of Langreth and Lundqvist et al. The vdW-DF proposed by Dion et al. is a non-local correlation functional that approximately accounts for dispersion interactions. In VASP the method is implemented using the algorithm of Roman-Perez and Soler which transforms the double real space integral to reciprocal space and reduces the computational effort. Several propsed versions of the method can be used: the original vdW-DF , the ``opt functionals (optPBE-vdW, optB88-vdW, and optB86b-vdW) where the exchange functionals were optimised for the correlation part , and the vdW-DF2 of Langreth and Lundqvist groups . This method is available since the 5.2.12.26May2011 version of VASP for the calculation of total energies and forces. The stress calculation for the cell optimisation (ISIF=3) is available since the VASP 5.2.12.11Nov2011 version for spin unpolarised systems and VASP 5.3.1 for spin polarised systems. N.B. : This feature has been implemented by J. Klimeš. If you make use of the vdW-DF functionals presented in this section, we ask you to cite Ref. . Please also cite the original vdW-DF paper of Dion et al. and the paper of Roman-Perez and Soler . In addtion, cite the paper of Lee et al. if you use the vdW-DF2 functional, the paper of Klimeš et al. if you use the optB88-vdW or optPBE-vdW functionals, and any other appropriate references, such as Ref. . Correlation functionals The method is invoked by setting LUSE_VDW = .TRUE. Moreover, the PBE correlation correction needs to be removed since only LDA correlation is used in the functionals. This is done by setting AGGAC = 0.0000 The two tags above need to be used for all of the following functionals. Exchange functionals The GGA tag is further used to choose the appropriate exchange functional. The original vdW-DF of Dion et al uses revPBE, therefore the vdW-DF can be chosen by setting GGA = RELUSE_VDW = .TRUE.AGGAC = 0.0000 More accurate exchange functionals for the vdW correlation functional have been proposed in and . To use these functionals set: GGA = ORLUSE_VDW = .TRUE.AGGAC = 0.0000 for optPBE-vdW, GGA = BOPARAM1 = 0.1833333333PARAM2 = 0.2200000000LUSE_VDW = .TRUE.AGGAC = 0.0000 for the optB88-vdW functional, or GGA = MK PARAM1 = 0.1234 PARAM2 = 1.0000LUSE_VDW = .TRUE.AGGAC = 0.0000 for the optB86b-vdW functional. In the vdW-DF2 functional the rPW86 exchange functional is used GGA = ML moreover, the vdW functional needs to be changed to the vdW2 correlation which requires only a change of a parameter: Zab_vdW = -1.8867 Therefore to use vdW-DF2, set: GGA = MLLUSE_VDW = .TRUE.Zab_vdW = -1.8867AGGAC = 0.0000 An overview of the performance of the different approaches can be found for example in for gas phase clusters and in for solids. Important remarks : The method needs a precalculated kernel which is distributed via the VASP download portal (VASP - src - vdw_kernel.bindat) and on the ftp server (vasp5/src/vdw_kernel.bindat). If VASP does not find this file, the kernel will be calculated. This, however, is rather demanding calculation. The kernel needs to be either copied to the VASP run directory for each calculation or can be stored in a central location and read from there. The location needs to be set in routine PHI_GENERATE. This does not work on some clusters and the kernel needs to be copied into the run directory in such cases. The distributed file uses little endian convention and won't be read on big endian machines. The big endian version of the file is available from the VASP team. There are no special POTCARs for the vdW-DF functionals and the PBE or LDA POTCARs can be used. Currently the evaluation of the vdW energy term is not done fully within the PAW method but the sum of the pseudo-valence density and partial core density is used. This approximation works rather well, as is discussed in , and the accuracy generally increases when the number of valence electrons is increased or when harder PAW datasets are used. For example, for adsorption it is recommended to compare the adsorption energy obtained with standard PAW datasets and more-electron POTCARs for both PBE calculation and vdW-DF calculation to assess the quality of the results. The spin polarised calculations are possible, but strictly speaking the non-local vdW correlation is not defined for spin-polarised systems. For spin-polarised calculation the non-local vdW correlation energy is evaluated on the sum of the spin-up and spin-down densities. The evaluation of the vdW energy requires some additional time. Most of it is spent on performing FFTs to evaluate the energy and potential. Thus the additional time is determined by the number of FFT grid points in the calculation, basically size of the simulation cell. It is almost independent on the number of the atoms in the cell. Thus the relative cost of the vdW-DF method depends on the ``filling of the cell and increases with the amount of vacuum in the cell. The relative increase is high for isolated molecules in large cells, but small for solids in smaller cells with many k-points. C( graphite ) +2I 2 = CI 4 我在没有考虑VDW的情况下算出反应的形成焓是正的，与预期相反。想着应该是CI 4 为分子，应考虑VDW的影响。同时，对于反应中的其它物质也应该考虑VDW。 ∆ H = + 2.7077 non ∆ H = + 2.75386 optPBE_VDW ∆ H = + 2.70642 optB88_ VDW 不 考虑 VDW(non) 和考虑时 ( 用两种算法 ), 单个结构的焓值有很大变化，但是反应的形成焓差别不大。总体来讲，这个反应的形成焓就是正值。
赝势相关 http://cms.mpi.univie.ac.at/vasp/vasp/pseudopotential_generation_package.html The pseudopotential generation package consists of two separate programs :the first one is called rhfsps and generates the l dependent pseudopotentials, the second one called fourpot3 prepares the pseudopotentials for VAMP and creates the POTCAR file, which can be used by VAMP. http://www.physics.rutgers.edu/gbrv/ This site hosts the GBRV pseudopotential library , a highly accurate and computationally inexpensive open source pseudopotential library which has been designed and optimized for use in high throughput DFT calculati-ons. https://www.vasp.at/index.php/news/36-highlights/100-new-release-paw-datasets After long and careful consideration we have decided to release a new set of POTCAR files covering the periodic table. This includes GW potentials for almost all elements. http://kitchingroup.cheme.cmu.edu/dft-book/dft.html 《 Modeling materials using density functional theory 》 http://blog.sciencenet.cn/blog-567091-796833.html GGA+U计算的时候除了在INCAR中设置LDAU LDAUTYPE LDAUL LDAUU LDAUJ等不需要设置其他，POTCAR中的赝势加和不加U的时候是一样的。
vdW-DF functional of Langreth and Lundqvist et al. http://cms.mpi.univie.ac.at/vasp/vasp/vdW_DF_functional_Langreth_Lundqvist_et_al.html The vdW-DF proposed by Dion et al. is a non-local correlation functional that approximately accounts for dispersion interactions. In VASP the method is implemented using the algorithm of Roman-Perez and Soler which transforms the double real space integral to reciprocal space and reduces the computational effort. Several propsed versions of the method can be used: the original vdW-DF , the ``opt functionals (optPBE-vdW, optB88-vdW, and optB86b-vdW) where the exchange functionals were optimised for the correlation part , and the vdW-DF2 of Langreth and Lundqvist groups . This method is available since the 5.2.12.26May2011 version of VASP for the calculation of total energies and forces. The stress calculation for the cell optimisation (ISIF=3) is available since the VASP 5.2.12.11Nov2011 version for spin unpolarised systems and VASP 5.3.1 for spin polarised systems. N.B. : This feature has been implemented by J. Klimeš. If you make use of the vdW-DF functionals presented in this section, we ask that you cite the following paper: J. Klimeš, D. R. Bowler, and A. Michaelides, Phys. Rev. B 83 , 195131 (2011). Correlation functionals The method is invoked by setting LUSE_VDW = .TRUE. Moreover, the PBE correlation correction needs to be removed since only LDA correlation is used in the functionals. This is done by setting AGGAC = 0.0000 The two tags above need to be used for all of the following functionals. Exchange functionals The GGA tag is further used to choose the appropriate exchange functional. The original vdW-DF of Dion et al uses revPBE, therefore the vdW-DF can be chosen by setting GGA = RE LUSE_VDW = .TRUE. AGGAC = 0.0000 More accurate exchange functionals for the vdW correlation functional have been proposed in and . To use these functionals set: GGA = OR LUSE_VDW = .TRUE. AGGAC = 0.0000 【1】 for optPBE-vdW , GGA = BO PARAM1 = 0.1833333333 PARAM2 = 0.2200000000 LUSE_VDW = .TRUE. AGGAC = 0.0000 【2】 for the optB88-vdW functional or for the optB86b-vdW functional： GGA = MK PARAM1 = 0.1234 PARAM2 = 1.0000 LUSE_VDW = .TRUE. AGGAC = 0.0000 【3】 In the vdW-DF2 functional the rPW86 exchange functional is used： GGA = ML LUSE_VDW = .TRUE. Zab_vdW = -1.8867 AGGAC = 0.0000 An overview of the performance of the different approaches can be found for example in for gas phase clusters and in for solids. Important remarks : The method needs a precalculated kernel which is distributed via the VASP download portal (VASP - src - vdw_kernel.bindat) and on the ftp server (vasp5/src/vdw_kernel.bindat). If VASP does not find this file, the kernel will be calculated. This, however, is rather demanding calculation. The kernel needs to be either copied to the VASP run directory for each calculation or can be stored in a central location and read from there. The location needs to be set in routine PHI_GENERATE. This does not work on some clusters and the kernel needs to be copied into the run directory in such cases. The distributed file uses little endian convention and won't be read on big endian machines. The big endian version of the file is available from the VASP team. There are no special POTCARs for the vdW DF functionals and the PBE or LDA POTCARs can be used. Currently the evaluation of the vdW energy term is not done fully within the PAW method but the sum of the pseudo-valence density and partial core density is used. This approximation works rather well, as is discussed in , and the accuracy generally increases when the number of valence electrons is increased or when harder PAW datasets are used. For example, for adsorption it is recommended to compare the adsorption energy obtained with standard PAW datasets and more electron POTCARs for both PBE calculation and vdW DF calculation to assess the quality of the results. The spin polarised calculations are possible, but strictly speaking the nonlocal vdW correlation is not defined for spinpolarised systems. For spinpolarised calculation the nonlocal vdW correlation energy is evaluated on the sum of the spin-up and spin-down densities. The evaluation of the vdW energy requires some additional time. Most of it is spent on performing FFTs to evaluate the energy and potential. Thus the additional time is determined by the number of FFT grid points in the calculation, basically size of the simulation cell. It is almost independent on the number of the atoms in the cell. Thus the relative cost of the vdW DF method depends on the“filling of the cell and increases with the amount of vacuum in the cell. The relative increase is high for isolated molecules in large cells, but small for solids in smaller cells with many k-points. This feature has been implemented by J. Klimeš. If you make use of the vdW-DF functionals presented in this section, we ask that you cite the following paper: J. Klimeš, D. R. Bowler, and A. Michaelides, Phys. Rev. B 83 , 195131 (2011).
研究高压结构时，经常需考察同一结构在不同压力下的键长、晶格参数、体积及相对能量的变化，此时需要对所研究的结构在不同压力下进行优化。在优化前务必对KPOINTS进行测试： 问题：低压下经测试使能量收敛的k点也试用于高压下的结构优化吗？ 对于下面的方式一产生的KPOINTS不行！原则上高压下经测试使能量收敛的k点可适用于低压。 因为，高压结构拥有更小的晶格参数，对应更长的k空间参数，若在高压下采用N1、N2、N3的数值分割k空间，且间隔足够密；则低压结构采用该k点进行优化时，对k空间的分割将更密(因为k空间长度小)，能量将更容易收敛。 对于方式二产生的KPOINTS，由于空间的k mesh密度始终维持不变，则有可能低压下测试的KPOINS可能适用于高压，但仍需具体考察。以第二种KPOINYS产生方式为例，考察结果见下面“KPOINTS设置考察”。 1. KPOINTS产生方式之一 如下面的KPOINTS文件： 4x4x4 ！ Comment 0 0 ！ automatic generation of k-points Monkhorst ！ M use Monkhorst Pack 4 4 4 ！ grid 4x4x4 0 0 0 ！ shift (usually 0 0 0) 如该采用奇数分割，则可将Gamma点包含进来。 2. KPOINTS产生方式之二 Automatic mesh 0 ! number of k-points = 0 -automatic generation scheme Auto ! fully automatic 10 ! length (l) KPOINTS文件最后一行即表示a linear k-point density of 20 A， 其含义为每1/A上的k点数， 注意其单位为angstrom(A=1/(1/A))，对于绝缘体通常取20，金属通常取40。【低压下测试的k密度试用于高压相】 例如正空间a=3A，倒空间2pai/3=~2/A，即分割数为20；若在压力下a变为1.5A，则倒空间变为 2pai/1.5 ~ 4/A，则分割数为40； 注意这是对倒空间上的三个矢量采用相等k密度进行分割，因此，视晶格参数不同，三个方向的分割数（ subdivisions ）可能不同。 该KPOINTS已经包含了Gama Point As before, the first line is treated as a comment. On the second line a number smaller or equal 0 must be specified. In the previous section, this value supplied the number of k-points, a zero in this line activates the automatic generation scheme. The fully automatic scheme, selected by the first character in the third line (’a’), generates G centered Monkhorst-Pack grids, where the numbers of subdivisions along each reciprocal lattice vector are given by Symmetry is used to map equivalent k-points to each other, which can reduce the total number of k-points significantly. Useful values for the length vary between 10 (large gap insulators) and 100 (d-metals). 3. KPOINTS产生方式之三 A slightly enhanced version, allows to supply the numbers for the subdivisions N1, N2 and N3 manually: Automatic mesh 0 ! number of k-points = 0 -automatic generation scheme Gamma ! generate a Gamma centered grid 4 4 4 ! subdivisions N_1, N_2 and N_3 along recipr. l. vectors 0. 0. 0. ! optional shift of the mesh (s_1, s_2, s_3) In this case, the third line (again, only the first character is significant) might start with ’G’ or ’g’ —for generating meshes with their origin at the G point (as above)— or ’M’ or ’m’, which selects the original Monkhorst-Pack scheme. In the latter case k-point grids【指， Monkhorst-Pack scheme？ 】, with even (mod(Ni,2) = 0) subdivisions are shifted off G: 4. hexagonal lattices We strongly recommend to use only Gamma centered grids for hexagonal lattices . Many tests we have performed indicate that the energy converges significantly faster with Gamma centered grids than with standard Monkhorst Pack grids. Grids generated with the “M” setting in the third line , in fact do not have full hexagonal symmerty. 5. 摘录 一般如非必要，可以先用自动模式生成 k 点， VASP 会自动生成一个简约化后的 k 点矩阵，存于 IBZKPT 文件，可以直接复制里面的数据到 KPOINTS 文件中使用 ，这也是该输入法的主要用途，可以减少重复自动生成格点的时间。另一个用途是为了做精确的 DOS （ Density of status ）计算，由于这类计算所需的 k 点数极大，通过全手动尽可能的优化 k 点也就必需了。 最后提醒一点，VASP的帮助文档特别提醒，对于六方晶系，不要用M来自 动生成格点，而要用G。 6. 不同KPOINTS设置下OUTCAR文件中k mesh查看 （1） Monkhorst 方式 K-Points 0 Monkhorst Pack 15 15 15 0 0 0 OUTCAR ..... KPOINTS: K-Points Automatic generation of k-mesh. Space group operators: ....... Dimension of arrays: k-points NKPTS = 120 k-points in BZ NKDIM = 120 number of bands NBANDS= 16 number of dos NEDOS = 301 number of ions NIONS = 4 non local maximal LDIM = 6 non local SUM 2l+1 LMDIM = 18 total plane-waves NPLWV = 8000 max r-space proj IRMAX = 1 max aug-charges IRDMAX= 19027 dimension x,y,z NGX = 20 NGY = 20 NGZ = 20 dimension x,y,z NGXF= 40 NGYF= 40 NGZF= 40 support grid NGXF= 40 NGYF= 40 NGZF= 40 ions per type = 3 1 NGX,Y,Z is equivalent to a cutoff of 11.15, 11.15, 11.15 a.u. NGXF,Y,Z is equivalent to a cutoff of 22.29, 22.29, 22.29 a.u. I would recommend the setting: dimension x,y,z NGX = 19 NGY = 19 NGZ = 19 （2） Auto 方式 Automatic mesh 0 ! number of k-points = 0 -automatic generation scheme Auto ! fully automatic 50 ! length (l) OUTCAR ..... KPOINTS: Automatic mesh Automatic generation of k-mesh. generate k-points for: 21 21 21 ..... Dimension of arrays: k-points NKPTS = 286 k-points in BZ NKDIM = 286 number of bands NBANDS= 16 number of dos NEDOS = 301 number of ions NIONS = 4 non local maximal LDIM = 6 non local SUM 2l+1 LMDIM = 18 total plane-waves NPLWV = 8000 max r-space proj IRMAX = 1 max aug-charges IRDMAX= 19106 dimension x,y,z NGX = 20 NGY = 20 NGZ = 20 dimension x,y,z NGXF= 40 NGYF= 40 NGZF= 40 support grid NGXF= 40 NGYF= 40 NGZF= 40 ions per type = 3 1 NGX,Y,Z is equivalent to a cutoff of 11.16, 11.16, 11.16 a.u. NGXF,Y,Z is equivalent to a cutoff of 22.32, 22.32, 22.32 a.u. I would recommend the setting: dimension x,y,z NGX = 19 NGY = 19 NGZ = 19 ....... （3） Gama 方式 --待补充 7. KPOINTS测试结果 KPOINTS 测试 Im-3m:600G （ 1 ）采用 Mateiralstudio 中的 fine 参数来设置 KPOINTS ./mskptest.x CONTCAR_6000 Thethree qualities of k-point separation for CASTEP (1/Angstrom) coarse.le..08 medium.le..05 fine.le..04 Please input thequality of Monkhorst-Pack grid (Default is .035) 0.03 Reciprocal lattice parameters 2.11841778 2.11841778 2.11841778 Meshparameters of Monkhorst-Pack grid 12 12 12 Actual spacing of Monkhorst-Pack grid 0.02809639 0.02809639 0.02809639 Reciprocal Lattice -1.49860797 -0.86886324 1.21940342 -0.00576644 1.72859416 1.22459088 -1.49727578 0.86297900 -1.22520461 Real lattice parameters 2.96598026 2.96598026 2.96598026 Cell Angles 90.00000000 90.00000000 90.00000000 Real Lattice -2.09818937 -1.21648867 1.70727724 -0.00807354 2.42019125 1.71454018 -2.09632419 1.20825019 -1.71539946 (2) 实际的 KPOINTS 测试过程，可见 KPOINTS 对能量的影响可达几个 meV ，特别是对于低压结构，若选择 30 Auto 则，则与选择 40 Auto 能量相差接近 20 meV 。 ../getEntafterkptest.sh OUTCAR_30 5987.85 87.484625 OUTCAR_40 5987.71 87.487188 OUTCAR_50 5995.14 87.491954 OUTCAR_60 5996.47 87.495100 OUTCAR_70 5994.86 87.494637 OUTCAR_75 5994.68 87.491468 OUTCAR_80 5994.63 87.493970 OUTCAR_90 5994.78 87.494334 OUTCAR_95 5995.41 87.494084 OUTCAR_30:Automatic generation of k-mesh. OUTCAR_30: generate k-points for: 10 10 10 OUTCAR_40:Automatic generation of k-mesh. OUTCAR_40: generate k-points for: 13 13 13 OUTCAR_50:Automatic generation of k-mesh. OUTCAR_50: generate k-points for: 17 17 17 OUTCAR_60:Automatic generation of k-mesh. OUTCAR_60: generate k-points for: 20 20 20 OUTCAR_70:Automatic generation of k-mesh. OUTCAR_70: generate k-points for: 24 24 24 OUTCAR_75:Automatic generation of k-mesh. OUTCAR_75: generate k-points for: 25 25 25 OUTCAR_80:Automatic generation of k-mesh. OUTCAR_80: generate k-points for: 27 27 27 OUTCAR_90:Automatic generation of k-mesh. OUTCAR_90: generate k-points for: 30 30 30 OUTCAR_95:Automatic generation of k-mesh. OUTCAR_95: generate k-points for: 32 32 32 grep NKPT OUTCAR_* OUTCAR_30: k-points NKPTS = 56 k-points in BZ OUTCAR_40: k-points NKPTS = 84 k-points in BZ OUTCAR_50: k-points NKPTS = 165 k-points in BZ OUTCAR_60: k-points NKPTS = 286 k-points in BZ OUTCAR_70: k-points NKPTS = 455 k-points in BZ OUTCAR_75: k-points NKPTS = 455 k-points in BZ OUTCAR_80: k-points NKPTS = 560 k-points in BZ OUTCAR_90: k-points NKPTS = 816 k-points in BZ OUTCAR_95: k-points NKPTS = 969 k-points in BZ Im-3m:0G Automatic mesh 0 Auto L (L=30,40,50,60) ../getEntafterkptest.sh OUTCAR_30 21.78 -48.119067 OUTCAR_40 22.89 -48.148892 OUTCAR_50 22.53 -48.142697 OUTCAR_60 23.06 -48.140251 OUTCAR_70 22.06 -48.144940 OUTCAR_75 22.96 -48.142481 OUTCAR_80 22.2 -48.144701 OUTCAR_90 22.98 -48.142411 OUTCAR_95 22.95 -48.143537 grep gene OUTCAR_* OUTCAR_30:Automatic generation of k-mesh. OUTCAR_30: generate k-points for: 7 7 7 OUTCAR_40:Automatic generation of k-mesh. OUTCAR_40: generate k-points for: 10 10 10 OUTCAR_50:Automatic generation of k-mesh. OUTCAR_50: generate k-points for: 12 12 12 OUTCAR_60:Automatic generation of k-mesh. OUTCAR_60: generate k-points for: 15 15 15 OUTCAR_70:Automatic generation of k-mesh. OUTCAR_70: generate k-points for: 17 17 17 OUTCAR_75:Automatic generation of k-mesh. OUTCAR_75: generate k-points for: 18 18 18 OUTCAR_80:Automatic generation of k-mesh. OUTCAR_80: generate k-points for: 20 20 20 OUTCAR_90:Automatic generation of k-mesh. OUTCAR_90: generate k-points for: 22 22 22 OUTCAR_95:Automatic generation of k-mesh. OUTCAR_95: generatek-points for: 23 23 23 grep NKPT OUTCAR_* OUTCAR_30: k-points NKPTS = 20 k-points in BZ OUTCAR_40: k-points NKPTS = 56 k-points in BZ OUTCAR_50: k-points NKPTS = 84 k-points in BZ OUTCAR_60: k-points NKPTS = 120 k-points in BZ OUTCAR_70: k-points NKPTS = 165 k-points in BZ OUTCAR_75: k-points NKPTS = 220 k-points in BZ OUTCAR_80: k-points NKPTS = 286 k-points in BZ OUTCAR_90: k-points NKPTS = 364 k-points in BZ OUTCAR_95: k-points NKPTS = 364 k-points in BZ 手册注释 NGX, NGY, NGZ FFT mesh for wavefunctions describes the number of grid points in x, y and z direction, and NPLWV the total number of points (i.e.NGX*NGY*NGZ) Most quantities (like charge densities) are defined on these 3-dimensional grids. In the sequential version NGX, NGY and NGZ were sufficient to perform a three dimensional FFT of quantities defined on these grids. NGXF,NGYF,NGZF FFT mesh for charge s NGX, NGY, NGZ controls the number of grid-points in the FFT-mesh into the direction of the three lattice-vectors. X corresponds to the first, Y to the second and Z to the third lattice-vector (X,Y and Z are not connected with cartesian coordinates, don’t be fooled by the historical naming conventions). NGXF, NGYF, NGZF controls the number of grid-points for a second finer FFT-mesh. On this mesh the localized augmentation charges are represented, if ultrasoft (US) Vanderbilt potentials or the PAW method are used. In addition, local potentials (exchange-correlation, Hartree-potential and ionic potentials) are also calculated on this second finer FFT-mesh if (and only if) US-pseudopotentials are used. 转载自： http://blog.sciencenet.cn/blog-567091-680538.html
VASP manual site http://cms.mpi.univie.ac.at/vasp/vasp/vasp.html LINK for VASP tips (update: 2015-04-27) David Karhánek's Homepage ( http://homepage.univie.ac.at/david.karhanek/downloads.html ) Nano Group from Budapest ( http://wiki.kfki.hu/nano/VASP_relevant_resources ) WaveTrans: Real-space wavefunctions from VASP WAVECAR file (Dr. R. M. Feenstra and Dr. M. Widom) ( http://www.andrew.cmu.edu/user/feenstra/wavetrans/ ) JR Kitchin, Modeling materials using density functional theory ( http://kitchingroup.cheme.cmu.edu/dft-book/dft.html ) Skelton's tips: Phonons Phonopy: -Pro Tips (2014) ( http://www.slideshare.net/jmskelton/phonons-phonopy-pro-tips-2014 ) -Phonons Phonopy: Pro Tips (2015) ( http://www.slideshare.net/jmskelton/phonons-phonopy-pro-tips-2015 ) -VASP: Some Accumulated Wisdom ( http://www.slideshare.net/jmskelton/vasp-some-accumulated-wisdom ) -VASP And Wannier90: A Quick Tutorial ( http://www.slideshare.net/jmskelton/vasp-and-wannier90-a-quick-tutorial ) -Wannier90: Band Structures, Tips and Tricks ( http://www.slideshare.net/jmskelton/wannier90-band-structures-tips-and-tricks ) Atomic structure optimization : to obtain optimized structures (update: 2013-08-24) - In many cases, the atomic positions are drifted if one use default option for structure relaxations. In other words, atomic forces are hardly reaching small values due to the uncertainty in calculation. - To obtain optimized geometry, it is suggested to use the following tag in INCAR. ADDGRID = T - Or you can increase calculation precision by lowering EDIFF or by increasing NELMIN It will increase the precision of the force and thereby you wil obtain the optimized structure. NELMIN = 6 Importance of k-space projection (LREAL=FALSE) - Not same energies for same structures, LREAL=Auto or LREAL=TRUE - Sometimes, the same structures can give different total energies if either LREAL=TRUE or LREAL=Auto are used. For example, one makes same structure but the atomic basis are equally displaced by constant vector. Many cases it is ok, because the energy difference is small, less than few meV/atom. However, when one calculates the formation energy, the error can be large if the supercell size is very big. If it happens, you should re-calculate the total energy using k-space projection scheme. (LREAL = False) - When you have to calculate very precise calculation for phonon mode or for small barrier calculations, it is strongly recommended to use the LREAL=F tag, with other high precision settings even though the supercell size is very large. Unit of density-of-states (DOS) in DOSCAR file (update: 2013-08-24) - The unit of DOS is the states per electrons. If you run the single interation, you will always obtain states/eV - However, sometimes you will get different unit in DOSCAR. In many cases, it happens becasue you run the molecular dynamics simulation or structure relaxation. Then the DOSCAR contatins the sum or average value of DOSs. Note that the NSW value is important while the number of interation is not. Density-of-states (DOS) on the vacant site (DOSCAR on the vacant site) (from here, http://cms.mpi.univie.ac.at/vasp-forum/forum_viewtopic.php?4.1692 site broken ) - One can project the electronic wavefunction on arbitrary positions. It is very helpful to find the integrated charge for given spheres without postprocessing of charge density files (CHG or CHGCAR) - To do that, you should add following lines in POSCAR file, with the use of proper RWIGS, NPAR, LORBIT tags in the INCAR file. Empty (spheres) 1 0.0000000000000000 0.000000000000000 0.0000000000000000 Drawing local density of states (LDOS) (update: 2013-10-27) - Local density of states is the real space distribution of the wavefunction square, i.e. charge density for given energy range. Following below steps to obtain the LDOS - Obtain wavefunction - obtain all partial charge density, rho(n,k,x,y,z) - You can obtain energy eigenvalue from EIGENVAL - Now you know rho(E(n,k),x,y,z) - Now average over rho for y,z. to obtain rho(E,x) - Now draw rho(E,x) with smearing. - Note that you need very many k-point. I know and you know that it is the time consuming job. Local-potential (LOCPOT) - Frequently, one should check the difference between LOCPOT files - In this case, it is a good idea to fix the number of charge-density or local-potential grids. - Related tags are NGX, NGY, NGZ. Maybe sometimes with NGXF, NGYF, NGZF Phonon calculations - Phonopy code is the very useful tool for phonon calculations. - Download from here: http://phonopy.sourceforge.net/ - Good texts to read: ftp://ftp.heanet.ie/.../introduction-phonon-calc.pdf (broken site) http://icms3.weebly.com/uploads/3/5/9/0/3590130/version1.pdf IR-active mode calculations - You can calculate the IR-active mode intensities using the following code. - Download from here: http://homepage.univie.ac.at/david.karhanek/downloads.html Dielectric function for metallic system - Dielectric functions are composed of interband-transition-term (bound electron term) plus intraband-transition-term (free electron term: Drude term) - Although you can obtain inerband-transition term from VASP, you should calculate intraband-transition-term by yourself especially for metallic system. - It means that you need plasma frequency from the band structure. There are two possible method to calculate the plasma freqeuncy using vasp code. Dielectric constant for dielectrics with clamped-ion (update: 2016-01-28) - Using the LEPSILON-tag, one can calculate the ion-clamped dielectric constant as written in VASP manual . - Optimizing structure - Calculating wavefunction - Turning on LEPSILON-tag (LEPSILON = TRUE) and commenting out NPAR-tag or no NPAR tag in the INCAR-file. Dielectric constant for dielectrics with relaxed-ion (update: 2016-01-28) - Using the LEPSILON-tag IBRION-tag, one can calculate the ion-clamped dielectric constant as written in VASP manual . - Optimizing structure. Note that ionic forces should be very small because of the force constant calculations. - Calculating wavefunction - LEPSILON = TRUE / no NPAR tag / IBRION = 5,6,7,8 (one of them) Prediction of plasma freqeuncy (update: 2013-10-27) - The plasma frequency I mention here is not the 'Screened plasma frequency'. You shoud distinguish between just plasma frequency and screened plasma frequency. The screened plasma frequency is the frequency where the real part of the dielctric function (inter-band term + intra-band term) becomes zero. And the screened plasma frequency is smaller than the plasma frequency. - Using BoltzTraP code, you can obtain electrical conductivity over electon scattering time (SIGMA/TAU) within constant relaxation time approximation. Then, SIGMA/TAU can be transformed to plasma frequency. You shoud be careful when calculating plasma frequency. As the plasma frequency is obtaiend from the intra-band transition, I mean the group velocity of the electron, the k-point sampling is very important. You need sufficiently many k-points to obtain correct band dispersion. It is recommended to test the energy convergence with increasing number of k-points. - You may also obtain the plasma frequency within vasp code. With the input tag of LOPTICS=TRUE you can find the plasma frequency in the OUTCAR file. dx:一定要注意这里的弛豫时间是电子的弛豫时间，涉及电子-电子相互作用和电子-声子相互作用，和shengBTE里面的声子-声子弛豫时间完全不是一回事。 Electrical conductivity using VASP code (update: 2013-11-13, 2015-04-17) - Within semiclassical theory, one may obtain the electrical conductivity from the band structure by calculating the intraband structure. However, when you use the semiclassical theory, you should get the relaxation time from outside or you should calculate the relaxation time. If you already know the electron scattering time or relaxation time, you can get the electrical conductivity using VASP. - There are two ways to get the electrical conductivity. One is using BoltzTraP code. And the other thing is to read the OUTCAR file from VASP code. Using LOPTICS = TRUE tag and RTIME = x.xx (femtotsec) in INCAR with sufficiently many kpts, you can finally found the reliable value of electrical conductivity tensor in OUTCAR from VASP. Unit of electrical conductivity per relaxation time is 10^22 . For example, if the tensor is ], then Sigma/Tau = 0.020 x 10^22 S/m/sec = 2 x 10^20 S/m/sec For recent VASP version of 5.3.3, you can set the relaxation time ( RTIME in unit of femtosec ) and obtain the electrical conductivity tensor in unit of mega S/m written in OUTCAR. // 2015-04-17 And the Energy dependent conductivity is also written in the vasprun.xml file. // 2015-04-17 MAGMOM tag in INCAR - You may reduce the sentence length of MAGMOM tag as shown below, MAGMOM = 1 1 1 1 1 == MAGMOM = 5*1 MAGMOM = 1 1 -3 -3 4 == MAGMOM = 2*1 2*-3 4 It really saves your time to make INCAR file for spin poloarized calculation of large supercell Many body perturbation theory (MBPT) : GW and BSE in VASP (I tested MBPT calculations with vasp version of 5.2.xx) - MBPT calculations are very time consuming and very memory consuming. - BSE tips in other site : http://cms.mpi.univie.ac.at/vasp-forum/forum_viewtopic.php?4.12605 (site broken) - Calculation flow: (1) DFT(Exc=LDA or GGA) calculations (or DFT+U, Hybrid-DFT) for wavefunctions (WAVECAR) and their derivatives (WAVEDER). You need sufficient many unoccupied bands. (2) GW calculations (ALGO=GW0, scGW0, scGW) (3) BSE calcualtoin (ALGO=BSE) : although you can not find the BSE tag in your manual, it works. And check the vasprun.xml which contains frequency dependent dielectric functions. - Related tags: LOPTICS = T for derivative of wavefunction (see manual) NBANDS number of total bands including occupied and unoccupied, convergence test is needed (see manual) ENCUT convergence test is needed (see manual) CSHIFT smearing value for optical response results (freq. dependent dielectric function) (see manual) ALGO (=GW0, scGW0, scGW, BSE) one shot G0W0 tag (GW0), W fixed GW tag (scGW0), fully self consistent GW (scGW), BSE calculation tag (see manual) NOMEGA (see manual) NBANDSO number of occupied bands you considered for BSE calculations NBANDSV number of unoccipied (virtual) bands you considered for BSE calculations Charge analysis (update: 2015-04-20) - Bader Charge You may analyze the atomic charge density distribution by using Bader Analysis. The concept of Bader Charge is based on the space separation based on the minimum charge density positions. - CHGCAR or CHG file You can directly analyze the atomic charge density using CHGCAR file. Or as mentioned in above (#6), you can get the atomic charge using DOSCAR or OUTCAR. Unit conversion for VASP phonon frequency from density functional perturbation theory (DFPT) results, Unit of phonopy (update: 2015-06-24) - When one run DFPT calculations or BSE calculations, one obtain the vibration frequency in OUTCAR or vasprun.xml files.VASP. freq (THz) : ang. moment w (2pi THz) : 1/wavelength (1/cm) : energy (meV) = 1 : 2 pi : 10^10/c(speed of light: c) : 10^15xh(planck const: h) = 1 : 6.238186 : 33.35641 : 4.135669 - the default unit for frequency is THz. Therefore you can obtain phonon energy by producting planckconstant (h) 1 THz = 4.135669 meV - example of PbTe phonon (VASP + phonopy, cubic lattice, lattice parameter of a=6.567843. First principles calculations of Lattice thermal conductivity - Thermal conductivity (kappa_tot) is a sum of electrical term and lattice term. - Electrical term (kappa_elec) = (Lorenz Number) x ( Elec. conductivitiy ) x (Temp.), where you can obtain Lorenz number and Elec. conductivity from the BoltzTrrap (electron boltzmann transport equation) - Lattice therm (kappa_latt) = (kappa_phonon) can be obtained by performing phono3py calculations, which calculates the transition rate of phonon-phonon scattering. ( phono3py LINK ) Using phono3py code with VASP (2015-12-21) (0) prepare the optimized supercell with the VASP parameter to be used in force calculations. The force should be less than 1.0e-08 eV/angstrom. Of course, you should optimized the structure considering the k-point mesh you will use in (1). WAVECAR file is needed to reduce computational resources for step (2). (1) generate POSCAR files with pair-distance configurations. There might be one to ten thousands POSCAR files. phono3py -d --dim=4 3 3 -c POSCAR-unitcell --cutoff_pair_distance=4.4 supercell size, unitcell file, cutoff pair-distance in angstrom (2) calculate force for pair-distance configurations using VASP. Make directories for configurations of POSCAR-00365. Copy INCAR, KPOINTS, POTCAR POSCAR, and WAVECAR. Using WAVECAR preapared from (0) will be useful to reduce the computational cost. In my case, I used the gamma point only but with the large supercell containing 200-400 atoms. Energy cutoff is also known to be sensitive to the lattice thermal conductivity tensor. For example, when the monoclinic unitcell volume is optimized within 400 eV and internal coordinate is optimized within 300 eV. then there could be non-zero off-digonal term in lattice thermal conductivity tensor. (3) collect vasprun.xml file phono3py --cf3 disp-{00001..00999}/vasprun.xml (or) phono3py --cf3 disp-*/vasprun.xml disp_fc3.yaml, POSCAR-unitcell files should be in same directory. (4) create fc2.hdf and fc3.hdf phono3py --dim=4 3 3 -c POSCAR-unitcell This step is not madatory, but you can avoid calculating fc2 and fc3 at every run tim. from http://phonopy.sourceforge.net/phono3py/workflow.html The file size of fc3.hdf is very large. It is about 1 to 10 GB. (5) calculate thermal conductivity (5.1) split thermal 3phonon process configurations phono3py --fc3 --fc2 --dim=4 3 3 -v --mesh=11 11 11 -c POSCAR-unitcell --br --thm ---wgp q-mesh grid, write grid point for split calculation The computational time for 3phonon-transition rates are very time consuming when q-mesh grid is fine. Thermal conductivity, especially for low temperature region, q-mesh grid affects the lattice thermal conductivity. Therefore it is recommended to split the calculation set and calculate the transition rate separately using many cpu cores. After run the commands in (5.1), you will get grid_address-m111111 and ir_grid_points.yaml (5.2) calculate 3phonon transition rate for each q grid point phono3py --fc3 --fc2 --dim=4 3 3 -v --mesh=11 11 11 -c POSCAR-unitcell --br --thm --write_gamma --gp=3 Check q grid point number (gp) from ir_grid_points.yaml and write correct gp. The calculation time is dependent on q-mesh size. Each gp-point calculation time is varying from 10 minutes to 10 days depending supercell size and q-mesh size. Note that there are many gp if q-mesh size is fine. (5.3) sum up thermal conductivity phono3py --fc3 --fc2 --dim=4 3 3 -v --mesh=11 11 11 -c POSCAR-unitcell --br --thm --read_gamma z_TC_def.txt Now you will get the lattice thermal conductivitiy file, z_TC_def.txt Type tail -n 150 z_TC_def.txt and find the lattice conductivity tensor as a function of temperature. (*) Useful tag: --isotope, --mv, --bmfp for isotope effect, mass variance effect, boundary scattering effect. (**) Auxilary tools kaccum and gaccum commands for accumulated lattice thermal conductivity and 3ph transition rate, respectively. Band structure using VASP+Wannier90 (2016-04-05) We can draw band structure by calculatinng Wannier Interpolation using Wannier90. You need wannier90 compiled VASP. Although I don't know well about wannier90, I want to share my knowledge to prevent other's failures Simple method (1-1) Preconverge wavefunction (1-2) Turn on LWANNIER90=T tag and run wan90 compiled VASP. Then Maximally Localized Wannier Functions (MLWFs) are generated with wannier90.win, wannier90.mmn, wannier90.eig, wannier90.out files. (1-3) Re-Run wannier90. wannier90.x wannier90 (1-4) Interpolate wannier functionals: first modify the wannier90.win file. For details, see the tutorial file. wannier90.x wannier90 GW band structure from wannier90 I will explain how to make GW band structure. You considered diamond Si and there are two Si atoms in primitive cell. Then there are 8 electrons. For GW calculations you may need many number of bands, about 96 or larger than it. Please note that NPAR should be tagged out or equal to number of nodes(cores). And you may know that NBANDS is dependent on NPAR tag. First check NBANDS without NPAR tag. Also you don't need too many wannier functions. It is a good idea to reduce the number of bands. Actually the number of bands for wannier90 is very sensitive to the band structure. So I recommend to check the suitable number of bands by just running simple DFT. Here I assume that the optimal number of wannier bands is 4 per Si atom (8 per primitive cell). In shortly, you need to determine NBANDS for VASP and NUM_WAN for wannier90. Here I will use NBANDS=96 (I have 12 core nodes, so NBANDS should be divived by 12) Here I will use NUM_WAN = 8 (the number of wannier functions are related to the symmetries) (1-1) DFT calculations with default NBANDS (1-2) Optical calculations (LOPTICS=TRUE, NBANDS=96 ) (1-3) GW calculations with wannier90 (ALGO=GW0, LWANNIER90=TRUE, no NPAR related tag) If you do not set wannier90.win file, then wannier90.win will be generated. But it is too heavy to play with it since there are too many BANDs. Make wannier90.win file (or modify wannier90.win file) as follows. (The default wannier90.win can be obtained by running without wannier90.win file. GW calculation is too time consuming. Therefore just run DFT with LWANNIER90=T tag) num_wann = 8 num_bands = 8 exclude_bands : 9-96 (1-4) Run wannier90 After getting wannier90.win, wannier90.eig, wannier90.mmn, rerun the wannier90.x cmd wannier90.x wannier90 Now you wannierize and obtain MLWFs. (1-5) Intepolate First modify or generate or make suitable wannier90.win or NAME.win. And run by typing wannier90.x wannier90 or wannier90.x NAME. Now you can interpolate for DOS, BANDS, or Fermi-Surface. IMPORTANT TIP for wannier90 with VASP (2016-11-02): you should comment out the NPAR tag when the calculation is crashed. When I ran the vasp calculation for Mg2Si primitive cell, my calculation is crashed with following error: internal error in GENERATE_KPOINTS_TRANS: G vector not found. By removing NPAR tag, I did succeed to calculate the wannier90 with VASP. Phonon thermal conductivity of low dimensional structure (2016/04/25) Please be careful that if there is vacuum, I mean your structure is 1d or 2d structures with vacuum, you should normalize the thermal conductivity tensors. Also be careful that there will be dipole interactions between adjacent supercell and the force from displacement supercell will be dependent on the lateral size of supercell. For example, I tested 1D carbon (C) atomic chain along z direction and the C-C distance is set to 1.28 angstrom. I considered 1x1x10 supercells containing 10 C atoms with various a=b values, 10, 14, 20, and 25 angstroms. I find that the thermal conductivity value of zz component times a square seems to be converged as a becoming larger from 10 to 25 (k_zz x a^2 = 2836, 2968, 3136, 3215 Wxm/K at 300K, 4560, 5181, 5281, 4685 Wxm/K at 500 K. It implies that the cross section area should be checked carefully because the phono3py code assumes that one's structure is bulk. The reason why the k_zz is not exactly same is that there is dipole when single C atom is displaced and it makes large cross sectional interactions between adjacent supercells. I think that if we consider the neutral atomic chain which has negligible dipole interaction between supercell, the k_zz will be converged more fastly. 原文链接：https://sites.google.com/site/cta4rbk/home/tips4vasp，因为google被墙了，所以在这里留个底。
对于一些磁性体系、镧系和锕系元素及相关化合物的静态计算（电子迭代），经常会遇到 “ 难收敛 ” 的问题。 下面给出几个相关 Flag 及设置方法： 1 、 LMAXMIX Default: LMAXMIX = 2 An additional flag controls up to which l quantum number the onsite PAW charge densities are passed through the charge density mixer. Higher l-quantum numbers are usually not handled by the mixer. In order to obtain fast convergence to the groundstate, you can try the following setting: LMAXMIX = 4 for d elements LMAXMIX = 6 for f elements 这个 FLAG 对于含 d 电子和 f 电子的体系是非常重要的，很大一部分体系的收敛问题可以通过设置合适的 LMAXMIX 值来解决。 2 、 ALGO, IALGO, LDIAG If the self-consistency loop does not converge within 40 steps, it will probably not converge at all. In this case you should reconsider the tags IALGO, LDIAG, and the mixing-parameters. 这是说明书上的建议。 一般情况下，或使用 IALGO=48 时遇到收敛问题的话，可以考虑设 IALGO 为 38 （ 4.5 以前的版本可设为 8 ），或设置 ALGO=Normal or Fast (in VAS P.4.5 and later versions) 。 3 、 NELMDL NELMDL gives the number of non-selfconsistent steps at the beginning; if one initializes the wave functions randomly the initial wave functions are far from anything reasonable. The resulting charge density is also 'nonsense'. Therefore it makes sense to keep the initial Hamiltonian, which corresponds to the superposition of atomic charge densities, fixed during the first few steps. Choosing a 'delay' for starting the charge density update becomes essential in all cases where the SC-convergence is very bad (e.g. surfaces or molecules/clusters chains). Without setting a delay VASP will probably not converge or at least the convergence speed is slowed down. NELMDL might be positive or negative. A positive number means that a delay is applied after each ionic movement — in general not a convenient option. A negative value results in a delay only for the start-configuration. 4 、 mixing-parameters 对于一些难收敛的体系，可以使用 “linear mixing”, 具体详见 VASP 说明书中的 “Mixing-tags” 。 For an initial linear mixing (BMIX ~ 0) an optimal setting for A(AMIX) can be found easily by setting Aopt=Acurrent*Γmean. For the Kerker scheme either A or q0(i.e. AMIX or BMIX) can be optimized, but we recommend to change only BMIX and keep AMIX fixed (you must decrease BMIX if the mean eigenvalue is larger than one, and increase BMIX if the mean eigenvalue is smaller than one). 尽管 VASP 说明书中给出了调节 AMIX 和 BMIX 的一些较为明确的建议，但是实际去调节的时候，还是挺难的，但原则上说，是可以通过调节这两个 Flag 来使得收敛问题得以解决的，只是得有耐心。 5 、 kmesh, SIGMA 收敛问题还跟 kmesh 及 SIGMA( 当使用 ISMEAR 不等于 -5 和 -4 时 ) 的设置有关。要达到同样的精度，较小的 SIGMA 则需要较大的 kmesh ；而且，当 SIGMA 较小时，若 kpoints 不够多，也会出现难收敛的情况。
1 常规DFT计算得到WAVCAR 2 HSE06自洽 LHFCALC = .TRUE. ; HFSCREEN = 0.2 ; AEXX = 0.25;PRECFOCK= F ALGO = D ; 比ALGO=N快点 TIME = 0.4 ; LDIAG = .TRUE. LPLANE=.TRUE. NPAR=96 NSIM=4 KPAR=4 几点经验： （用96个核在cray上计算，0.02的K-points，半小时左右就计算完了。） cray 上的aprun并行效率会比openmpi快2倍以上（HSE06的计算上）。 2D层状材料，15与25的真空层，HSE06计算时间会小2倍多。 HSE06自洽完算能带， ICHARG 注释掉，如设置为11，半导体可能算出金属的结果来。 用常规的DFT的波函数来计算hse06的能带，虽然结果一样，但计算时间会多3倍 vasp.5.4.1，the speed would be accelerated greatly via the following parallel setting. NPAR=24 KPAR=8 warning ----------------------------------------------------------------------------- | | | W W AA RRRRR N N II N N GGGG !!! | | W W A A R R NN N II NN N G G !!! | | W W A A R R N N N II N N N G !!! | | W WW W AAAAAA RRRRR N N N II N N N G GGG ! | | WW WW A A R R N NN II N NN G G | | W W A A R R N N II N N GGGG !!! | | | | ALGO=A and IALGO=5X tend to fail with the tetrahedron method | | (e.g. Bloechls method ISMEAR=-5 is not variational) | | please switch to IMSEAR=0-n, except for DOS calculations | | For DOS calculations use IALGO=53 after preconverging with ISMEAR=0 | | I HOPE YOU KNOW, WHAT YOU ARE DOING | | incar.hse NPAR=24 KPAR=8 PREC = high ENCUT = 500 eV NELMIN= 5 LREAL = F # ALGO = Fast EDIFF = 1E-5 ISMEAR = 0 ISPIN = 1 # MAGMOM = 3*1 SIGMA = 0.05 #hse06 ICHARG = 2 ISTART = 1 NELM = 100 NSW = 0 IBRION = -1 LHFCALC = .TRUE. HFSCREEN = 0.2 ; AEXX = 0.25;PRECFOCK= F ALGO = D TIME = 0.4 LDIAG = .TRUE. LCHARG = F LWAVE = F LPLANE=.TRUE.
主要参考以下博文 http://cms.mpi.univie.ac.at/wiki/index.php/Bandstructure_of_SrVO3_in_GW#The_dielectric_function http://blog.sciencenet.cn/blog-1460330-807089.html step1 ALGO = Exact NBANDS = 64 LOPTICS = .TRUE. NEDOS = 2000 ## you might try #LPEAD = .TRUE. ISMEAR = 0 SIGMA = 0.05 GGA = PE NPAR=4 step2: ## Frequency dependent dielectric tensor including ## local field effects within the RPA (default) or ## including changes in the DFT xc-potential (LRPA=.FALSE.). ## N.B.: beware one first has to have done a ## calculation with ALGO=Exact and LOPTICS=.TRUE. ## and a reasonable number of virtual states (see above) ALGO = GW0 ; LSPECTRAL = .TRUE. ; NOMEGA = 50 #LRPA = .FALSE. ## be sure to take the same number of bands as for ## the LOPTICS=.TRUE. calculation, otherwise the ## WAVEDER file is not read correctly NBANDS = 64 step3 ## Frequency dependent dielectric tensor including ## local field effects within the RPA (default) or ## including changes in the DFT xc-potential (LRPA=.FALSE.). ## N.B.: beware one first has to have done a ## calculation with ALGO=Exact and LOPTICS=.TRUE. ## and a reasonable number of virtual states (see above) #ALGO = GW0 ; LSPECTRAL = .TRUE. ; NOMEGA = 50 #LRPA = .FALSE. ## be sure to take the same number of bands as for ## the LOPTICS=.TRUE. calculation, otherwise the ## WAVEDER file is not read correctly NBANDS = 64 #bse ALGO=BSE NOMEGA=72 ENCUTGW = 500 NBANDSO = 4 NBANDSV = 4 LSPECTRAL = .TRUE. 附相关参数说明：(from http://blog.sciencenet.cn/blog-1460330-807089.html ) GW NOMEGA= integer ！ spcifies the number of frequency grid points NOMEGA= integer ！ spcifies the number of frequency points along real axis default： NOMEGA= 50 NOMEGA= NOMEGA for GW。 NOMEGA的取值一般为50-100之间，且必须为CPU核心数数整数倍， 同时NBANDS也必须为CPU核心数整数倍。 如果不是的话，vasp会随机产生波函数，补成整数倍，这些随机产生的波函数会导致结果出现很大误差。 NEDOS ！number of energy points 可决定介电函数的取值密度 默认NEDOS=2000， 在OUTCAR里介电函数的取值密度是 能量范围除以2000 在vasprun.xml里有10倍精细度的介电函数 NBANDSGW= ！determines how many QP energies are calculated and updated in GW calculations. This value usually needs to be increased somewhat for partially or fully selfconsistent calculations. Very accurate results are only obtained when NBANDSGW approaches NBANDS, although this dramatically increases the computational requirements. BSE ALGO=BSE 用GW算介电常数，一般是结合BSE一起计算的，算完BSE后在vasprun.xml中有一个精细的介电矩阵。 NBANDSO = 4 ! number of bands for electron-hole treatment (occupied) NBANDSV = 6 ! number of bands for electron-hole treatment (virtual) OMEGAMAX=6 ! max frequency yunhailiseu： when vasp reports： ———————————————————————————— “BSE diagnonalizing matrix （****） BSE calculating oscillator strength” ———————————————————————————— the whole calculation is almost over. And if no error message is given after that, the job just finished sucessfully . The BSE oscillator strength and dielectric tensor is written in vasprun.xml, not OUTCAR。 ———————————————————————————— dielectricfunction imag array dimension dim=1gridpoints/dimension fieldenergy/field fieldxx/field fieldyy/field fieldzz/field fieldxy/field fieldyz/field fieldzx/field set r 0.0000 0.2951 0.2951 0.0068 0.0000 0.0000 0.0000 /r r 0.0120 0.2952 0.2952 0.0068 0.0000 0.0000 0.0000 /r r 0.0239 0.2952 0.2952 0.0068 0.0000 0.0000 0.0000 /r r 0.0359 0.2953 0.2953 0.0068 0.0000 0.0000 0.0000 /r ———————————————————————————— Then in the following BSE calculation the excitation energies can be determined ( in vasprun.xml just before the dielectric tensor, together with transition matrix elements ). If one exciton has a corresponding excitation energy lower than Eg, it is a bound exciton（束缚激子）. And if the excitation energy is larger than Eg, it is a resonant exciton（共振激子）. The type of exciton may also be determined from absorption spectrum. Nano Lett. 2010, 10, 426-431, Nano Lett. 2007, 7, 3112-3115 and PRB 83, 085405 (2011).
HSE+GW+BSE几类计算的关键词小结 已有 394 次阅读 2014-6-27 15:52 | 系统分类: 科研笔记 GW NOMEGA= integer ！ spcifies the number of frequency grid points NOMEGA= integer ！ spcifies the number of frequency points along real axis default： NOMEGA= 50 NOMEGA= NOMEGA for GW。 NOMEGA的取值一般为50-100之间，且必须为CPU核心数数整数倍， 同时NBANDS也必须为CPU核心数整数倍。 如果不是的话，vasp会随机产生波函数，补成整数倍，这些随机产生的波函数会导致结果出现很大误差。 NEDOS ！number of energy points 可决定介电函数的取值密度 默认NEDOS=2000， 在OUTCAR里介电函数的取值密度是 能量范围除以2000 在vasprun.xml里有10倍精细度的介电函数 NBANDSGW= ！determines how many QP energies are calculated and updated in GW calculations. This value usually needs to be increased somewhat for partially or fully selfconsistent calculations. Very accurate results are only obtained when NBANDSGW approaches NBANDS, although this dramatically increases the computational requirements. BSE ALGO=BSE 用GW算介电常数，一般是结合BSE一起计算的，算完BSE后在vasprun.xml中有一个精细的介电矩阵。 NBANDSO = 4 ! number of bands for electron-hole treatment (occupied) NBANDSV = 6 ! number of bands for electron-hole treatment (virtual) OMEGAMAX=6 ! max frequency yunhailiseu： when vasp reports： ———————————————————————————— “BSE diagnonalizing matrix （****） BSE calculating oscillator strength” ———————————————————————————— the whole calculation is almost over. And if no error message is given after that, the job just finished sucessfully . The BSE oscillator strength and dielectric tensor is written in vasprun.xml, not OUTCAR。 ———————————————————————————— dielectricfunction imag array dimension dim=1gridpoints/dimension fieldenergy/field fieldxx/field fieldyy/field fieldzz/field fieldxy/field fieldyz/field fieldzx/field set r 0.0000 0.2951 0.2951 0.0068 0.0000 0.0000 0.0000 /r r 0.0120 0.2952 0.2952 0.0068 0.0000 0.0000 0.0000 /r r 0.0239 0.2952 0.2952 0.0068 0.0000 0.0000 0.0000 /r r 0.0359 0.2953 0.2953 0.0068 0.0000 0.0000 0.0000 /r ———————————————————————————— Then in the following BSE calculation the excitation energies can be determined ( in vasprun.xml just before the dielectric tensor, together with transition matrix elements ). If one exciton has a corresponding excitation energy lower than Eg, it is a bound exciton（束缚激子）. And if the excitation energy is larger than Eg, it is a resonant exciton（共振激子）. The type of exciton may also be determined from absorption spectrum. Nano Lett. 2010, 10, 426-431, Nano Lett. 2007, 7, 3112-3115 and PRB 83, 085405 (2011).
vasp运行中常见错误的解决 已有 725 次阅读 2014-6-25 10:12 | 系统分类: 科研笔记 1. forrtl: severe (174): SIGSEGV, segmentation fault occurred 在FFLAGS 选项的赋值中增加 -heap-arrays 64 运行的时候加 ulimit -s unlimited 再make得到vasp,运行即可 Point 2 运行vasp的计算速度慢 对openmpi试试用export OMP_NUM_THREADS=1 或者换用intel-mpi 2. VERY BAD NEWS! internal error in subroutine INVGRP: inverse of rotation matrix was not found (increase SYMPREC) 一，取消对称性 ISYM=0 二，减小对称性判断标准 SYMPREC默认是10的负4次，你把它改得稍微小点。如：有5E-4 或 1E-3 3. 计算dos时 ismear=-5 提示 VERY BAD NEWS! internal error in subroutine IBZKPT: Routine TETIRR needs special values for the k-mesh shifts! -3 1) for slab calculations, use ONE k-point along z only. 2) the tetrahedron integration method needs the Gamma point to be included for a subdivision of the BZ into tertahedra, probably the even mesh did not have the Gamma point in the sample. Please try 11x11x1 or 9x9x1 or give 'Gamma' or 'gamma' in the third line of your KPOINTS file 4. When calculate the Macroscopic dielectric properties and Born effective charge tensors(LCALCEPS=.TRUE.),errors: internal error in GENERATE_KPOINTS_TRANS: number of G-vector changed in star 3199 3197 internal error in GENERATE_KPOINTS_TRANS: number of G-vector changed in star 3199 3197 internal error in GENERATE_KPOINTS_TRANS: number of G-vector changed in star 3199 3197 Usually it is because the NPAR is not set proper, remove NPAR from INCAR. 5. EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 http://cms.mpi.univie.ac.at/vasp-forum/forum_viewtopic.php?3.12829 http://cms.mpi.univie.ac.at/vasp-forum/forum_viewtopic.php?3.12158 First: If the vasp was compiled by intel compiler 2013, then change to 2012 or lower; Seceond: change the the optimization flag from -O3 to -O1 . 6. VASP在电子自洽计算的中间步骤中会出现如下的错误 WARNING: DENTET: can't reach specified precision Number of Electrons is NELECT = 196.0137087990377 RMM: 7 -0.461353114525E+03 0.15540E+03 -0.29356E+02 6562 0.456E+01BRMI X: very serious problems the old and the new charge density differ old charge density: 195.99999 new 196.01370 0.758E+01 RMM: 8 -0.228026134405E+03 0.23333E+03 -0.10404E+02 4963 0.286E+01BRMI X: very serious problems the old and the new charge density differ old charge density: 196.01370 new 195.99999 0.376E+01 出现此警告(DENTET)的原因是因为无法通过tetrahedron方法得到足够精确的费米能级。也就是将态密度积分到费米面的电子数和体系的价电子数目不一致。可以尝试采用 选择另一种布里渊区内的积分方法 ( 改变 ISMEAR) 以解决此问题。 VASP计算中Sub-Space-Matrix is not hermitian in DAV的错误 我在计算界面体系时候，其他计算条件不变，仅改变了一些k格点数，就一直提示如下的错误： DAV: 13 -0.242323773333E+03 0.98155E+02 -0.87140E+01 48832 0.949E+01BRMIX: very serious problems the old and the new charge density differ old charge density: 252.00012 new 252.29979 0.809E+01 DAV: 14 -0.392866843695E+03 -0.15054E+03 -0.76122E+01 50857 0.731E+01BRMIX: very serious problems the old and the new charge density differ old charge density: 252.29979 new 252.48257 0.484E+01 WARNING: Sub-Space-Matrix is not hermitian in DAV 9 0.133520549894753 WARNING: Sub-Space-Matrix is not hermitian in DAV 17 495.153990161108 WARNING: Sub-Space-Matrix is not hermitian in DAV 6 0.250235927490523 WARNING: Sub-Space-Matrix is not hermitian in DAV 9 1876.75162244581 解决办法只需调整 AMIX, BMIX的值，把他们设置小一些。
1. forrtl: severe (174): SIGSEGV, segmentation fault occurred 在FFLAGS 选项的赋值中增加 -heap-arrays 64 运行的时候加 ulimit -s unlimited 再make得到vasp,运行即可 Point 2 运行vasp的计算速度慢 对openmpi试试用export OMP_NUM_THREADS=1 或者换用intel-mpi 2. VERY BAD NEWS! internal error in subroutine INVGRP: inverse of rotation matrix was not found (increase SYMPREC) 一，取消对称性 ISYM=0 二，减小对称性判断标准 SYMPREC默认是10的负4次，你把它改得稍微小点。如：有5E-4 或 1E-3 3. 计算dos时 ismear=-5 提示 VERY BAD NEWS! internal error in subroutine IBZKPT: Routine TETIRR needs special values for the k-mesh shifts! -3 1) for slab calculations, use ONE k-point along z only. 2) the tetrahedron integration method needs the Gamma point to be included for a subdivision of the BZ into tertahedra, probably the even mesh did not have the Gamma point in the sample. Please try 11x11x1 or 9x9x1 or give 'Gamma' or 'gamma' in the third line of your KPOINTS file 4. When calculate the Macroscopic dielectric properties and Born effective charge tensors(LCALCEPS=.TRUE.),errors: internal error in GENERATE_KPOINTS_TRANS: number of G-vector changed in star 3199 3197 internal error in GENERATE_KPOINTS_TRANS: number of G-vector changed in star 3199 3197 internal error in GENERATE_KPOINTS_TRANS: number of G-vector changed in star 3199 3197 Usually it is because the NPAR is not set proper, remove NPAR from INCAR. 5. EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 EDWAV: internal error, the gradient is not orthogonal 1 1 http://cms.mpi.univie.ac.at/vasp-forum/forum_viewtopic.php?3.12829 http://cms.mpi.univie.ac.at/vasp-forum/forum_viewtopic.php?3.12158 First: If the vasp was compiled by intel compiler 2013, then change to 2012 or lower; Seceond: change the the optimization flag from -O3 to -O1 . 6. VASP在电子自洽计算的中间步骤中会出现如下的错误 WARNING: DENTET: can't reach specified precision Number of Electrons is NELECT = 196.0137087990377 RMM: 7 -0.461353114525E+03 0.15540E+03 -0.29356E+02 6562 0.456E+01BRMI X: very serious problems the old and the new charge density differ old charge density: 195.99999 new 196.01370 0.758E+01 RMM: 8 -0.228026134405E+03 0.23333E+03 -0.10404E+02 4963 0.286E+01BRMI X: very serious problems the old and the new charge density differ old charge density: 196.01370 new 195.99999 0.376E+01 出现此警告(DENTET)的原因是因为无法通过tetrahedron方法得到足够精确的费米能级。也就是将态密度积分到费米面的电子数和体系的价电子数目不一致。可以尝试采用 选择另一种布里渊区内的积分方法 ( 改变 ISMEAR) 以解决此问题。 VASP计算中Sub-Space-Matrix is not hermitian in DAV的错误 我在计算界面体系时候，其他计算条件不变，仅改变了一些k格点数，就一直提示如下的错误： DAV: 13 -0.242323773333E+03 0.98155E+02 -0.87140E+01 48832 0.949E+01BRMIX: very serious problems the old and the new charge density differ old charge density: 252.00012 new 252.29979 0.809E+01 DAV: 14 -0.392866843695E+03 -0.15054E+03 -0.76122E+01 50857 0.731E+01BRMIX: very serious problems the old and the new charge density differ old charge density: 252.29979 new 252.48257 0.484E+01 WARNING: Sub-Space-Matrix is not hermitian in DAV 9 0.133520549894753 WARNING: Sub-Space-Matrix is not hermitian in DAV 17 495.153990161108 WARNING: Sub-Space-Matrix is not hermitian in DAV 6 0.250235927490523 WARNING: Sub-Space-Matrix is not hermitian in DAV 9 1876.75162244581 解决办法只需调整 AMIX, BMIX的值，把他们设置小一些。
单机4核 32位vasp并行安装 采用IFC编译器，MKL数据库mpich2-1.0.8对VASP编译的过程 1． 准备 系统为 suse linux enterprise Desktop service Pack2 For x86 VASP源代码（vasp.4.6.tar.gz和vasp.4.lib.tar.gz）,mkl数据库（l_mkl_p_9.1.023.tar),ifc编译器（l_fc__pl_9.1.036.tar.gz），mpich2-1.0.8。我们将以上安装所需文件都放在/root/vasp目录下并解压。测试成功效率在95%以上 2． Ifc编译器安装 先解压tar –zxvf l_fc_c_9.1.036.tar.gz 得到l_fc_c_9.1.036.文件夹 进入l_fc__pl_9.1.036文件夹找到install.sh文件 执行./install.sh开始安装ifc，安装过程都选用默认的路径(/opt/intel/fc/9.1.036)安装 安装完毕，进入/opt/intel/fc/9.1.036/bin目录 执行cp ifort /bin ifc安装完成 3． mkl的安装 进入l_mkl_p_9.1.023.tar所在目录 tar –zxvf l_mkl_p_9.1.023.tar 进入解压得到的目录l_mkl_p_9.1.023 进入install文件夹，可以看到一个可执行的文件install ./install.sh 默认安装即可，默认目录为/opt/intel/mkl/9.1.023 三 ，设置环境变量 编辑root下的 .bashrc # .bashrc # User specific aliases and functions alias rm='rm -i' alias cp='cp -i' alias mv='mv -i' PATH=$PATH:/usr/local/bin （安装mpich2要添加的路径） # Source global definitions if ; then . /etc/bashrc fi export LD_LIBRARY_PATH=/opt/intel/mkl/9.1.023/lib/32:/opt/intel/fc/9.1.036/lib . /opt/intel/fc/9.1.036/bin/ifortvars.sh （新添加路径） export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/opt/intel/mkl/9.1.023/lib/32 （新添加路径） 二、安装MPICH2(在节点root目录下) 1、解压缩 #tar -zxvf mpich2-1.0.1.tar.gz 或者 #gunzip -c mpich2-1.0.1.tar.gz|tar xf mpich2-1.0.1.tar 2、创建安装目录 #mkdir /usr/MPICH-instsll 3、进入mpich2解压目录 #cd mpich2-1.0.1 4、默认安装目录 #./configure 5、编译 #make 6、安装 #make install 7、退出到root目录 #cd .. 8、通过编辑.bashrc文件修改环境变量 #vi .bashrc 修改后的.bashrc文件如下： # .bashrc # User specific aliases and functions alias rm='rm -i' alias cp='cp -i' alias mv='mv -i' PATH=$PATH:/usr/local/bin （新增加的mpich2的路径） #Source global definitions if ; then . /etc/bashrc fi export LD_LIBRARY_PATH=/opt/intel/mkl/9.1.023/lib/32:/opt/intel/fc/9.1.036/lib . /opt/intel/fc/9.1.036/bin/ifortvars.sh （新添加路径） export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/opt/intel/mkl/9.1.023/lib/32 （新添加路径） 9、测试环境变量设置 #which mpd #which mpicc #which mpiexec #which mpirun 10、修改/etc/mpd.conf文件，添加内容为secretword=myword #vi /etc/mpd.conf 添加 secretword=myword 设置文件读取权限和修改时间 #touch /etc/mpd.conf #chmod 600 /etc/mpd.conf 1、本地测试 #mpd 启动 #mpdtrace 观看启动机器 #mpdallexit 退出 三。 6．进入vasp.4.lib所在目录 cd /root/vasp/vasp.4.lib cp makefile.linux_ifc_P4 makefile 根据自己机子的情况选择合适的makefile 编辑makfile 19行的 FC＝ifc 修改为 FC＝ifort make 如果编译通过，说明前面安装的数学库和编译器等都是正确的 7．进入vasp.4.6所在目录 cd /root/vasp/vasp.4.6 cp makefile.linux_ifc_P4 makefile 编辑 makefile 对makefile做如下修改： 50行， 52行 前加 ＃ 80行前加 ＃ 128行 的 #BLAS=-L/opt/intel/mkl/lib/32 -lmkl_p4 -lpthread 修改为 BLAS=-L/opt/intel/mkl/9.1.023/lib/32 -lmkl_p4 -lsvml -lvml -lguide -lpthread 136行前加 # 145 行 的 #LAPACK= -lmkl_lapack 修改為LAPACK= -lmkl_lapack 或 LAPACK=-L/opt/intel/mkl/9.1.023/lib/32 -lmkl_lapack -lsvml -lvml -lguide –lpthread 也可以 149行 前 加 # 166行 前 加 # 201行 的 #FC=mpif77 修改為 FC=mpif90 202行的 #FCL=$(FC) 修改為 FCL=$(FC) 211-214行 前面的 # 去掉 修改為 CPP = $(CPP_) -DMPI -DHOST=\LinuxIFC\ -DIFC \ -Dkind8 -DNGZhalf -DCACHE_SIZE=4000 -DPGF90 -Davoidalloc \ -DMPI_BLOCK=500 \ -DRPROMU_DGEMV -DRACCMU_DGEMV 224行 前 加 # 227行 前 加 # 233-235行 前 的 # 去掉 238行 前 的 # 去掉 343行 的 -e95 去掉 保存退出后 编译make 编译通过则 cp vasp /bin 在任何目录下vasp命令都可以调用了 这时编译的vasp只能在root用户下使用 四．重启电脑，以普通用户(XXX代表普通用户名)的身份登陆 1. 进入 /home/XXX目录找到隐藏的 .bashrc文件（在窗口的工具栏中点击“查看”,选择显示隐藏文件） .bashrc 文件如下 # Sample .bashrc for SuSE Linux # Copyright (c) SuSE GmbH Nuernberg # There are 3 different types of shells in bash: the login shell, normal shell # and interactive shell. Login shells read ~/.profile and interactive shells # read ~/.bashrc; in our setup, /etc/profile sources ~/.bashrc - thus all # settings made here will also take effect in a login shell. # # NOTE: It is recommended to make language settings in ~/.profile rather than # here, since multilingual X sessions would not work properly if LANG is over- # ridden in every subshell. # Some applications read the EDITOR variable to determine your favourite text # editor. So uncomment the line below and enter the editor of your choice :-) #export EDITOR=/usr/bin/vim #export EDITOR=/usr/bin/mcedit # For some news readers it makes sense to specify the NEWSSERVER variable here #export NEWSSERVER=your.news.server # If you want to use a Palm device with Linux, uncomment the two lines below. # For some (older) Palm Pilots, you might need to set a lower baud rate # e.g. 57600 or 38400; lowest is 9600 (very slow!) # #export PILOTPORT=/dev/pilot #export PILOTRATE=115200 test -s ~/.alias . ~/.alias || true 修改后变为 # Sample .bashrc for SuSE Linux # Copyright (c) SuSE GmbH Nuernberg # There are 3 different types of shells in bash: the login shell, normal shell # and interactive shell. Login shells read ~/.profile and interactive shells # read ~/.bashrc; in our setup, /etc/profile sources ~/.bashrc - thus all # settings made here will also take effect in a login shell. # # NOTE: It is recommended to make language settings in ~/.profile rather than # here, since multilingual X sessions would not work properly if LANG is over- # ridden in every subshell. # Some applications read the EDITOR variable to determine your favourite text # editor. So uncomment the line below and enter the editor of your choice :-) #export EDITOR=/usr/bin/vim #export EDITOR=/usr/bin/mcedit # For some news readers it makes sense to specify the NEWSSERVER variable here #export NEWSSERVER=your.news.server # If you want to use a Palm device with Linux, uncomment the two lines below. # For some (older) Palm Pilots, you might need to set a lower baud rate # e.g. 57600 or 38400; lowest is 9600 (very slow!) # #export PILOTPORT=/dev/pilot #export PILOTRATE=115200 test -s ~/.alias . ~/.alias || true 以下为增加内容 export LD_LIBRARY_PATH=/opt/intel/mkl/9.1.023/lib/32:/opt/intel/fc/9.1.036/lib . /opt/intel/fc/9.1.036/bin/ifortvars.sh （新添加路径） export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/opt/intel/mkl/9.1.023/lib/32 （新添加路径） 保存， 2. 配置mpd 在普通用户下 执行 cd $HOME vi .mpd.conf 添加 MPD_SECRETWORD=mr45-j9z 保存， 设置文件读取权限和修改时间 touch .mpd.conf chmod 600 .mpd.conf 设置完成。 上述编译通过的vasp还不能运行，会提示找不到目标文件libsvml.so 1. locate libsvml.so可以找到libsvml.so所在目录/opt/intel_fc_80/lib中 2. 在root目录下找到.bash_profile和.bashrc文件（ls –a 显示当前目录下所有内容 3．在root目录下运行env看环境变量LD_LIBRARY_PATH是否为 LD_LIBRARY_PATH=/opt/intel/mkl/8.0.2/lib/32: /opt/intel/fce/9.1.036/lib 若不是则执行如下命令 export LD_LIBRARY_PATH=/opt/intel/mkl/8.0.2/lib/32: /opt/intel/fce/9.1.036/lib 4．reboot后即可以使用vasp 计算了
INCAR文件中LREAL参数 ----------------------------------------------------------------------------- | | | ADVICE TO THIS USER RUNNING 'VASP/VAMP' (HEAR YOUR MASTER'S VOICE ...): | | | | You have a (more or less) 'large supercell' and for larger cells | | it might be more efficient to use real space projection operators | | So try LREAL= Auto in the INCAR file. | | Mind: If you want to do a very accurate calculations keep the | | reciprocal projection scheme (i.e. LREAL=.FALSE.) | | | ----------------------------------------------------------------------------- INCAR文件中LREAL参数是什么意思，应该如何设置。 现将这个参数和不同选择项意思列于下，供参考。 LREAL: Default= .FALSE. LREAL是赝势的非局域部分用到的一个积分，其在倒格空间或者实空间都可以求值。这个 选项就是决定是在哪个空间里求。在倒格空间里，采用平面波基组求解，在实空间里，采 用积 分球求解。缺省是.FALSE，即不在实空间求。但效率会低一些。 其他选项是 O or On,A or Auto 和.True.。 On和.TRUE.的差别在于是否使用King-Smith算法优化，Auto则自动选择，推荐。 如在OUTCAR中出现如下字样 ----------------------------------------------------------------------------- | | | ADVICE TO THIS USER RUNNING 'VASP/VAMP' (HEAR YOUR MASTER'S VOICE ...): | | | | You havea (more or less) 'large supercell' and for larger cells | | it mightbe more efficient to use real space projection operators | | So tryLREAL= Auto in the INCAR file. | | Mind: Ifyou want to do an extremely accurate calculations keep the | | reciprocal projection scheme (i.e. LREAL=.FALSE.) | | | ----------------------------------------------------------------------------- 表明，选择LREAL选择为在倒格空间求解，这样的求解积分的方法对于大的晶包体系，会 影响计算效率，但不影响计算进度。可以根据需求选择LREAL参数项。 LREAL-tag (and ROPT-tag)LREAL = .TRUE. | .FALSE. ROPT = Default LREAL = .FALSE. .FALSE. projection done in reciprocal space .TRUE. projection done in real space, (old, superseded by LREAL=O) On or O projection done in real space, projection operators are re-optimized Auto or A projection done in real space, fully automatic optimization of projection operators no user interference required Determines whether the projection operators are evaluated in real-space or in reciprocal space:The non local part of the pseudopotential requires the evaluation of an expression .The ``projected wavefunction character'' is defined as: This expression can be evaluated in reciprocal or real space:In reciprocal space (second line)the number of operations scales with the size of the basis set (i.e.number of plane-waves). In real space (first line) theprojection-operators are confined to spheres around each atom. Therefore the number of operations necessary to evaluate one does not increase with the system size (usually the number of grid points within the cut-off-sphere is between 500 and 2000). One of the major obstacles of the method working in real space is that the projection operatorsmust be optimized, i.e. all high frequency components must be removed fromthe projection operators. If this is not done 'aliasing' can happen(i.e. the high frequency components of the projection operators are aliasedto low frequency components and a random noise is introduced). Currently VASP supports three different schemes to remove the high frequency components from the projectors. LREAL = .TRUE. is the simplest one. If LREAL = .TRUE. is selected,the real space projectors which have been generated by the pseudopotential generationcode are used. This requires no user interference.For LREAL = On the real space projectors are optimized by VASP using an algorithm proposed by King-Smith et al. . For LREAL = Auto a new scheme is used which isconsiderably better (resulting in more localized) projector functionsthan the King-Smith et al. method. To fine-tune the optimization procedure the flag ROPT can be usedif LREAL = Auto or LREAL = On is used. We recommend to use the real-space projection scheme for systems containingmore than 20 atoms. We also recommend to use only LREAL = Auto (for version VASP.4.4 and newer releases) and LREAL = On (for all other versions).Version 4.4 also supports the old mode LREAL= O toallow calculations that are fully compatible to VASP.4.3 (and VASP.3.2).The best performance is generally achieved with LREAL= Auto, butif performance is not that important you can also use LREAL=.TRUE. whichgenerally requires less user interference. You can skip therest of the paragraph, if you use only LREAL=.TRUE.. For LREAL = O and LREAL = Athe projection operators are optimized by VASPon the fly (i.e. on startup). Several flags influence the optimization ENCUT (i.e. the energy cutoff), components beyond the energy cutoff are 'removed' from the projection operators. PREC tag specifies how precise the real space projectorsshould be, and sets the variables ROPT accordingly to the followingvalues: For LREAL = On PREC = Low 700 points in the real space sphere (ROPT=0.67) PREC = Med 1000 points in the real space sphere (ROPT=1.0) PREC = High 1500 points in the real space sphere (ROPT=1.5) For LREAL = Auto PREC = Low accuracy (ROPT=0.01) PREC = Med accuracy 2 (ROPT=0.002) PREC = High accuracy 2 (ROPT=2E-4) These defaults can be superseded by the line ROPT = one_number_for_each_species in the INCAR file. For instance ROPT = 0.7 1.5 will set the number of real space points within the cutoff sphere for the first speciesto approximately 700, and that for the second species to 1500.In VASP.4.4 alternatively the ``precision'' of the operators can be specifiedwriting i.e.ROPT = 1E-3 1E-3 In that case the real space operators will be optimized foran accuracy of approximately 1meV/atom ( ). The ``precision'' mode works both for LREAL=On and LREAL=Auto (but to maintaincompatibility with older VASP versions it is only selected if LREAL = Auto is specified in the INCAR file).The precision mode is generally switched on if the value for ROPTis smaller than 0.1. The ``precision'' mode and the conventional mode can be intermixed,i.e. it is possible to specifyROPT = 0.7 1E-3 in that case the number of real space points within the cutoff sphere for the first specieswill be approximately 700, whereas the real space projector functions for the secondspecies are optimized for an accuracy of approximately 1 meV.We recommend to use the ``precision'' mode with a target accuracy of around eV/atom if your version supports this. If you use the mode in which the number of grid points in the real space projection sphere is specified, you have to selectROPT carefully, especially if a hard species is mixed with a softspecies. In that case the following lines in the OUTCAR file mustbe checked (here is the output for LREAL = On, but that one forLREAL = Auto is quite similar ) Optimization of the real space projectors maximal supplied Q-value = 12.85 optimization between = = Ry Optimized for a Real-space Cutoff 2.30 Angstroem l X(QCUT) X(cont) X(QGAM) max X(q) W(q)/X(q) e(spline) 0 9.518 9.484 -.004 18.582 .11E-03 .16E-06 0 -2.149 -2.145 .001 3.059 .17E-03 .25E-06 1 8.957 8.942 .003 9.950 .14E-03 .34E-06 1 1.870 1.870 .001 1.837 .95E-03 .51E-06 2 3.874 3.866 .000 4.764 .15E-03 .68E-07 The meaning of QCUT and QGAM is explained in Sec. 11.5.6 . The most important information is given in thecolumn W(q)/X(q) (respectively the column W(low)/X(q) forLREAL = Auto).The values in these columns must be as small as possible.If these values are too large, increase the ROPT tag from the default value.As a rule of thumb the maximum allowed value in this column is for PREC = Med.(For PREC = Low errors might be around and for PREC = High errorsshould be smaller than ). If W(q)/X(q) is larger than theerrors introduced by the real space projections can be substantial.In this case ROPT must be specified in the INCAR file to avoidincorrect results. If the new precision mode is used in VASP.4.4 (ROPT 0.1) thecode automatically selects the real-space cutoff so that therequired precision is reached. A few comments for non-experts and experts:Real space optimization (LREAL = .TRUE., LREAL = On or LREAL = Auto)always results in a small (not necessarily negligible)error (the error is usually a constant energy shift for each atom). If you are interested in energy differences of a fewmeV use only calculations with the same setup (i.e. same ENCUT,PREC, LREAL and ROPT setting) forall calculations.For example, if you want to calculatesurface energies recalculate the bulk groundstate energy withexactly the same setting you are going to use for the surface. Another possibility is to relax the surface with real space projection, and to do one final total energy calculation with LREAL = .FALSE. to get exact energies. Anyway, for PREC = Med, the errors introduced by the real spaceprojection are usually of the same order magnitude as those introduced by the wrap around errors.For PREC = High errors are usually less than meV. PREC = Lowshould be used only for high speed MD's,if computer resources are really a problem.
HSE hybrid functional：Hartree-Fock (HF) type and hybrid functional calculations Available only in VASP.5.X. http://cms.mpi.univie.ac.at/vasp/vasp/Hartree_Fock_HF_type_hybrid_functional_calculations.html INCAR FOR HSE # output options LWAVE = .FALSE. # write or don't write WAVECAR LCHARG = .FALSE. # write or don't write CHG and CHGCAR LELF = .FALSE. # write ELF # ionic relaxation NSW = 100 # number of ionic steps IBRION = 1 # 2=conjucate gradient, 1=Newton like ISIF = 3 # 3=relax everything, 2=relax ions only, 4=keep volume fixed # precision parameters EDIFF = 1E-7 # 1E-3 very low precision for pre-relaxation, use 1E-5 next EDIFFG = -1E-3 # usually: 10 * EDIFF PREC = high # precision low, med, high, accurate # electronic relaxation ISMEAR = 0 # -5 = tetraedon, 1..N = Methfessel SIGMA = 0.1 ENCUT = 600 # cutoff energy PSTRESS = 0 #ISYM=0 # Choose DFT functional - HSE06 ISTART = 1 LHFCALC = .TRUE . ; HFSCREEN = 0.2 NBANDS = 16 ALGO = All ; TIME = 0.4 PRECFOCK = Fast ! used PRECFOCK = Normal for high quality calculations #NKRED = 2 ! omit flag for high quality calculations HSE hybrid functional：Hartree-Fock (HF) type and hybrid functional calculations Available only in VASP.5.X. http://cms.mpi.univie.ac.at/vasp/vasp/Hartree_Fock_HF_type_hybrid_functional_calculations.html Subsections Introduction: Hartree-Fock LHFCALC-tag Amount of exact/DFT exchange and correlation: AEXX, AGGAX, AGGAC and ALDAC tags ENCUTFOCK: FFT grid in the Hartree-Fock related routines PRECFOCK: FFT grid in the Hartree-Fock and GW related routines LMAXFOCK (or old HFLMAXF ) LMAXFOCKAE HFSCREEN and LTHOMAS NKRED, NKREDX, NKREDY, NKREDZ and EVENONLY, ODDONLY When NKRED should not be applied Typical hybrid functional and Hartree-Fock calculations (1) Typical hybrid functional and Hartree-Fock calculations It is strongly recommended to perform standard DFT calculations first , and to start Hartree-Fock type calculations from a preconverged WAVECAR file . (a) A typical INCAR file for a Hartree-Fock or hybrid HF/DFT calculation for an insulator or semiconductor has the following input lines: ISTART = 1 LHFCALC = .TRUE. ; HFSCREEN = 0.2 NBANDS = number of occupied bands ALGO = All ; TIME = 0.4 PRECFOCK = Fast #! used PRECFOCK = Normal for high quality calculations NKRED = 2 #! omit flag for high quality calculationsFor (b) For metals and small gap semiconductors it is recommended to use . ISTART = 1 LHFCALC = .TRUE. ; HFSCREEN = 0.2 ALGO = Damped ; TIME = 0.4 PRECFOCK = Fast ! used PRECFOCK = Normal for high quality calculations NKRED = 2 ! omit flag for high quality calculations These input files select the HSE06 functional, which tends to yield very similar thermochemistry as the PBE0 functional, but converges more rapidly with respect to the number of k-points . We thus recommend to apply and use this functional instead of the more demanding PBE0 functional. The NKRED flag is applicable, if and only if the number of k-points is dividable by NKRED (see Sec. 6.71.9 ). PRECFOCK= fast selects a smaller FFT grid for the fast-Fourier-transforms (see Sec. 6.71.5 ). For high accuracy NKRED and in particular PRECFOCK= fast should be ommited, but we recommend to do this only after preconverging the orbitals and atomic positions with the flags specified above. Mind, that the parameter TIME defaults to 0.4, and for the present algorithm this hardly ever needs to be changed. If divergence is observed, simply decrease TIME until the damped or conjugate gradient algorithm become stable (see Sec. 6.47 and 6.51 ). Standard Hartree-Fock type calculations require one to set the flag AEXX = 1.0 to switch on full non-local exchange (local exchange and correlation are automatically switched off): ISTART = 1 LHFCALC = .TRUE. ; AEXX = 1.0 ; NBANDS = number of occupied bands ALGO = All ; TIME = 0.4 PRECFOCK = Fast ! used PRECFOCK = Normal for high quality calculations NKRED = 2 ! omit flag for high quality calculations Concerning NKRED and PRECFOCK the same considerations as above apply. Matter of fact, it is also possible to try to converge using the ``metallic'' setup given above. (2)Notes （1）LHFCALC-tag LHFCALC= .TRUE. | .FALSE. Default: LHFCALC=.FALSE. The flag specifies, whether Hartree-Fock type calculations are performed. At the moment, it is recommended to select an all bands simultaneous algorithm , i.e. ALGO =Damped ( IALGO =53) or ALGO =All ( IALGO =58) in the INCAR file (see Sec. 6.46 6.47 ). The blocked Davidson algorithm ALGO =Normal is, with certain caveat, also supported, whereas calculations for the other algorithms (ALGO=Fast) are not currently supported (note: no warning is printed). The blocked Davidson algorithm ALGO=Normal is generally rather slow, and in many cases the Pulay mixer will be unable to determine the proper ground-state. We hence recommend to select the blocked Davidson algorithm only in combination with straight mixing or a Kerker like mixing. The following combination have been successfully applied for small and medium sized systems LHFCALC = .TRUE. ; ALGO = Normal ; IMIX = 1 ; AMIX = a Decrease the parameter a until convergence is reached. In most cases, however, it is recommended to use the damped algorithm with suitably chosen timestep . The following setup for the electronic optimization works reliably in most cases: LHFCALC = .TRUE. ; ALGO = Damped ; TIME = 0.4 If convergence is not obtained, it is recommended to reduce the timestep TIME . （2） HFSCREEN HFSCREEN determines the range separation parameter in range separated hybrid functionals. In combination with PBE potentials, attributing a value to HFSCREEN will switch from the PBE0 functional (in case LHFCALC=.TRUE.) to the closely related HSE03 or HSE06 functional . Note: A comprehensive study of the performance of the HSE03/HSE06 functional compared to the PBE and PBE0 functionals can be found in Ref. . The B3LYP functional was investigated in Ref. . Further applications of hybrid functionals to selected materials can be found in the following references: Ceria (Ref. ), lead chalcogenides (Ref. ), CO adsorption on metals (Refs. ), defects in ZnO (Ref. ), excitonic properties (Ref. ), SrTiO and BaTiO (Ref. ). LTHOMAS= .TRUE. | .FALSE. Default: LTHOMAS=.FALSE. If the flag LTHOMAS is set, a similar decomposition of the exchange functional into a long range and a short range part is used. This time, it is more convenient to write the decomposition in reciprocal space: (6.69) where is the Thomas-Fermi screening length. HFSCREEN is used to specify the parameter . For typical semi-conductors, the Thomas-Fermi screening length is about 1.8 Å , and setting HFSCREEN to this value yields reasonable band gap s for most materials. In principle, however, the Thomas-Fermi screening length depends on the valence electron density; VASP determines this parameter from the number of valence electrons (POTCAR) and the volume and writes the corresponding value to the OUTCAR file: Thomas-Fermi vector in A = 2.00000Since, VASP counts the semi-core states and -states as valence electrons, although these states do not contribute to the screening, the values reported by VASP are, however, often incorrect. Details can be found in literature . Another important detail concerns that implementation of the density functional part in the screened exchange case. Literature suggests that a global enhancement factor (see Equ. (3.15) in Ref. ) should be used, whereas VASP implements a local density dependent enhancement factor , where is the Fermi wave vector corresponding to the local density (and not the average density as suggested in Ref. ). The VASP implementation is in the spirit of the local density approximation. （3） ALGO ALGO-tag ALGO = Normal | VeryFast | Fast | Conjugate | All | Damped | Subrot | Eigenval | None | Nothing | Exact | Diag Default ALGO = Normal The ALGO tag is a convenient option to specify the electronic minimisation algorithm in VASP.4.5 and later versions. Except for ``None'' and ``Nothing'', ``Exact'' and ``Diag'' (which must be spelled out), the first letter determines the applied algorithm. Conjugate, Subrot, Eigenval, Exact, None and Nothing are only supported by VASP.5.2.12 and newer versions. ALGO = Normal selects IALGO = 38 (blocked Davidson iteration scheme), whereas ALGO = Very_Fast selects IALGO = 48 (RMM-DIIS). A faily robust mixture of both algorithm is selected for ALGO = Fast. In this case, Davidson (IALGO = 38) is used for the initial phase, and then VASP switches to RMM-DIIS (IALGO = 48). Subsequencly, for each ionic update, one IALGO = 38 sweep is performed for each ionic step (except the first one). The ``all band simultaneous update of orbitals'' can be selected using ALGO = Conjugate or ALGO = All (IALGO = 58, in both cases the same conjugate gradient algorithm is used). A damped velocity friction algorithm is selected using ALGO = Damped (IALGO = 53). ALGO = Subrot selects subspace rotation or diagonalization in the sub-space spanned by the calculated NBANDS orbitals (IALGO = 4). ALGO = Exact or ALGO = Diag performs an exact diagonalization (IALGO = 90), and we recommend to use this if more than 30-50 % of the states are calculated (e.g. for or RPA calculations). ALGO = Eigenval allows to recalculate one electron energies, density of state and perform selected postprocessing using the current orbitals (IALGO = 3) e.g. read from WAVECAR. ALGO = None or ALGO = Nothing allows to recalculate the density of states (eigenvalues from WAVECAR, e.g. using different smearing or tetrahedron method) or perform other selected postprocessing using the current orbitals and one electron energies (IALGO = 2) e.g. read from WAVECAR. （4）PRECFOCK PRECFOCK: FFT grid in the Hartree-Fock and GW related routines PRECFOCK= Low | Medium | Fast | Normal | Accurate Default: PRECFOCK=Normal The PRECFOCK parameter controls the FFT grid for the exact exchange (Hartree-Fock) routines, i.e. it is possible to chose a different grid for the exact exchange part, and for the local Hartree and DFT potentials. In fact, the exchange is rather insensitive to the FFT grids, and in many cases a rather coarse grid can be used to calculate the overlap density and the potentials. Since the exact exchange requires the evaluation of an overlap density (compare 6.59 ) errors in the convolution (aliasing errors) are only avoided, if the FFT grid contains all Fourier components up to twice the plane wave with the largest wave vector ( ). For Low and Fast, however, the smallest possible FFT grid, which just encloses the cutoff sphere ( ) determined by the plane wave cutoff (ENCUT), is used. This accelerates the calculations by roughly a factor two to three, but causes slight changes in the total energies and some noise in the calculated forces. The corresponding FFT grid that is used in the Hartree Fock routines is written to the OUTCAR file after the lines FFT grid for exact exchange (Hartree Fock)For PRECFOCK=Normal, the FFT grid for the exact exchange is identical to the FFT grid used for the orbitals for PREC=Normal in the DFT part. For PRECFOCK=Accurate, the FFT grid for the exact exchange is identical to the FFT grid used for the orbitals for PREC=Accurate in the DFT part (any combination of PREC and PRECFOCK is allowed). For PRECFOCK=Fast, Normal and Accurate, the augmentation charges--which are required to restore the norm and dipoles of the overlap density on the plane wave grid --are made soft, such that an accurate presentation on the plane wave grid is possible even for relatively coarse FFT grids. The sphere size is printed out after Radii for the augmentation spheres in the non-local exchangeThe following table summarises the possible setting: PRECFOCK FFT grid augmentation charge advantage/disatvantage VASP.5.2.2 compatible, not recommended Low identical to standard DFT large noise in forces/energy errors Medium identical to std. FFT identical to standard DFT some noise in forces/good energy VASP.5.2.4 and newer, recommended Fast very soft augmentation charge some noise in forces/good energy Normal 3/2 soft augmentation charge accurate forces and energy Accurate 2 soft augmentation charge very accurate forces and energy soft augmentation charge: radius for augmentation sphere is increased by factor 1.25 compared to default very soft augmentation charge: radius is increased by factor 1.35 compared to default except for like charge, for the channel the radius of the augmentation sphere is increased by a factor 1.25 Even PRECFOCK=Fast yields fairly low noise in the forces and virtually no egg-box effects (aliasing errors). In the forces, the noise is usually below 0.01 eV/Å. For PRECFOCK=N and PRECFOCK=A, noise is usually not an issue, and the accuracy is sufficient even for phonon calculations in large supercells.
VASP Tools Peter’s collection of small, but useful, VASP scripts. Some of these may be found on NSC’s computers by loading the “vasptools” module. vaspcheck The “vaspcheck” script scans your input files and looks for common errors such as: forgetting to copy CONTCAR to POSCAR, misspelled INCAR tags, and various inconsistent configurations such as having N atoms but only N-x magnetic moments specified. If you supply it with information about the intended number of cpu cores and nodes that you will run on, it will also check the NPAR and KPAR values and warn if they are unreasonable. An alpha version is available in the “vasptools” module on Triolith, and once the script has become more robust, I will make a version available for download here. The idea is to integrate it with the queue system in the future, so that your VASP jobs will be automatically checked when running the “sbatch” command. Install: Click here to download to the Python script. Install it by putting it in a folder in your $PATH (such as~/bin) and then do: mv vaspcheck-standalone.py vaspcheckchmod u+x vaspcheck This should allow you to write “vaspcheck” in any job folder. grad2 The gradient program is used to quickly get an overview of VASP geometry optimization runs. It also works for molecular dynamics. It parses the OUTCAR file and writes a one-line summary per ionic step, including average forces and run-time per step. The output looks like this: $ grad2 OUTCARStep 1 Energy: -42.794090 Log|dE|: +1.631 Avg|F|: 1.720176 Max|F|: 2.710237 SCF: 15 Time: 0.52mStep 2 Energy: -43.336931 Log|dE|: -0.265 Avg|F|: 0.513256 Max|F|: 0.696803 SCF: 9 Time: 0.32mStep 3 Energy: -43.384337 Log|dE|: -1.324 Avg|F|: 0.095940 Max|F|: 0.178813 SCF: 7 Time: 0.24mStep 4 Energy: -43.387267 Log|dE|: -2.533 Avg|F|: 0.050399 Max|F|: 0.071963 SCF: 6 Time: 0.20mStep 5 Energy: -43.388250 Log|dE|: -3.007 Avg|F|: 0.043184 Max|F|: 0.061866 SCF: 4 Time: 0.13mStep 6 Energy: -43.390384 Log|dE|: -2.671 Avg|F|: 0.030630 Max|F|: 0.047413 SCF: 5 Time: 0.18mStep 7 Energy: -43.392498 Log|dE|: -2.675 Avg|F|: 0.003760 Max|F|: 0.005849 SCF: 6 Time: 0.20mStep 8 Energy: -43.392509 Log|dE|: -4.959 Avg|F|: 0.001285 Max|F|: 0.001861 SCF: 3 Time: 0.09m The Log\|dE\| is the base-10 logarithm of the absolute value of the energy difference between steps, essentially the number of “stable” decimals in the total energy, i.e. for a VASP run with EDIFF=1.0E-5, Log\|dE\| should approach -5 for convergence. If your terminal has colors, each line will be color-coded with red color signalling convergence problems, and green color corresponding to energy convergence having been reached (as judged by the EDIFF tag in the INCAR file.) Install: Click here to download to the Python script (version 2012-01-07). Install it by putting it in a folder in your $PATH (such as~/bin) and then do: mv grad2.py grad2chmod u+x grad2 This should allow you to write “grad2 OUTCAR” in any job folder. vasp2cif The vasp2cif program creates a CIF-file from a POSCAR/CONTCAR file, which can be used for visualization in graphical applications, such as VESTA. Note that it will not preserve symmetries from VASP – the output CIF is always in P1 symmetry. If there is no POTCAR file available, the script has a flag (“–elements”) whereby you can specify the atomic species. Install: Click here to download to the Python script (version 2012-12-15). Install it by putting it in a folder in your $PATH (such as~/bin) and then do: mv vasp2cif.py vasp2cifchmod u+x vasp2cif This should allow you to write e.g. “vasp2cif CONTCAR” in any job folder. cif2vasp For CIF to VASP input conversion, a recommend a program by Torbjrn Bjrkman called cif2cell , which can generate input from CIF files for many ab initio programs, including VASP. If anyone is interested in doing high-performance CIF to VASP conversion with optimized C code, I can provide source code to a CIF lexer written with the Ragel state machine compiler and CIF parser written with Lemon . Unfortunately, the tool I built to do CIF to VASP with this framework does not compile correctly anymore on recent platforms, so I cannot put up any binaries for download. You can email me at pla@nsc.liu.se if you want to give it a try. supersize Yet another script to make supercells out of VASP’s POSCAR/CONTCAR files. It is basically a wrapper around the corresponding feature in ASE . It works like this: $ supersize POSCAR 4x4x4 And you get a new POSCAR file called “POSCAR.4x4x4” containing the supercell repeated 4 times in a,b,c directions. You can achieve the same effect with three lines of Python code in ASE: import ase.io.vaspcell = ase.io.vasp.read_vasp(POSCAR)ase.io.vasp.write_vasp(POSCAR.4x4x4,cell*(4,4,4), label='444supercell',direct=True,sort=True) The script is a little bit more elaborate, however, since it includes error checks.
B ADER C HARGE A NALYSIS News 07/12/12 - Version 0.28a Released Fixed a sign problem for PWscf cube files. Introduction Richard Bader , from McMaster University, developed an intuitive way of dividing molecules into atoms. His definition of an atom is based purely on the electronic charge density. Bader uses what are called zero flux surfaces to divide atoms. A zero flux surface is a 2-D surface on which the charge density is a minimum perpendicular to the surface. Typically in molecular systems, the charge density reaches a minimum between atoms and this is a natural place to separate atoms from each other. Besides being an intuitive scheme for visualizing atoms in molecules, Bader's definition is often useful for charge analysis. For example, the charge enclosed within the Bader volume is a good approximation to the total electronic charge of an atom. The charge distribution can be used to determine multipole moments of interacting atoms or molecules. Bader's analysis has also been used to define the hardness of atoms, which can be used to quantify the cost of removing charge from an atom. Program Overview We have developed a fast algorithm for doing Bader's analysis on a charge density grid. The program (see below) can read in charge densities in the VASP CHGCAR format, or the Gaussian CUBE format. The program outputs the total charge associated with each atom, and the zero flux surfaces defining the Bader volumes. Note for VASP users Information about generating and analyzing charge densities from vasp . For users using our DOS projection code in Bader volumes , the -vac flag is recommended. Note for CASTEP users Aaron Hopkinson and Dr Matt Probert from the University of York have provided a den2cube.tar.gz utility to convert from the CASTEP charge density to the cube format so that it can be read in by the Bader analysis program. Download Select the appropriate platform to download a binary of the Bader analysis program: Linux x86-64 (ifort) Windows (gfortran) Mac OS X, Intel (ifort) The F90 source code is also available: Source Code (v0.28 06/26/11) Version history Running the Program The program can be run with the command bader chargefile It will automatically determine if the chargefile is a VASP CHGCAR file or a Gaussian CUBE file. The only required input argument is the name of the charge density file. Command line arguments and output files The following options can be used when running the Bader analysis program. bader chargefile To get a description of the options, run 'bader -h'. Output files The following output files are generated: ACF.dat, BCF.dat, AtomVolumes.dat. ACF.dat contains the coordinates of each atom, the charge associated with it according to Bader partitioning, percentage of the whole according to Bader partitioning and the minimum distance to the surface. This distance should be compared to maximum cut-off radius for the core region if pseudo potentials have been used. BCF.dat contains the coordinates of each Bader maxima, the charge within that volume, the nearest atom and the distance to that atom. AtomVolumes.dat contains the number of each volume that has been assigned to each atom. These numbers correspond to the number of the BvAtxxxx.dat files. The Bader volumes can be written using the print options. bader chargefile bader chargefile bader chargefile bader chargefile -p none The default is to write no charge density files. -p all_atom Combine all volumes associated with an atom and write to file. This is done for all atoms and written to files named BvAtxxxx.dat. The volumes associated with atoms are those for which the maximum in charge density within the volume is closest to the atom. -p all_bader Write all Bader volumes (containing charge above threshold of 0.0001) to a file. The charge distribution in each volume is written to a separate file, named Bvolxxxx.dat. It will either be of a CHGCAR format or a CUBE file format, depending on the format of the initial charge density file. These files can be quite large, so this option should be used with caution. -p sel_atom Write the selected atomic volumes, read from the subsequent list or range of volumes. -p sel_bader Write the selected Bader volumes, read from the subsequent list or range of volumes. -p sum_atom Write the sum of selected atomic volumes, read from the subsequent list of volumes. -p sum_bader Write the sum of selected Bader volumes, read from the subsequent list of volumes. -p atom_index Write the atomic volume index to a charge density file. -p bader_index Write the Bader volume index to a charge density file. Visualization The Bader volumes can be written and visualized with the VASP Data Viewer , VMD , or a cube file viewer (such as GaussView) for Gaussian cube files. Examples NaCl crystal (vasp chgcar) NaCl crystal including core charges (vasp chgcar) C 2 H 4 molecule, orientation 1 (vasp chgcar) C 2 H 4 molecule, orientation 2 (vasp chgcar) H 2 O molecule (gaussian cube) Discussion Forum We have a discussion forum to address issues related to the code and running the program. References W. Tang, E. Sanville, and G. Henkelman A grid-based Bader analysis algorithm without lattice bias , J. Phys.: Condens. Matter 21 , 084204 (2009). E. Sanville, S. D. Kenny, R. Smith, and G. Henkelman An improved grid-based algorithm for Bader charge allocation , J. Comp. Chem. 28 , 899-908 (2007). G. Henkelman, A. Arnaldsson, and H. Jónsson, A fast and robust algorithm for Bader decomposition of charge density , Comput. Mater. Sci. 36 , 254-360 (2006). This program was written by Andri Arnaldsson , Wenjie Tang , and Sam Chill , Graeme Henkelman . Improvements to the original algorithm were developed by Ed Sanville (Loughborough University, UK). With contributions from: Johannes Voss (DTU), Erik McNellis (FHI), and Matthew Dyer (Liverpool) Voronoi Polyhedra The electron density contains cusps at the nuclei and decays exponentially away from the nuclei. This suggests that the placement of quadrature points should be dependent on the nuclear positions. This can be achieved by partitioning space into regions (Voronoi Polyhedra), each containing only one atom. The total integral is the sum of single atom integrals. There are a number of such partitioning schemes ; however, the scheme used by the package produced in our research group, Q-C HEM , is that developed by Becke . The integral is split up via weight functions where the weight functions obey These weights are constructed to be almost unity if A is the closest atom, and almost zero in the vicinity of other atoms. This is achieved through the variable , defined by where r A and r B are the distance from atoms A and B , while R AB is the interatomic distance between atoms A and B . The weight function is then described by where the function is defined as Becke then removed the discontinuity at by redefining as a function with Above is the implementation in Q-CHEM. Becke's original implementation also includes a correction for atomic size. This is accomplished by a change in variable, working with instead of , where with a AB defined by where R A and R B are the Bragg-Slater radii . Each of these single-center integrals is then calculated with a spherical polar quadrature grid (requiring the insertion of more weights, w i ), making the final expression for E xc :
转自： http://www.nsc.liu.se/~pla/blog/2012/09/26/vaspkpar/ Testing the K-point Parallelization in VASP SEP 26 TH , 2012 VASP 5.3.2 finally introduced official support for k-point parallelization. What can we expect from this new feature in terms of performance? In general, you only need many k-points in relatively small cells, so up front we would expect k-point parallelization to improve time-to-solution for small cells with hundreds or thousands of k-points. We do have a subset of users at NSC, running big batches of these jobs, so this may be a real advantage in the prototyping stage of simulations, when the jobs are set up. In terms of actual job throughput for production calculations, however, k-point parallelization should not help much, as the peak efficiency is reached already with 8-16 cores on a single node. So let’s put this theory to test. Previously, I benchmarked the 8-atom FeH system with 400 k-points for this scenario. The maximum throughput was achieved with two 8-core jobs running on the same node, and the time-to-solution peaked at 3 minutes (20 jobs/h) with 16 cores on one compute node. What can k-point parallelization do here? KPAR is the new parameter which controls the number of k-point parallelized groups. KPAR=1 means no k-point parallelization, i.e. the default behavior of VASP. For each bar in the chart, the NPAR value has been individually optimized (and is thereby different for each number of cores). Previously, this calculation did not scale at all beyond one compute node (blue bars), but with KPAR=8 (purple bars), we can get close to linear (1.8x) speed-up going from 1 to 2 nodes, cutting the time-to-solution in half. As suspected, in terms of efficiency, the current k-point parallelization is not more efficient than the old scheme when running on a single node, which means that peak throughput remains the same at roughly 24 jobs/h per compute node. This is a little surprising, given that there should be overhead associated with running two job simultaneously on a node, compared to using k-point parallelization. What must be remembered, though, is that it is considerably easier to handle the file and job management for several sequential KPAR runs vs juggling several jobs per node with many directories, so in this sense, KPAR seems like a great addition with respect to workflow optimization. Posted by Peter Larsson Sep 26 th , 2012 转自： http://www.nsc.liu.se/~pla/blog/2012/11/26/vaspkpar2/ K-point Parallelization in VASP, Part 2 NOV 26 TH , 2012 Previously, I tested the k-point parallelization scheme in VASP 5.3 for a small system with hundreds of k-points. The outcome was acceptable, but less than stellar. Paul Kent (who implemented the scheme in VASP) suggested that it would be more instructive to benchmark medium to large hybrid calculations with just a few k-points, since this was the original use case, and consequently where you would be able to see the most benefit. To investigate this, I ran a 63-atom MgO cell with HSE06 functional and 4 k-points over 4 to 24 nodes: A suitable number of bands here is 192, so the maximum number of nodes we could expect to use with standard parallelization is 12, due to the fact that 12 nodes x 16 cores/node = 192 cores. And we do see that KPAR=1 flattens out at 1.8 jobs/h on 12 nodes. But with k-point parallelization, the calculation can be split into “independent” groups, each running on 192 cores. This enables us, for example, to run the job on 24 nodes using KPAR=2, which in this case translates into a doubling of speed (4.0 jobs/h), compared to the best case scenario without k-point parallelization. So there is indeed a real benefit for hybrid calculations of cells that are small enough to need a few k-points. And remember that in order for the k-point parallelization to work correctly with hybrids, you should set: NPAR = total number of cores / KPAR.
转自： http://cms.mpi.univie.ac.at/vasp/guide/node138.html The optimum setting of NPAR and LPLANE depends very much on the type of machine you are running. Here are a few guidelines SGI power challenge: Usually one is running on a relatively small number of nodes, so that load balancing is no problem. Also the communication band width is reasonably good on SGI power challenge machines. Best performance is often achived with LPLANE = .TRUE. NPAR = 1 NSIM = 1Increasing NPAR usually worsens performance. For NPAR=1 we have in fact observed a superlinear scaling w.r.t. the number of nodes in many cases. This is due to the fact that the cache on the SGI power challenge machines is relatively large (4 Mbytes); if the number of nodes is increased the real space projectors (or reciprocal projectors) can be kept in the cache and therefore cache misses decrease significantly if the number of nodes are increased. SGI Origin: The SGI Origin behaves quite different from the SGI Power Challenge. Mainly because the memory bandwidth is a factor of three faster than on the SGI Power Challenge. The following setting seems to be optimal when running on 4-16 nodes: LPLANE = .TRUE. NPAR = 4 NSIM = 4Contrary to the SGI Power Challenge superlinear scaling could not be observed, obviously because data locality and cache reusage is no only of minor importance on the Origin 2000. T3D, T3E On many T3D, T3E platforms one is forces to use a huge number of nodes. In that case load balancing problems and problems with the communication bandwidth are likely to be experienced. In addition the cache is fairly small on T3E and T3D machines so that it is impossible to keep the real space projectors in the cache with any setting. Therefore, we recommend to set NPAR on these machines to (explicit timing can be helpful to find the optimum value). The use of LPLANE = .TRUE. is only recommend if the number of nodes is significantly smaller than NGX, NGY and NGZ. In summary the following setting is recommended LPLANE = .FALSE. NPAR = sqrt(number of nodes) NSIM = 1 转自： http://www.nsc.liu.se/~pla/blog/2012/02/22/nparnsim/ Optimizing NPAR and NSIM FEB 22 ND , 2012 Our VASP users at NSC are sometimes asking about how to set NSIM and NPAR for good parallel performance in VASP. I wrote about NPAR before. But what about the NSIM parameter? The VASP manual says that NSIM=4 by default. It means that 4 bands are optimized in parallel in the RMM-DIIS algorithm, which allows VASP to exploit matrix-matrix BLAS operations instead of matrix-vector operations. But the NSIM/NPAR parameters should be adjusted based on actual underlying hardware (network, typ of processor, caches etc). Here are some results for the 24-atom PbSO4 cell running on a single Kappa compute node. Each bar in the chart below represents the average of three runs. It looks like NPAR and NSIM are largely indepedent factors, with NPAR being the most important one. Varying NPAR can give a performance boost of up to ca 50%, while varying NSIM gives about 10%. The internal variability between runs is less than 1% in this case, so the differences are real. We can conclude that NPAR=1 is optimal for a single-node job, as expected, and that NSIM=2 might be beneficial, instead of keeping the default NSIM. A more realistic example is a 128-atom Li2FeSiO4 supercell. This one we run on 4 nodes (32 cores) on Matter. It is a highly symmetric case with 512 bands. Like before, 3 runs were made for each data point. We find the best performance for NPAR=2/4, in line with previous results. But here, the default NSIM=4 setting seems to produce the worst performance, and the influence of NSIM is higher (up to +20% speed). The optimal choice seems to be NSIM=16. It is tempting to conjecture that NSIM should be increased even more for larger jobs. To investigate the upper limit of VASP jobs, let us look at the NiSi supercell with 504-atoms. It takes about 23 minutes to run a full SCF cycle on 32 Matter nodes. The outcome is not so encouraging, however: NSIM=16 does not deliver an increase in performance, and the influence of smaller NSIM values is dwarfed by other measurement errors. In this case, a likely culprit is the variation in MPI performance when running over many Infiniband switches. So for large jobs, NSIM seems to make less difference. You can leave it at default value. In conclusion: Use NPAR = 1 and NSIM = 2 for single-node jobs Use NPAR = nodes/2 (or nodes) and NSIM=2 for medium jobs. If you have enough bands per core and want to optimize, you can try NSIM=8/16 and see if it helps. Use NPAR = sqrt(nodes) and NSIM=4 for large jobs. 转自； http://www.nsc.liu.se/~pla/blog/2013/03/25/vaspabisko/ Test setup Here, we will be looking at the Li2FeSiO4 supercell test case with 128 atoms. I am running a standard spin-polarized DFT calculation (no hybrid), which I run to self-consistency with ALGO=fast. I adjusted the number of bands to 480 to better match the number of cores per node. Naive single node test A first (naive) test is to characterize the parallel scaling in a single compute node, without doing anything special such as process binding. This produced an intra-node scaling that looks like this after having tuned the NPAR values: Basically, this is what you get when you ask the queue system for 12,16,24,36,48 cores on 1 node with exclusively running rights (no other jobs on the same node), and you just launch VASP with srun $VASPin the job script. We see that we get nice intra-node scaling. In fact, it is much better than expected, but we will see in the next section that this is an illusion. The optimal choice of NPARturned out to be: 12 cores NPAR=1 16 cores NPAR=1 24 cores NPAR=3 36 cores NPAR=3 48 cores NPAR=6 This was also surprising, since I had expected NPAR=8 to be optimal. With these settings, there would be MPI process groups of 6 ranks which exactly fit in a NUMA zone. Unexpectedly, NPAR=6 seems optimal when using all 48 cores, and either NPAR=1 or NPAR=3 for the other cases. This does not fit the original hypothesis, but a weakness in our analysis is that we don’t actually know were the processes end up in this scenario, since there is no binding. The only way that you can get a symmetric communication pattern with NPAR=6 is to place ranks in a round robin scheme around each NUMA zone or socket. Perhaps this is what the Linux kernel is doing? An alternative hypothesis is that the unconventional choice of NPAR creates a load imbalance that may actually be beneficial because it allows for better utilization of the second core in each module. To explore this, I decided to test different binding scenarios. The importance of process binding To bind MPI processes to a physical core and prevent the operating system from moving them on around inside the compute node, you need to give extra flags to either srun or your MPI launching command such as mpirun. On Abisko, we use srun, where binding is controlled through SLURM by setting e.g. in the job script: srun --cpu_bind=rank ... This binds the MPI rank #1 to core #1, and so on in a sequential manner. It is also possible to explicitly specify where each rank should go. The following example binds 24 ranks to alternating cores, so that there is one rank running per Interlagos module: srun --cpu_bind=map_cpu=0,2,4,6,8,10,12,14,16,18,20,22,24,26,28,30,32,34,36,38,40,42,44,46 ... In this scheme, neighboring ranks are close to each other: i.e. rank #1 and #2 are in the same NUMA zone. The aim is to maximize the NUMA locality. The third type of binding I tried was to distribute the ranks in a round-robin scheme in steps of 6 cores. The aim is to minimize NUMA locality, since neighboring ranks are far apart from each other, i.e. rank #1 and #2 are in different NUMA zones. srun --cpu_bind=map_cpu=0,6,12,18,24,30,36,42,2,8,14,20,26,32,38,44,4,10,16,22,28,34,40,46 ... Below are the results when comparing the speed of running with 48 cores and 24 cores with different kinds of process binding. The 48 core runs are with NPAR=6 and the 24 cores with NPAR=3. It turns out that you can get all of the performance, and even more, by running with 24 cores in the “fat” mode. The trick is, however, that we need to enable the process binding ourselves. It does not happen by default when you run with half the number of cores per node (the “None” section in the graph). We can further observe that straight sequential process binding actually worsens performance in the 48 core scenario. Only in the round-robin NUMA scheme (“RR-NUMA”) can we reproduce the performance of the unbound case. This leads me to believe that running with no binding gets you in similar situation with broken NUMA locality, which explains why NPAR=3/6 is optimal, and not NPAR=4. The most surprising finding,however, is that the top speed was achieved not with the “alternate” binding scheme, which emphasizes NUMA memory locality, but rather with the round-robin scheme, which breaks memory locality of NPAR groups. The difference in speed is small (about 3%), but statistically significant. There are few scenarios where this kind of interleaving over NUMA zones is beneficial, so I suspect that it is not actually a NUMA issue, but rather related to memory caches. The L3 cache is shared between all cores in a NUMA zone, so perhaps the L3 cache is being trashed when all the ranks in an NPAR group are accessing it? It would be interesting to try to measure this effect with hardware counters… NSIM Finally, I also made some tests with varying NSIM: NSIM=4 is the default setting in VASP. It usually gives good performance in many different scenarios. NPAR=4 works on Abisko too, but I gained about 7% by using NPAR=8 or 12. An odd finding was that NPAR=16 completely crippled the performance, doubling the wall time compared to NPAR=4. I have no good explanation, but it obviously seems that one should be careful with too high NPAR values on Abisko. Conclusion and overall speed In terms of absolute speed, we can compare with Triolith, where one node with 16 cores can run this example in 380s (9.5 jobs/h) with 512 bands, using the optimal settings of NPAR=2 and NSIM=1. So the overall conclusion is that one Abisko node is roughly 20% faster than one Triolith node. You can easily become disappointed by this when comparing the performance per core, which is 2.5x higher on Triolith, but I think it is not a fair comparison. In reality, the performance difference per FPU is more like 1.3x, and if you compensate for the fact that the Triolith processors in reality run at much higher frequency than the listed 2.2 Ghz, the true difference in efficiency per core-GHz is closer to 1.2x. Hopefully, I can make some multi-node tests later and determine whether running with 24 cores per node and round-robin binding is the best thing there as well. Posted by Peter Larsson Mar 25 th , 2013
Q1: 最近学习USPEX, 对T0进行测试，运行 nohup matlab USPEX.m log ，大概几秒钟后，就显示如下的错误。整个过程完全没有调用vasp进行优化计算，你能帮我分析一下吗，或者要做什么修改？谢谢。 M A T L A B (R) Copyright 1984-2012 The MathWorks, Inc. R2012a (7.14.0.739) 64-bit (glnxa64) February 9, 2012 To get started, type one of these: helpwin, helpdesk, or demo. For product information, visit www.mathworks.com. ans = 0 status = Local optimisation finished status = Too many structures have errors or failed the constraints after optimization. Please check the input files. The calculation has to stop. status1 = Possible reasons: badly tuned optimization parameters or unreasonable contraints. 1, 上面的提示，完全不是INPUT.txt的设置问题，真是KD，害我不断调试了好多参数。 2. matlab普通用户的执行参考http://wenku.baidu.com/view/01f08cd126fff705cc170ad7.html 3 我的root帐户没有vasp.5.2的执行权限，而多次测试是在root帐户下进行的（CalcFold1 ） 4. 重新提交任务时，必须删除上次的Current_EXE.mat Current_ORG.mat Current_POP.mat三个文件，不然仍然会显示上面的错误提示。 Q2:moab 作用的作用提交系统与uspex结合时， submitJob_local.m需要做红色的修改，不然 function jobNumber = submitJob_local() %------------------------------------------------------------- %This routine is to check if the submitted job is done or not %One needs to do a little edit based on your own case. %1 : whichCluster (default 0, 1: local submission, 2: remote submission) %------------------------------------------------------------- %Step 1: to prepare the job script which is required by your supercomputer fp = fopen('myrun', 'w'); fprintf(fp, '#!/bin/sh\n'); fprintf(fp, '# MSUB -l nodes=1:ppn=16,walltime=01:00:00\n'); fprintf(fp, '#MSUB -N USPEX\n'); fprintf(fp, '#MSUB -q jsc\n'); fprintf(fp, '#MSUB -j oe\n'); fprintf(fp, '#MSUB -V tpt=1\n'); fprintf(fp, 'cd $PBS_O_WORKDIR \n'); fprintf(fp, 'mpiexec -np 16 /lustre/jhome15/hhb17/hhb172/bin/vasp.st.5.3 vasp.out\n'); fclose(fp); %Step 2: to submit the job with the command like qsub, bsub, llsubmit, .etc. %It will output some message on the screen like '2350873.nano.cfn.bnl.local' =unix( ) %Step 3: to get the jobID from the screen message end_marker = findstr(b,'.'); jobNumber = b(2:end-1); disp( ) %jobNumber = b(1:end_marker(1)-1); (moab提交作业时，msub myrun会空一行再显示进程号，所以这里要改为2，不然显示 提交uspex后，显示这样的提示 ??? Attempted to access end_marker(1); index out of bounds because numel(end_marker)=0. Error in == submitJob_local at 27 jobNumber = b(1:end_marker(1)-1); Error in == submitJob at 26 jobNumber = submitJob_local(); Error in == SubmitJobs_M200 at 26 POP_STRUC.POPULATION(DO_NOW).JobID = submitJob(DO_NOW); Error in == EA_M200 at 12 SubmitJobs_M200(); Error in == Start at 73 eval( ); ) Q3:USPEX调用matlab提交vasp的优化计算， 能够成功提交任务，但会显示如下的错误： 而我 单独在alcfold*目录里面msub myrun却能成功运行vasp．后通过 ldd /lustre/jhome15/hhb17/hhb172/bin/vasp.2d.5.3与在matlab里面!ldd /lustre/jhome15/hhb17/hhb172/bin/vasp.2d.5.3发现vasp依赖的库不一样，需在环境变量里面把依赖的库加到export LD_LIBRARY_PATH里
对于一些磁性体系、镧系和锕系元素及相关化合物的静态计算（电子迭代），经常会遇到“难收敛”的问题。 下面给出几个相关Flag及设置方法：（有些话从说明书上摘录的，不大好翻译，就copy了，请见谅） 1、LMAXMIX Default: LMAXMIX = 2 An additional flag controls up to which l quantum number the onsite PAW charge densities are passed through the charge density mixer. Higher l-quantum numbers are usually not handled by the mixer. In order to obtain fast convergence to the groundstate, you can try the following setting: LMAXMIX = 4 for d elements LMAXMIX = 6 for f elements 这个FLAG对于含d电子和f电子的体系是非常重要的，很大一部分体系的收敛问题可以通过设置合适的LMAXMIX值来解决。 2、ALGO, IALGO, LDIAG If the self-consistency loop does not converge within 40 steps, it will probably not converge at all. In this case you should reconsider the tags IALGO, LDIAG, and the mixing-parameters. 这是说明书上的建议。 一般情况下，或使用IALGO=48时遇到收敛问题的话，可以考虑设IALGO为38（4.5以前的版本可设为8），或设置ALGO=Normal or Fast (in VASP.4.5 and later versions)。 3、NELMDL NELMDL gives the number of non-selfconsistent steps at the beginning; if one initializes the wave functions randomly the initial wave functions are far from anything reasonable. The resulting charge density is also 'nonsense'. Therefore it makes sense to keep the initial Hamiltonian, which corresponds to the superposition of atomic charge densities, fixed during the first few steps. Choosing a 'delay' for starting the charge density update becomes essential in all cases where the SC-convergence is very bad (e.g. surfaces or molecules/clusters chains). Without setting a delay VASP will probably not converge or at least the convergence speed is slowed down. NELMDL might be positive or negative. A positive number means that a delay is applied after each ionic movement — in general not a convenient option. A negative value results in a delay only for the start-configuration. 4、mixing-parameters 对于一些难收敛的体系，可以使用“linear mixing”,具体详见VASP说明书中的“Mixing-tags”。 For an initial linear mixing (BMIX ~ 0) an optimal setting for A(AMIX) can be found easily by setting Aopt=Acurrent*Γmean. For the Kerker scheme either A or q0(i.e. AMIX or BMIX) can be optimized, but we recommend to change only BMIX and keep AMIX fixed (you must decrease BMIX if the mean eigenvalue is larger than one, and increase BMIX if the mean eigenvalue is smaller than one). 尽管VASP说明书中给出了调节AMIX和BMIX的一些较为明确的建议，但是实际去调节的时候，还是挺难的，但原则上说，是可以通过调节这两个Flag来使得收敛问题得以解决的，只是得有耐心。 5、kmesh, SIGMA 收敛问题还跟kmesh及SIGMA(当使用ISMEAR不等于-5 和-4时)的设置有关。要达到同样的精度，较小的SIGMA则需要较大的kmesh；而且，当SIGMA较小时，若kpoints不够多，也会出现难收敛的情况。 P.S. 对于不同系统的计算，问题的原因不一定一样，因而可能解决之一问题的方法也不一定会一样。其实以上大部分都来自VASP的说明书，只是我在遇到收敛问题时试过这些方法. 写得比较仓促，以后再修改及补充。 希望有这方面经验的朋友也把自己的见解提出来一起讨论。 金属表面体系。换算法，改Mixing都不管用。 收敛和k-mesh有关系吗？谢谢了 别改Mixing。 你把IBRION设为3，然后把SMASS设为0.6试试看 不好意思，我没说清楚 IBRION是对应离子步的算法 我的问题是电子步不收敛 dE 一直0.1左右反复, 200步也不收敛 改IALGO为38, 增加NBANDS都不行:( 那你试试不要采用自恰计算。但是仍然采用damping驰豫。 是不是你的初始构型离平衡态太远了。 还没遇到过这种问题，可能比较运气吧。 调一下cutoff试试。 IALGO我一直用48，5x5x1 Kpoints。我做金属的表面体系都是在100步内就OK了。 如果愿意，可以把INCAR贴出来，让大家出出注意，这样更可能有的放矢。 应该不会 有的时候甚至弛豫出来的结构再进行static的计算都不收敛 收敛跟k-mesh是有一定关系的。 你INCAR中没有设置ISMEAR，应该就是使用默认的ISMEAR=1。此时SIGMA应该也是默认值0.2。 你看看OUTCAR中的EENTRO的值，看是不是符合要求（一般要求小于1 meV/atom，这与你想要的精度有关，如果你要求精度不高的话，这个标准可以设大一些）。 有几个方向你可以试试，或许会有用： 1、如果EENTRO的值够小，你可以尝试加大SIGMA。（SIGMA较大时，相对收敛会快些） 2、增加k点多。 其实就是说，要调整k-mesh和SIGMA。一般SIGMA跟k点数会有一定的关系。 若想达到同一精度，SIGMA越小所需要的K点数就越多。如果SIGMA较小，而k点数不足，有可能会出现难收敛的情况，也可能得到不准确的结果。 3、调节AMIX与BMIX 我之前也遇到过难收敛的问题，当时调这两个参数，还是有点用的。具体调节的方法VASP上有讲到一些，但可能还是要多试才知道。我现在一时想不起来细节了，你先试试，如果有情况，再一起讨论。 VASP的说明书也有提到，如果电子迭代超过40次仍未收敛的话，一般就很难收敛了。当然，我想这只是VASP说明书的一个建议，并不是绝对的。一般设置NELM=60（默认）就够了，如果60次都不能收敛，那可能就要找找问题或是改相应的设置了。 有不对或不妥的地方，大家帮忙指正：） 从你的INCAR来看，貌似没什么问题。如果还是不收敛的话，你可以增加一点截断能，比如400 eV，并把POTIM减小至0.1。另外，添加三个tag:NELM=100； NELMDL=6； WEIMIN=0. 大家遇到收敛问题怎么解决啊？？？ 小弟仔细读过使用手册，并结合一些帖子，可是还是解决不了问题。恳请高手指点。 以下是我所用的方法和遇到的问题： 解决方法一：设置scf参数 a1. Iterations. 加大循环次数。计算过程中默认50次，我的体系无论加多大次数就是不收敛。 可能这不是主要原因。我现在一般设置一个中等大小的值400。 a2. mixing.对强阻尼体系，降低mixing值到0.01-0.05.甚至用mixing=0.001,diis n=0.保存TAPE21文件， 再用mixing=0.03,diis n=10.我都用了，都不收敛。分析结果计算过程中有能量跳跃，因此我分析我的体系就 是较强的阻尼，因此我设置mixing=0.03，这个因该可以吧？ a3. lshift, vshift.手册上未解释怎么用。 这个我不懂。请问在什么情况下用？怎么用？大家一般设置值为多少？ 解决方法二：用occupations指定电子占据。 b1. 直接指定分子占据，由于小弟的量化基础薄弱，弄不懂怎么指定。 我是按论坛上的方法之一，先算一个单点能，然后得到所有不可约表示的电子占据数。但是， 我计算单点能虽然scf收敛了，但是计算得到的能量明显不对，估计是收敛到激发态上去了。请问怎么指定到基态？ 另外，对电子占据大家都有什么具体的着？哪位有这方面的文献和资料？我算单点的输出文件在附件中，请高手帮我 指定一下基态的电子占据。先谢谢了！！！ 这是最重要，最有用的方法，可是我不懂，郁闷啊！！！ b2. 用keeporbitals或freeze控制能量跳跃，但是我的体系用这个关键词还是无法得到几何结构收敛。我只用过 keeporbitals=20,10,5. 但是scf有几个收敛，可是最后得不到几何结构收敛。 大家一般设置值为多少？freeze又该如何设置？ b3. smearq. 手册上说先用smearq=一个值。几何收敛后再用smearq=0, restart 前面的TAPE21文件。 我设置 smearq=0.2,几何结构收敛了，但是我用smearq=0,restart不收敛。我用smearq=0.1计算连几何结构 都不收敛。请问大家设置smearq=多少合适？用smearq且令smearq=0之后的结果可信度有多大啊？我看过几篇文献， 但是都没有什么有用的东西。 我已收到你的邮件了,个人觉得在这么多方法中用smearq 再 occupations比较有效: 先用启用smearing (例如:smearq=0.02), 再根据结果文件指认电子占据一般都能解决问题, 例如: occupations A1 a//b B1 a//b E a//b . . . end. Dear Vasp Community: I'm trying to do some very basic defect calculations in a class of complex oxides with the stoichiometry AB6O12, where A=U,W,Re and B=Y,Yb,Lu. In most cases, especially when A=W, I haven't had any problems. But, I have had problems with convergence in two cases that I haven't been able to make any progress on. When A=U and B=Y, I can converge the perfect structure very nicely. But, when I look at an antisite defect, by replacing one Y with U, the calculation doesn't converge. The electronic iterations get stuck, changing by +/- the same amount. A number of colleagues have suggested I change the smearing, set ISPIN=2, set LMIXMAX=6, etc. Everything I've tried, either doesn't change the behavior or does fix the electronic iterations, but after 3 or 4 ionic iterations, the forces on the ions blow up (become unphyiscally large, like 100 eV/angstrom). The second problem, when A=Re and B=Lu, is that I don't get any kind of convergence, even in the perfect case. I get a warning at the top VERY BAD NEWS! internal error in subroutine IBZKPT: Reciprocal lattice and k-lattice belong to different class of lattices. 168 But all I've done compared to the A=U, B=Y case is change the POTCAR file. If I let it go, I get similar convergence problems as before. Any suggestions on what I might try to fix these problems would be greatly appreciated. I'm using ISYM=0 and a Monkhorst-Pack grid of 2x2x2 for these calculations, along with the GGA PAW pseudopotentials. Thank you, Blas Back to top admin Tue Jan 08 2008, 01:27PM posts 1815 1) for the electronic convergence --) LMAXMIX=6 should be set if you have rare-earth elements in the cell, for d-elements, please set it to 4 --) it may help to decrease the mixing parameters (AMIX, BMIX, and in addition AMIX_MAX and BMIX_MAX for the spin-polarized runs to improve convergence ) --) also, it may help to increase NELMDL (the number of non-selconsistent steps at the beginning) to improve the pre-convergence of the wavefunction --) maybe it also helps to use L(S)DA+U ( please check the electronic structure and compare to experiment, if results are available) 2) for the forces which suddenly increase: --) please check (in OSZICAR) if each ionic step was fully converged electronically up to the requested convergence limit (EDIFF) before the forces are calculated. (i.e. whether the forces were calculated after el. convergence is reached or because the number of electronic convergence steps (NELM=60, by default) was reached without actually being converged. The latter yields unreasonable forces. (in that case, please either increase NELM and/or try one of the above mentioned to improve the electronic convergence) --) if each of the ionic steps was converged, please check the relaxation history (XDATCAR file). Sometimes, especially the first ionic steps 'overshoot' a little. If that is the case, please choose a different relaxation algorithm (IBRION) and/or decrease the step width (POTIM) 3) the VERY BAD NEWS message is due to different lattice types obtained for the unit cell and the symmetrized k-mesh. (triclinic and simple monoclinic as it seems). maybe the axes ratios of the unit cells are not 1:1:1 ? please check. for the symmetry analysis, the Bravais Matrix , atoms' positions, and, if set, magnetic moments and atoms' velocities are taken into account. If nothing has been changed except the atoms in the PP, it seems strange that this warning shows up for some calculations whereas it does not for others. In any case, if ISYM=0 is set, symmetrization is not taken into account in the calculation later on, the only point to be careful with is the axis' ratio then. LDAU = .TRUE. LDAUL = LDAUU = LDAUJ =
在这篇文章中，我将首先介绍Partial Charge的概念，以及如何用VASP具体的计算Partial Charge。首先，所谓的Partial Charge是针对与Total Charge来说的，指的是某个能量范围、某个K点或者某个特定的态所对应的电荷密度。在文献中最常见的是价带顶部，导带底部，表面态或者局域态所对应的Partial Charge。通过分析这些态所对应的Partial Charge，可以得到体系的一些性质，比如局域态具体的是局域在哪个原子上等。我将通过具体的例子说明如何用VASP进行Partial Charge Analysis。 进行Partial Charge Analysis的第一步是进行自洽的计算，得到体系的电子结构。这一步的计算采用通常的INCAR和KPOINTS文件。在自洽计算结束后，我们需要保存WAVECAR文件。（通过在INCAR文件中设置LWAVE=TRUE实现）在这个例子中，假设我们需要计算一个硅纳米线的导带和价带的Partial Charge。硅纳米线的结构如下： http://www.quantumchemistry.net/Experience/UploadFiles_6872/200603/20060331155154521.jpg 第二步是画出能带结构，以决定你需要画哪条能带的那个K点的态所对应的Partial Charge。关于具体如何用VASP画能带，请参见用VASP4.6计算晶体硅能带实例一文。我们得到硅纳米线的能带结构如下： 画能带时有些小技巧。你可以用一些支持列模块的编辑器，如UltraEdit，将OUTCAR里的各个K点所对应的本征值粘贴到Origin中。这一步完成后，在Origin中做一个矩阵转置，然后将K点坐标贴到第一列，并将其设为X坐标。如此画出来的基本上就是能带图了。在Origin中可以通过设置纵轴范围来更加清楚的区分费米能级附近的各条能带。如上的硅纳米线所对应的能带结构图如下： http://www.quantumchemistry.net/Experience/UploadFiles_6872/200603/20060331155540648.jpg 决定画哪条能带，或者那些感兴趣的K点之后，有如下几种方法计算不同的Partial Charge。如果你希望计算价带顶端的Partial Charge，则需要首先通过能带结构图确定价带的能带标号。需要注意，进行Partial Charge分析必须要保留有自洽计算的WAVECAR才可以。 第一种Partial Charge分析的INCAR ISTART = 1 job : 0-new 1-cont 2-samecut ICHARG = 1 charge: 1-file 2-atom 10-const LPARD=.TRUE. IBAND= 20 21 22 23 KPUSE= 1 2 3 4 LSEPB=.TRUE. LSEPK=.TRUE. 这样的INCAR给出的是指定能带，指定K点所对应的Partial Charge。分析导带、价带等的Partial Charge特性，通常采用的都是这种模式。 第二种Partial Charge分析的INCAR ISTART = 1 job : 0-new 1-cont 2-samecut ICHARG = 1 charge: 1-file 2-atom 10-const LPARD=.TRUE. EINT = -10.3 -5.1 LSEPB=.FALSE. LSEPK=.FALSE. 这样的INCAR给出的是在 能量之间的Partial Charge。这种模式适合于分析某个能量区间内的波函数的性质。 第三种Partial Charge分析的INCAR ISTART = 1 job : 0-new 1-cont 2-samecut ICHARG = 1 charge: 1-file 2-atom 10-const LPARD=.TRUE. NBMOD=-3 EINT = -1 LSEPB=.FALSE. LSEPK=.FALSE. 这样的INCAR给出的是从 能量之间的Partial Charge。这种模式最利于分析费米面附近的波函数的性质。 用第一种方法，我们可以得到硅纳米线价带顶部和导带底部的Partial Charge如下： http://www.quantumchemistry.net/Experience/UploadFiles_6872/200603/20060331155753136.jpg • LPARD: Evaluate partial (band and/or k-point) decomposed charge density. We want to stress again, that the wavefunctions read from WAVECAR must be converged in a separate prior run. If only LPARD is set (and none of the tags discussed below), the total charge density is evaluated from the wavefunctions and written to CHGCAR. • There are several ways how to specify for which bands the charge density is evaluated: In general the input lines with IBAND, EINT and NBMOD control this respect of the routine: • IBAND: Calculate the partial charge density for all bands specified in the array IBAND. If IBAND is specified in the INCAR file and NBMOD is not given, NBMOD is set automatically to the size of the array. If IBAND is for instance IBAND= 20 21 22 23 the charge density will be calculated for bands 20 to 23. • EINT: Specifies the energy range of the bands that are used for the evaluation of the partial charge density. Two real values should be given, if only one value is specified, the second one is set to ǫf . If EINT is given and NBMOD is not specified, NBMOD is set automatically to -2. • NBMOD: This integer variable can take the following values 0 Number of values in the array IBAND. If IBAND is specified, NBMOD is set automatically to the correct value (in that case NBMOD should not be set manually in the INCAR file) 0 Take all bands to calculate the charge density, even unoccupied bands are taken into account. -1 Calculate the total charge density as usual. This is the default value if nothing else is given. -2 Calculate the partial charge density for electrons with there eigenvalues in the range specified by EINT. -3 The same as before, but the energy range is given vs. the Fermi energy. • KPUSE: KPUSE specifies which k-points are used in the evaluation of the partial dos. KPUSE is an array of integer values. KPUSE= 1 2 3 4 means that the charge density is evaluated and summed for the first four k-points. Be careful: VASP changes the kpoint weights if KPUSE is specified. • LSEPB: Specifies whether the charge density is calculated for every band separately and written to a file PARCHG.nb.⋆ (TRUE) or whether charge density is merged for all selected bands and write to the file PARCHG.ALLB.⋆ or PARCHG. Default is FALSE. • LSEPK: Specifies whether the charge density of every k-point is write to the files PARCHG.⋆.nk (TRUE) or whether it is merged (FALSE) to a single file. If the merged file is written, then the weight of each k-point is determined from the KPOINTS file, otherwise the kpoints weights of one are chosen.
转自： http://cms.mpi.univie.ac.at/vasp/vasp/vdW_DF_functional_Langreth_Lundqvist_et_al.html vdW-DF functional of Langreth and Lundqvist et al.The vdW-DF proposed by Dion et al. is a non-local correlation functional that approximately accounts for dispersion interactions . In VASP the method is implemented using the algorithm of Roman-Perez and Soler which transforms the double real space integral to reciprocal space and reduces the computational effort. Several propsed versions of the method can be used : the original vdW-DF , vdW-DF with exchange functionals optimised for the correlation part , and the vdW-DF2 of Langreth and Lundqvist groups . N.B. : This feature has been implemented by J. Klimeš. If you make use of the vdW-DF functionals presented in this section, we ask that you cite the following paper: J. Klimeš, D. R. Bowler, and A. Michaelides, Phys. Rev. B 83 , 195131 (2011). Correlation functionals The method is invoked by setting LUSE_VDW = .TRUE. Moreover, the PBE correlation correction needs to be removed since only LDA correlation is used in the functionals. This is done by setting AGGAC = 0.0000 The two tags above need to be used for all the following functionals. Exchange functionals To use the different exchange functionals, the GGA tag needs to be set appropriately. The original version of Dion et al uses revPBE which can be set by GGA = RE More accurate exchange functionals for the vdW correlation functional have been proposed in and . They can be used by setting GGA = OR for optPBE, GGA = BO PARAM1 = 0.1833333333 PARAM2 = 0.2200000000 for the optB88 functional , or GGA = MK PARAM1 = 0.1234 PARAM2 = 1.0000 for the optB86b functional . For the vdW-DF2 functional the rPW86 exchange functional is used: GGA = ML moreover, the vdW functional needs to be changed to the vdW2 correlation which requires only a change of a parameter: Zab_vdW = -1.8867 An overview of the performance of the different approaches can be found for example in for gas phase clusters and in for solids. Important remarks : The method needs a precalculated kernel which is distributed via the VASP download portal ( VASP - src - vdw_kernel.bindat ) and on the ftp server ( vasp5/src/vdw_kernel.bindat ). If VASP does not find this file, the kernel will be calculated. This, however, is rather demanding calculation. The kernel needs to be either copied to the VASP run directory for each calculation or can be stored in a central location and read from there. The location needs to be set in routine PHI_GENERATE. This does not work on some clusters and the kernel needs to be copied into the run directory in such cases. Currently the evaluation of the vdW energy term is not done fully within the PAW method but the sum of the pseudo-valence density and partial core density is used. This approximation works rather well, as is discussed in , and the accuracy generally increases when the number of valence electrons is increased or when harder PAW datasets are used . For example, for adsorption it is recommended to compare the adsorption energy obtained with standard PAW datasets and more-electron POTCARs for both PBE calculation and vdW-DF calculation to asses the quality of the results. The optimisation of the cell (ISIF=3 and higher) is currently not possible. The spin polarised calculations are possible, but strictly speaking the non-local vdW correlation is not defined for spin-polarised systems. For spin-polarised calculation the non-local vdW correlation energy is evaluated on the sum of the spin-up and spin-down densities. The evaluation of the vdW energy requires some additional time. Most of it is spent on performing FFTs to evaluate the energy and potential. Thus the additional time is determined by the number of FFT grid points in the calculation, basically size of the simulation cell. It is almost independent on the number of the atoms in the cell. Thus the relative cost of the vdW-DF method depends on the ``filling" of the cell and increases with the amount of vacuum in the cell. The relative increase is high for isolated molecules in large cells, but small for solids in smaller cells with many k-points. This feature has been implemented by J. Klimeš. If you make use of the vdW-DF functionals presented in this section, we ask that you cite the following paper: J. Klimeš, D. R. Bowler, and A. Michaelides, Phys. Rev. B 83 , 195131 (2011). http://cms.mpi.univie.ac.at/vasp/vasp/DFT_D2_method_Grimme.html#tab:grimme DFT-D2 method of Grimme LVDW= .TRUE. | .FALSE. (Available as of VASP.5.2.11) Default: LVDW=.FALSE. Popular density functionals are unable to describe correctly van der Waals interactions resulting from dynamical correlations between fluctuating charge distributions. A pragmatic method to work around this problem has been given by the DFT-D approach , which consists in adding a semi-empirical dispersion potential to the conventional Kohn-Sham DFT energy: (90) In the DFT-D2 method of Grimme , the van der Waals interactions are described via a simple pair-wise force field, which is optimized for several popular DFT functionals. The dispersion energy for periodic systems is defined as (91) where the summations are over all atoms and all translations of the unit cell , the prime indicates that for , is a global scaling factor, denotes the dispersion coefficient for the atom pair , is a position vector of atom after performing translations of the unit cell along lattice vectors. In practice, terms corresponding to interactions over distances longer than a certain suitably chosen cutoff radius contribute only negligibly to and can be ignored. The term is a damping function (92) whose role is to scale the force field such as to minimize contributions from interactions within typical bonding distances. Combination rules for dispersion coefficients and vdW radii are (93) and (94) respectively. The global scaling parameter has been optimized for several different DFT functionals such as PBE ( ), BLYP ( ), and B3LYP ( ). The DFT-D2 method can be activated by setting LVDW=.TRUE. Optionally, the forcefield parameters can be controlled using the following flags (the default values are listed): VDW_RADIUS = 30.0 cutoff radius () for pair interactions VDW_SCALING = 0.75 global scaling factor VDW_D = 20.0 damping parameter VDW_C6 = ,... parameters ( ) for each species defined in POSCAR VDW_R0 = ,... parameters () for each species defined in POSCAR The default values for VDW_C6 and VDW_R0 are compiled in Tab. 2 . As the potential energy, interatomic forces as well as stress tensor are corrected by adding contribution from the forcefield, simulations such as the atomic and lattice relaxations, molecular dynamics, and vibrational analysis can be performed. The number of atomic pairs contributing to and the estimated vdW energy are written in OUTCAR (check lines following the expression 'Grimme's potential'). The forces and stresses written in OUTCAR contain the vdW correction but the corrected energy should be read from OSZICAR (energies in OUTCAR do not contain the vdW term). IMPORTANT NOTE: The defaults for VDW_C6 and VDW_R0 are defined only for the first five rows of periodic table of elements (see Tab. 2 ) - if the system contains other elements the user must provide the corresponding parameters. Table 2: Parameters used in the empirical force-field of Grimme . Element C R Element C R Jnm mol Jnm mol H 0.14 1.001 K 10.80 1.485 He 0.08 1.012 Ca 10.80 1.474 Li 1.61 0.825 Sc-Zn 10.80 1.562 Be 1.61 1.408 Ga 16.99 1.650 B 3.13 1.485 Ge 17.10 1.727 C 1.75 1.452 As 16.37 1.760 N 1.23 1.397 Se 12.64 1.771 O 0.70 1.342 Br 12.47 1.749 F 0.75 1.287 Kr 12.01 1.727 Ne 0.63 1.243 Rb 24.67 1.628 Na 5.71 1.144 Sr 24.67 1.606 Mg 5.71 1.364 Y-Cd 24.67 1.639 Al 10.79 1.716 In 37.32 1.672 Si 9.23 1.716 Sn 38.71 1.804 P 7.84 1.705 Sb 38.44 1.881 S 5.57 1.683 Te 31.74 1.892 Cl 5.07 1.639 I 31.50 1.892 Ar 4.61 1.595 Xe 29.99 1.881
Platform: Dell T7500 workstation; vasp 4.638; Ubuntu 12.04 (also ubuntu 11.04) software for compilation: gcc; gfortran(maybe necessary); ifort MKL (l_fcompxe_intel64_2011.9.293); openmpi-1.60 for parallelization version Procedures as follows 1. install ifort + mkl 2. make libintel64 3. install openmpi 4. make vasp.4.lib 5. make vasp Please reference to attachment file for detailed procedures and Makefiles both for series and parallel version compilaton process series.zip ; parallelization.zip
I successfully installed VASP 5.2.12 on my DELL desktop, 2 Intel Xeon Quad cores, 8 CPUs in total. I spend a lot of time for it, I'd like to share the experience and I wish you to spend less. :) System: Ubuntu 11.04, 64 bit C compiler: gcc (you can also use icc, i installed Intel C++ Composer XE) Fortran compiler: Intel Fortran Composer XE 2011 GotoBLAS2-1.13 ( http://www.tacc.utexas.edu/tacc-projects/gotoblas2/downloads) FFTW 3.2.2 ( http://www.fftw.org/) , you can use newer if available Openmpi 1.4.3 ( http://www.open-mpi.org/) , you can use newer if available Since I use ubuntu I almost always use sudo command like sudo make TARGET=CORE2 USE_THREAD=0 I don't write sudo thing further in this post since it's feature of ubuntu. Parallel version 1) installation of GotoBLAS2 Basically it's suggested to type just make but this does not always work. I advise to build like make TARGET=CORE2 USE_THREAD=0 be careful with TARGET flag, don't use default or at list check, you could have problems because of that. Second one has to specify USE_THREAD=0 to have not threaded version, it's faster when used with OpenMPI. You can also specify the compiler: make CC=gcc (or icc) FC=ifort TARGET=CORE2 USE_THREAD=0 I had problems to compile with icc, so I used gcc + ifort 2) Installation of FFTW To compile read supported readme files. ./configure --prefix=/folder/for/fftw make make install I didn't specify any CC=gcc for ./configure, so I used default for me gcc compiler . 3) Installation of OpenMPI To compile OpenMPI read supported readme files. ./configure CC=icc CXX=icc F77=ifort F90=ifort --prefix=/folder/for/openmpi make all install Again I had problems with icc, so I compiled it with gcc, again I used gcc + ifort 4) VASP lib I used makefile.linux_ifc_P4 file the only change you need, it's FC=ifort 5) Compilation of VASP Again I started from makefile.linux_ifc_P4 file. Things I changed: 1) FFLAGS = -FR -names lowercase -assume byterecl 2) If you use GotoBLAS BLAS= /home/alex/VASP/GotoBLAS2/libgoto2.so or you can use MKL, but MKL is slow for parallel version for my case BLAS=-L/opt/intel/composerxe-2011.4.191/mkl/lib/intel64 -lmkl_intel_lp64 -lmkl_intel_thread -lmkl_core -liomp5 -lpthread 3) LAPACK LAPACK= ../vasp.5.lib/lapack_double.o 4) fortran compiler/linker for mpi FC=/home/alex/VASP/openmpi/bin/mpif77 5) FFTW for mpi FFT3D = fftmpiw.o fftmpi_map.o fftw3d.o fft3dlib.o /home/alex/VASP/fft/lib/libfftw3.a INCS = -I/home/alex/VASP/fft/include 6) I tried to set -O3 optimization flag in line OFLAG=-O2 -ip -ftz but I didn't find any advantage and I left it as it is. 7) Probably I forgot something but anyway you can see the makefile. Makefile itself: .SUFFIXES: .inc .f .f90 .F #----------------------------------------------------------------------- # Makefile for Intel Fortran compiler for Pentium/Athlon/Opteron # bases systems # we recommend this makefile for both Intel as well as AMD systems # for AMD based systems appropriate BLAS and fftw libraries are # however mandatory (whereas they are optional for Intel platforms) # # The makefile was tested only under Linux on Intel and AMD platforms # the following compiler versions have been tested: # - ifc.7.1 works stable somewhat slow but reliably # - ifc.8.1 fails to compile the code properly # - ifc.9.1 recommended (both for 32 and 64 bit) # - ifc.10.1 partially recommended (both for 32 and 64 bit) # tested build 20080312 Package ID: l_fc_p_10.1.015 # the gamma only mpi version can not be compiles # using ifc.10.1 # # it might be required to change some of library pathes, since # LINUX installation vary a lot # Hence check ***ALL*** options in this makefile very carefully #----------------------------------------------------------------------- # # BLAS must be installed on the machine # there are several options: # 1) very slow but works: # retrieve the lapackage from ftp.netlib.org # and compile the blas routines (BLAS/SRC directory) # please use g77 or f77 for the compilation. When I tried to # use pgf77 or pgf90 for BLAS, VASP hang up when calling # ZHEEV (however this was with lapack 1.1 now I use lapack 2.0) # 2) more desirable: get an optimized BLAS # # the two most reliable packages around are presently: # 2a) Intels own optimised BLAS (PIII, P4, PD, PC2, Itanium) # http://developer.intel.com/software/products/mkl/ # this is really excellent, if you use Intel CPU's # # 2b) probably fastest SSE2 (4 GFlops on P4, 2.53 GHz, 16 GFlops PD, # around 30 GFlops on Quad core) # Kazushige Goto's BLAS # http://www.cs.utexas.edu/users/kgoto/signup_first.html # http://www.tacc.utexas.edu/resources/software/ # #----------------------------------------------------------------------- # all CPP processed fortran files have the extension .f90 SUFFIX=.f90 #----------------------------------------------------------------------- # fortran compiler and linker #----------------------------------------------------------------------- #FC=ifc # fortran linker #FCL=$(FC) #----------------------------------------------------------------------- # whereis CPP ?? (I need CPP, can't use gcc with proper options) # that's the location of gcc for SUSE 5.3 # # CPP_ = /usr/lib/gcc-lib/i486-linux/2.7.2/cpp -P -C # # that's probably the right line for some Red Hat distribution: # # CPP_ = /usr/lib/gcc-lib/i386-redhat-linux/2.7.2.3/cpp -P -C # # SUSE X.X, maybe some Red Hat distributions: CPP_ = ./preprocess $*.F | /usr/bin/cpp -P -C -traditional $*$(SUFFIX) #----------------------------------------------------------------------- # possible options for CPP: # NGXhalf charge density reduced in X direction # wNGXhalf gamma point only reduced in X direction # avoidalloc avoid ALLOCATE if possible # PGF90 work around some for some PGF90 / IFC bugs # CACHE_SIZE 1000 for PII,PIII, 5000 for Athlon, 8000-12000 P4, PD # RPROMU_DGEMV use DGEMV instead of DGEMM in RPRO (depends on used BLAS) # RACCMU_DGEMV use DGEMV instead of DGEMM in RACC (depends on used BLAS) # tbdyn MD package of Tomas Bucko #----------------------------------------------------------------------- CPP = $(CPP_) -DHOST=\"LinuxIFC\" \ -DCACHE_SIZE=12000 -DPGF90 -Davoidalloc -DNGXhalf \ # -DRPROMU_DGEMV -DRACCMU_DGEMV #----------------------------------------------------------------------- # general fortran flags (there must a trailing blank on this line) # byterecl is strictly required for ifc, since otherwise # the WAVECAR file becomes huge #----------------------------------------------------------------------- FFLAGS = -FR -names lowercase -assume byterecl #----------------------------------------------------------------------- # optimization # we have tested whether higher optimisation improves performance # -axK SSE1 optimization, but also generate code executable on all mach. # xK improves performance somewhat on XP, and a is required in order # to run the code on older Athlons as well # -xW SSE2 optimization # -axW SSE2 optimization, but also generate code executable on all mach. # -tpp6 P3 optimization # -tpp7 P4 optimization #----------------------------------------------------------------------- # ifc.9.1, ifc.10.1 recommended OFLAG=-O2 -ip -ftz OFLAG_HIGH = $(OFLAG) OBJ_HIGH = OBJ_NOOPT = DEBUG = -FR -O0 INLINE = $(OFLAG) #----------------------------------------------------------------------- # the following lines specify the position of BLAS and LAPACK # VASP works fastest with the libgoto library # so that's what we recommend #----------------------------------------------------------------------- # mkl.10.0 # set -DRPROMU_DGEMV -DRACCMU_DGEMV in the CPP lines #BLAS=-L/opt/intel/mkl100/lib/em64t -lmkl -lpthread # even faster for VASP Kazushige Goto's BLAS # http://www.cs.utexas.edu/users/kgoto/signup_first.html # parallel goto version requires sometimes -libverbs BLAS= /home/alex/VASP/GotoBLAS2-1.13_bsd/GotoBLAS2/libgoto2.so # LAPACK, simplest use vasp.5.lib/lapack_double LAPACK= ../vasp.5.lib/lapack_double.o # use the mkl Intel lapack #LAPACK= -lmkl_lapack #----------------------------------------------------------------------- LIB = -L../vasp.5.lib -ldmy \ ../vasp.5.lib/linpack_double.o $(LAPACK) \ $(BLAS) # options for linking, nothing is required (usually) LINK = #----------------------------------------------------------------------- # fft libraries: # VASP.5.2 can use fftw.3.1.X ( http://www.fftw.org) # since this version is faster on P4 machines, we recommend to use it #----------------------------------------------------------------------- #FFT3D = fft3dfurth.o fft3dlib.o # alternatively: fftw.3.1.X is slighly faster and should be used if available #FFT3D = fftw3d.o fft3dlib.o /opt/libs/fftw-3.1.2/lib/libfftw3.a #======================================================================= # MPI section, uncomment the following lines until # general rules and compile lines # presently we recommend OPENMPI, since it seems to offer better # performance than lam or mpich # # !!! Please do not send me any queries on how to install MPI, I will # certainly not answer them !!!! #======================================================================= #----------------------------------------------------------------------- # fortran linker for mpi #----------------------------------------------------------------------- FC=/home/alex/VASP/openmpi/bin/mpif77 FCL=$(FC) #----------------------------------------------------------------------- # additional options for CPP in parallel version (see also above): # NGZhalf charge density reduced in Z direction # wNGZhalf gamma point only reduced in Z direction # scaLAPACK use scaLAPACK (usually slower on 100 Mbit Net) # avoidalloc avoid ALLOCATE if possible # PGF90 work around some for some PGF90 / IFC bugs # CACHE_SIZE 1000 for PII,PIII, 5000 for Athlon, 8000-12000 P4, PD # RPROMU_DGEMV use DGEMV instead of DGEMM in RPRO (depends on used BLAS) # RACCMU_DGEMV use DGEMV instead of DGEMM in RACC (depends on used BLAS) # tbdyn MD package of Tomas Bucko #----------------------------------------------------------------------- #----------------------------------------------------------------------- CPP = $(CPP_) -DMPI -DHOST=\"LinuxIFC\" -DIFC \ -DCACHE_SIZE=4000 -DPGF90 -Davoidalloc -DNGZhalf \ -DMPI_BLOCK=8000 # -DRPROMU_DGEMV -DRACCMU_DGEMV #----------------------------------------------------------------------- # location of SCALAPACK # if you do not use SCALAPACK simply leave that section commented out #----------------------------------------------------------------------- #BLACS=$(HOME)/archives/SCALAPACK/BLACS/ #SCA_=$(HOME)/archives/SCALAPACK/SCALAPACK #SCA= $(SCA_)/libscalapack.a \ # $(BLACS)/LIB/blacsF77init_MPI-LINUX-0.a $(BLACS)/LIB/blacs_MPI-LINUX-0.a $(BLACS)/LIB/blacsF77init_MPI-LINUX-0.a SCA= #----------------------------------------------------------------------- # libraries for mpi #----------------------------------------------------------------------- #LIB = -L../vasp.5.lib -ldmy \ # ../vasp.5.lib/linpack_double.o $(LAPACK) \ # $(SCA) $(BLAS) # FFT: fftmpi.o with fft3dlib of Juergen Furthmueller #FFT3D = fftmpi.o fftmpi_map.o fft3dfurth.o fft3dlib.o # alternatively: fftw.3.1.X is slighly faster and should be used if available FFT3D = fftmpiw.o fftmpi_map.o fftw3d.o fft3dlib.o /home/alex/VASP/fft/lib/libfftw3.a INCS = -I/home/alex/VASP/fft/include #----------------------------------------------------------------------- # general rules and compile lines #----------------------------------------------------------------------- BASIC= symmetry.o symlib.o lattlib.o random.o SOURCE= base.o mpi.o smart_allocate.o xml.o \ constant.o jacobi.o main_mpi.o scala.o \ asa.o lattice.o poscar.o ini.o mgrid.o xclib.o vdw_nl.o xclib_grad.o \ radial.o pseudo.o gridq.o ebs.o \ mkpoints.o wave.o wave_mpi.o wave_high.o \ $(BASIC) nonl.o nonlr.o nonl_high.o dfast.o choleski2.o \ mix.o hamil.o xcgrad.o xcspin.o potex1.o potex2.o \ constrmag.o cl_shift.o relativistic.o LDApU.o \ paw_base.o metagga.o egrad.o pawsym.o pawfock.o pawlhf.o rhfatm.o paw.o \ mkpoints_full.o charge.o Lebedev-Laikov.o stockholder.o dipol.o pot.o \ dos.o elf.o tet.o tetweight.o hamil_rot.o \ steep.o chain.o dyna.o sphpro.o us.o core_rel.o \ aedens.o wavpre.o wavpre_noio.o broyden.o \ dynbr.o rmm-diis.o reader.o writer.o tutor.o xml_writer.o \ brent.o stufak.o fileio.o opergrid.o stepver.o \ chgloc.o fast_aug.o fock.o mkpoints_change.o sym_grad.o \ mymath.o internals.o dynconstr.o dimer_heyden.o dvvtrajectory.o vdwforcefield.o \ hamil_high.o nmr.o pead.o mlwf.o subrot.o subrot_scf.o \ force.o pwlhf.o gw_model.o optreal.o davidson.o david_inner.o \ electron.o rot.o electron_all.o shm.o pardens.o paircorrection.o \ optics.o constr_cell_relax.o stm.o finite_diff.o elpol.o \ hamil_lr.o rmm-diis_lr.o subrot_cluster.o subrot_lr.o \ lr_helper.o hamil_lrf.o elinear_response.o ilinear_response.o \ linear_optics.o linear_response.o \ setlocalpp.o wannier.o electron_OEP.o electron_lhf.o twoelectron4o.o \ ratpol.o screened_2e.o wave_cacher.o chi_base.o wpot.o local_field.o \ ump2.o bse_te.o bse.o acfdt.o chi.o sydmat.o dmft.o \ rmm-diis_mlr.o linear_response_NMR.o vasp: $(SOURCE) $(FFT3D) $(INC) main.o rm -f vasp $(FCL) -o vasp main.o $(SOURCE) $(FFT3D) $(LIB) $(LINK) makeparam: $(SOURCE) $(FFT3D) makeparam.o main.F $(INC) $(FCL) -o makeparam $(LINK) makeparam.o $(SOURCE) $(FFT3D) $(LIB) zgemmtest: zgemmtest.o base.o random.o $(INC) $(FCL) -o zgemmtest $(LINK) zgemmtest.o random.o base.o $(LIB) dgemmtest: dgemmtest.o base.o random.o $(INC) $(FCL) -o dgemmtest $(LINK) dgemmtest.o random.o base.o $(LIB) ffttest: base.o smart_allocate.o mpi.o mgrid.o random.o ffttest.o $(FFT3D) $(INC) $(FCL) -o ffttest $(LINK) ffttest.o mpi.o mgrid.o random.o smart_allocate.o base.o $(FFT3D) $(LIB) kpoints: $(SOURCE) $(FFT3D) makekpoints.o main.F $(INC) $(FCL) -o kpoints $(LINK) makekpoints.o $(SOURCE) $(FFT3D) $(LIB) clean: -rm -f *.g *.f90 *.o *.L *.mod ; touch *.F main.o: main$(SUFFIX) $(FC) $(FFLAGS)$(DEBUG) $(INCS) -c main$(SUFFIX) xcgrad.o: xcgrad$(SUFFIX) $(FC) $(FFLAGS) $(INLINE) $(INCS) -c xcgrad$(SUFFIX) xcspin.o: xcspin$(SUFFIX) $(FC) $(FFLAGS) $(INLINE) $(INCS) -c xcspin$(SUFFIX) makeparam.o: makeparam$(SUFFIX) $(FC) $(FFLAGS)$(DEBUG) $(INCS) -c makeparam$(SUFFIX) makeparam$(SUFFIX): makeparam.F main.F # # MIND: I do not have a full dependency list for the include # and MODULES: here are only the minimal basic dependencies # if one strucuture is changed then touch_dep must be called # with the corresponding name of the structure # base.o: base.inc base.F mgrid.o: mgrid.inc mgrid.F constant.o: constant.inc constant.F lattice.o: lattice.inc lattice.F setex.o: setexm.inc setex.F pseudo.o: pseudo.inc pseudo.F poscar.o: poscar.inc poscar.F mkpoints.o: mkpoints.inc mkpoints.F wave.o: wave.F nonl.o: nonl.inc nonl.F nonlr.o: nonlr.inc nonlr.F $(OBJ_HIGH): $(CPP) $(FC) $(FFLAGS) $(OFLAG_HIGH) $(INCS) -c $*$(SUFFIX) $(OBJ_NOOPT): $(CPP) $(FC) $(FFLAGS) $(INCS) -c $*$(SUFFIX) fft3dlib_f77.o: fft3dlib_f77.F $(CPP) $(F77) $(FFLAGS_F77) -c $*$(SUFFIX) .F.o: $(CPP) $(FC) $(FFLAGS) $(OFLAG) $(INCS) -c $*$(SUFFIX) .F$(SUFFIX): $(CPP) $(SUFFIX).o: $(FC) $(FFLAGS) $(OFLAG) $(INCS) -c $*$(SUFFIX) # special rules #----------------------------------------------------------------------- # these special rules are cummulative (that is once failed # in one compiler version, stays in the list forever) # -tpp5|6|7 P, PII-PIII, PIV # -xW use SIMD (does not pay of on PII, since fft3d uses double prec) # all other options do no affect the code performance since -O1 is used fft3dlib.o : fft3dlib.F $(CPP) $(FC) -FR -names lowercase -O2 -c $*$(SUFFIX) fft3dfurth.o : fft3dfurth.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) fftw3d.o : fftw3d.F $(CPP) $(FC) -FR -names lowercase -O1 $(INCS) -c $*$(SUFFIX) wave_high.o : wave_high.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) radial.o : radial.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) symlib.o : symlib.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) symmetry.o : symmetry.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) wave_mpi.o : wave_mpi.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) wave.o : wave.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) dynbr.o : dynbr.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) asa.o : asa.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) broyden.o : broyden.F $(CPP) $(FC) -FR -names lowercase -O2 -c $*$(SUFFIX) us.o : us.F $(CPP) $(FC) -FR -names lowercase -O1 -c $*$(SUFFIX) LDApU.o : LDApU.F $(CPP) $(FC) -FR -names lowercase -O2 -c $*$(SUFFIX)
内容：什么是VASP、优点、主要输入文件 一 VASP is acomplex package for performing ab-initio quantum-mechanical molecular dynamics (MD) simulations using pseudopotentials or the projector-augmented wave method and aplane wave basis set. The approach imple-mentedin VASP is based on the (?nite-temperature) local-density approximation with the free energy as variationa l quantity and anexact evaluationof the instantaneous electronic ground state at each MD time step. （译文：） 二 计算的体系十余第一行元素和过渡金属；体系计算快；基于Linux/Unix等平行计算；自动对称分析；加速收敛 三 一个最简单的VASP计算体系，应该包括四个输入文件： INCAR（计算细节）；POSCAR（体系坐标）；POTCAR（赝势）；KPONITS（k空间描述）。 吉宗威
标签: VASP
![[转载]VASP与电荷密度差](http://image.sciencenet.cn/album/201210/10/200538sy105p1hqpf88k1s.png.thumb.jpg)
