||
关注:
1) GW方法的物理含义;
2) VASP 计算过程及参数设置;
4) 手册是最好的老师,很多问题手册上都有答案:
http://cms.mpi.univie.ac.at/vasp/vasp/Contents.html
1、To obtain band structure plots from the GW method, first you have to perform a DFT calculation.
The INCAR file I used for DFT on Si looks like this:
ISTART = 0
ISYM = 0
ICHARG = 2
ISMEAR = 0
SIGMA = 0.1
LORBIT = 11
NSW = 200
IBRION = 2
ISIF = 3
EDIFF = 0.0001
EDIFFG = -0.01
KPOINTS file for the run is:
Monkhorst Pack
0
Gamma
4 4 4
0 0 0
POSCAR file is:
Si 5.38936000000000
0.5000000000000000 0.5000000000000000 0.0000000000000000
0.0000000000000000 0.5000000000000000 0.5000000000000000
0.5000000000000000 0.0000000000000000 0.5000000000000000
2
Direct
0.0000000000000000 0.0000000000000000 0.0000000000000000
0.2500000000000000 0.2500000000000000 0.2500000000000000
POTCAR file is supplied with VASP 5.2, remember to use paw pseudopotentials so that in the band structure calculation the Wigner-Seitz radius for Si need not be directly specified in INCAR.
2.After this run completes you are ready to perform the HF-type (HSE06) calculation. Change
the INCAR file to:
ISTART=1
ICHARG=2
EDIFF=0.0001
EDIFFG=-0.01
ENCUT=400
ENAUG=800
LREAL=.FALSE.
LWAVE=.TURE.
NELM=200
NSE=0
IBRION=-1
LMAXMIX=4
ISMEAR=0
SIGMA=0.1
NSIM=4
IALGO=48
SIYM=0
LHFCALC=.TRUE.
HFSCREEN=0.2
ALGO=Damped
TIME=0.4
ENCUTFOCK=0NKRED=2
AEXX=0.25
keep the other files constant for this run (e.g. KPOINTS, POSCAR, POTCAR need not be changed)
3. once this run has completed, add to the bottom of the INCAR file from the HSE06 calculation:
LOPTICS=.TRUE.
This step calculats the dielectric matrix
4. Once this run has completed ,change
ALGO=Damped
to
ALGO=GWo
and
ICHARG=2
to ICHARG=11
and add:
NOMEGA=50
ENCUTGW=100
LORBIT=11
and change the KPOINTS file to:
k-points along high symmetry lines
20 !20 intersections
Line-mode
rec
0 0 0 !gamma
0.5 0.0 0.0 !X
.....
These values are the points of high symmetry along the lines of high symmetry within the first Brillouin Zone of the Si unit
cell structures in reciprocal space. And that's it !
Once complete you can plot the PROCAR, EIGENVAL, and DOSCAR files from the output of this run as described in the previous posts.
The band structure plot produce from the plotting utility is shown below.
附手册解读 :Frequency dependent GW calculations
http://cms.mpi.univie.ac.at/vasp/vasp/Frequency_dependent_GW_calculations.html
Frequency dependent GW calculations
Available as of VASP.5.X. For details on the implementation and use of the GW routines we recommend the following references: Ref. [111,112,113,114].
ODDONLYGW and EVENONLYGW and NKRED: reducing the -grid for the response functions
Using the GW routines for the determination of frequency dependent dielectric matrix
112M. Shishkin and G. Kresse, ``Self-consistent calculations for semiconductors and insulators'', Phys. Rev. B 75, 235102 (2007).
113F. Fuchs, J. Furthmüller, F. Bechstedt, M. Shishkin, and G. Kresse, ``Quasiparticle band structure based on a generalized Kohn-Sham scheme'', Phys. Rev. B 76, 115109-1-8 (2007).
114M. Shishkin, M. Marsman, and G. Kresse, ``Accurate quasiparticle spectra from self-consistent GW with vertex corrections'', Phys. Rev. Lett. 99, 246403 (2007).
115F. Bruneval, N. Vast, and L. Reining, Phys. Rev. B 74, 45102 (2006).
(1)Recipe for selfconsistent GW calculations
Selfconsistent GW calculations are only supported in a QP picture. As for GW0, it is possible to update the eigenvalues only (ALGO=GW), or the eigenvalues and one-electron orbitals (ALGO=scGW). In all cases, a quasiparticle picture is maintained, i.e. satellite peaks (shake ups and shake downs) can not be accounted for in the selfconsistency cycle. Selfconsistent GW calculations can be performed by simply repeatedly calling VASP using:
System = Si NBANDS = 96 ISMEAR = 0 ; SIGMA = 0.05 ALGO = GW # eigenvalues only or alternatively ALGO = scGW # eigenvalues and one electron orbitalsFor scGW0 or scGW non diagonal terms in the Hamiltonian are accounted for, e.g. the linearized QP equation is diagonalized, and the one electron orbitals are updated.[114]
Alternatively (and preferably), the user can specify an electronic iteration counter using NELM:
System = Si NBANDS = 96 ISMEAR = 0 ; SIGMA = 0.05 ALGO = GW # or ALGO = scGW NELM = 3In this case, the one electron energies (=QP energies) are updated 3 times (starting from the DFT eigenvalues) in both G and W. For ALGO = scGW (or ALGO = QPGW in VASP.5.2.13), the one electron energies and one electron orbitals are updated 3 times.[114] As for ALGO = scGW0, the ``static'' COHSEX approximation can be selected by setting NOMEGA = 1 [115]. To improve convergence to the groundstate, the charge density is mixed using a Kerker type mixing starting with VASP.5.2.13 (see Sec. 6.49). The mixing parameters AMIX, BMIX, AMIX_MAG, BMIX_MAG, AMIN can be adjusted, if convergence problems are encountered.
Alternatively the mixing may be switched off by setting IMIX=0 and controlling the step width for the orbitals using the parameter TIME (which defaults to 0.4). This selects a fairly sophisticated damped MD algorithm, that is also used for DFT methods when ALGO is set the ``Damped''. In general, this method is more reliable for metals and materials with strong charge sloshing.
附--网言摘录
Recipe for G0W0 calculations
GW calculations always require the calculation of a standard DFT WAVECAR file in an initial step,
using for instance the following INCAR file:
System = Si
NBANDS = 96
ISMEAR = 0 ; SIGMA = 0.05 ! small sigma is required to avoid partial occupancies
Note, that the a significant number of empty bands is required for GW calculations, so that it might be better to perform the calculations in two steps: first a standard grounstate calculation with few unoccupied orbitals only,
System = Si groundstate occupied orbitals
ISMEAR = 0 ; SIGMA = 0.05
! small sigma is required to avoid partial occupancies EDIFF = 1E-8 ! required tight tolerance for groundstate orbitalsand second a calculation of a large number of unoccupied orbitals
System = Si unoccupied orbitals
ALGO = Exact ! use exact diagonalization of the Hamiltonian
NELM = 1 ! since we are already converged stop after one step
NBANDS = 96
ISMEAR = 0 ; SIGMA = 0.05 ! small sigma is required to avoid partial occupancies
LOPTICS = .TRUE.Furthermore note that the flag
LOPTICS=.TRUE. is required in order to write the file WAVEDER, which contains the derivative of the orbitals with respect to the k-points
more precissely the matrix [compare (6.76)]
Calculation of this matrix requires the knowledge of the Hamiltonian, and therefore needs to be done in the preparatory DFT or hybrid functional run. The actual GW calculations are performed in a second step using an INCAR file such as (it is convenient to add a single line):System = Si NBANDS = 96 ISMEAR = 0 ; SIGMA = 0.05 LOPTICS = .TRUE. ALGO = GW0 ; NOMEGA = 50The head and wings of the dielectric matrix are constructed using k.p perturbation theory (this requires the file WAVEDER). In the present release the interaction between the core and the valence electrons is always treated on the Hartree Fock level [104].
For hybrid functionals, the three step procedure will accordingly involve the following INCAR files. In the first two steps, converged HSE03 orbitals are determined (usually HSE03 calculations should be preceeded by standard DFT calculations, we have not documented this step here, see Sec. 6.69.11):
System = Si groundstate occupied orbitals ISMEAR = 0 ; SIGMA = 0.05 ALGO = Damped ; TIME = 0.5 ! or ALGO = Conjugate LHFCALC = .TRUE. ; AEXX = 0.25 ; HFSCREEN = 0.3 EDIFF = 1E-6 ! required tight tolerance for groundstate orbitalsSecond determine the HSE03 orbitals for unoccupied states:System = Si unoccupied orbitals NBANDS = 96 ALGO = Exact ! perform exact diagonalization NELM = 1 ! since we are already converged stop after one step ISMEAR = 0 ; SIGMA = 0.05 LHFCALC = .TRUE. ; AEXX = 0.25 ; HFSCREEN = 0.3 LOPTICS = .TRUE.As before, in the GW step, the head and the wings of the response matrix are determined by reading the required data from the WAVEDER file.System = Si NBANDS = 96 ISMEAR = 0 ; SIGMA = 0.05 ALGO = GW0 ; NOMEGA = 50Convergence with respect to the number of empty bands NBANDS and with respect to the number of frequencies NOMEGA must be checked carfully.
附. 参数说明
ENCUTFOCK: FFT grid in the HF related routines
ENCUTFOCK = real number (energy cutoff determining the FFT grids in the HF related routines)
default: none
The ENCUTFOCK parameter controls the FFT grid for the HF routines. The only sensible value for ENCUTFOCK is
ENCUTFOCK=0. This implies that the smallest possible FFT grid, which just encloses the cutoff sphere corresponding to the
plane wave cutoff, is used. This accelerates the calculations by roughly a factor two to three, but causes slight changes in the
total energies and a small noise in the calculated forces. The FFT grid used internally in the Hartree Fock routines is written
to the OUTCAR file. Simply search for lines starting with
FFT grid for exact exchange (Hartree Fock)
In many cases, a sensible approach is to determine the electronic and ionic groundstate using ENCUTFOCK = 0, and to
make one final total energy calculation without the flag ENCUTFOCK.
网络摘录:
http://emuch.net/html/201312/6800622.html
NOMEGA, NOMEGAR number of frequency points
NOMEGA= [integer] (number of frequency points)
NOMEGAR= [integer] (number of frequency points along real axis)
Default: | ||
NOMEGA | =50 | for calculations |
谢谢你的回帖。我依照你的建议再做一次Si的测试。
我在计算过程中注意了CPU的要求,NBAND和NOMEGA和CPU核心数目是一样的。Si的赝势,Si_GW,ENMAX只有245.345 eV,我使用的300 eV应该已经足够了,不过我会调大再看看结果。
现在看来,您提到的空带数,K网格,频率采样密度是最可能的问题原因。可惜,这三个参数都和计算量有直接的关系。现在的参数设置下,GW的计算已经非常慢了,实际考察的复杂结构所需的时间估计也会很长。我看了Abinit GW计算部分的介绍,使用Eective Energy Technique (EET) 可以避免使用上百条的空带也能得到很好的结果。如果VASP速度实在不可接受,最后也考虑试一下Abinit
对于Si这样的元素,可以用Plasmon-Pole模型来算介电函数,这样就不用针对每个频率点都算一次介电函数了,可以节约不少计算量。
但是vasp只能算全频率GW,因此计算量很大。abinit支持4种PPA模型,对这四种模型的测试可以看下PRB 88, 125205(2013)这篇文章。
重新优化结构计算了Si的能隙,NBAND=240,ENCUT=300,KPOINTS 8*8*8,NOMEGA = 80, GW计算出来的能隙仍然是在1.8 eV。我搜索到有些Si GW计算的结构,即使使用高得多的收敛参数,仍然能得到接近实验值的结果。我的计算结构差距太大了,现在很怀疑是不是哪里做错了。
这是在GW计算上一步中计算介电矩阵的INCAR。计算介电矩阵这一步使用了静态计算的结果。
PREC = Accurate
ISTART = 1
ICHARG = 2
ISPIN = 0
ENCUT = 400
ISYM = 2
GGA = PS
NBANDS = 240
!GW0
LOPTICS = .TRUE.
NEDOS = 2000
LPEAD = .TRUE.
!Electronic Relaxation
ALGO = Exact
EDIFF = 1E-05 stopping-criterion for ELM
! LREAL = AUTO real-space projection
!Ionic relaxation
EDIFFG = -1E-03 stopping-criterion for IOM
ADDGRID= .TRUE. reduce the noise in the forces
!DOS related values:
ISMEAR = 0 gaussian smear 0 tetrahedron method -5
SIGMA = 0.05
# NOMEGA specifies the number of frequency grid points
# ENCUTGW energy cutoff for response function ENCUTGW= [real] (energy cutoff for response function)
# Default: ENCUTGW=ENCUT
ALGO for response functions and GW calculations
ALGO= CHI | GW0 | GW | scGW | scGW0 | QPGW | QPGW0
For ALGO=GW and ALGO=GW0 the orbitals of the previous groundstate calculations are maintained, and single shot GWcalculations are performed. If NELM is set as well, several iterations are performed, and the eigenvalues are updated in the calculation of G (ALGO=GW0) or W and G (ALGO=GW). A full update of the orbitals can be performed by specifying ALGO=SCGW and ALGO=SCGW0 (QPGW and QPGW0 are synonymous to these setting, and available in VASP.5.2.13). In the former case, the orbitals and eigenvalues are updated in G and W, whereas in the latter case the orbitals and eigenvalues are only updated in G. Convergence of the eingevalues with respect to ENCUTGW (see Sec. 6.73.7 and the number of unoccupied bands NBANDS is usually fairly slow and should be checked carefully.
We strongly recommend to read the following literature before performing calculations using VASP [111,112,113,114].
## 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
Archiver|手机版|科学网 ( 京ICP备07017567号-12 )
GMT+8, 2024-9-15 21:20
Powered by ScienceNet.cn
Copyright © 2007- 中国科学报社