FHI-gap: A GW Code Based on the All-electron Augmented Plane Wave Method

Table of Contents

• Introduction
• Installation
• Using FHI-gap
• Input parameters
• Postprocessing
• Citing FHI-gap
• Obtain FHI-gap
• Publications using FHI-GAP

Introduction

FHI-gap is an all-electron GW implementation based on the full-potential linearized augmented planewave plus local orbital ((L)APW+lo) method, which handle core, semicore, and valence states on the same footing and therefore allows for a correct treatment of core-valence interaction. It is, therefore, able to handle a wide range of materials, irrespective of their composition. At present, FHI-gap is interfaced to the WIEN2k code, and an implementation into the EXCITING code is in progress.

The main features of FHI-gap are the following:

• The products of Kohn-Sham wavefunctions are expanded into a dual basis, consisting of planewaves in the interstitial region and spherical harmonics inside the atom-centered non-overlapping spheres.
• The numerical Brillouin zone (BZ) integrations are carried out by a generalized tetrahedron method.
• The Gamma point (q=0), at which both the bare and screened Coulomb interaction are singular, is included in the BZ integrations, and the divergences are accurately treated within k.p perturbation theory.
• Frequency integrations are carried out on the imaginary axis.
• The code is parallelized using the message-passing interface library.

An advantage of our code is the capability to explore d- and f-electron systems, materials traditionally categorized as strongly correlated. For such materials, a full-potential all-electron treatment is highly desirable. DFT with LDA or GGA exchange-correlation functionals fails dramatically for many such systems, and so may G₀W₀ carried out on top of LDA or GGA. Here, much of the problem may stem from the LDA/GGA starting point. As a first step towards establishing G₀W₀ for d and f-electron systems, we have implemented G₀W₀ based on LDA+U. This simple and effective approach has been applied to a series of prototype d and f-electron systems and shown to overcome the major shortcomings of LDA/GGA.

Warning: The treatment of metallic systems is still in the experimental state, and has not been well tested.

Return to ToC

Installation

All system-specific setups are defined in make.inc-$SYSTEM, which is used in makefile.
Currently it supports SYSTEM=Linux/SunOS/AIX

You are strongly encouraged to use the INTEL MKL library.

First generate the sequential version by "make seq" and then parallel version with "make para" (only subroutines related to parallelization are recompiled)

In some systems like SunOS, awk does not support -v < variable = value > but nawk does. In that case, one should change alias for awk in gap_init

In more recent versions of WIEN2k (maybe since version 11), core states wave functions are not written to files (case.corewf) by default. It is possible to enable writing case.corewf by setting IPRINT in case.inc manually. In a more straightforward way, one can also modify SRC_lcore/hfsd.f slightly so that case.corewf is always written. For example, in WIEN2k version 13, the following changes can be made:
Original

		   ......
 3999 if(iprint.ne.0) then  
         WRITE(29,2003)  NQN(J),TITRE(J),DEN(j)*2.0d0+shift
         write(29,iform1) (dp(i),i=1,np-npxx)
         write(29,iform1) (dq(i),i=1,np-npxx)
      endif  
	  ...... 
     

1. Define the environmental variable SYSTEM.

2. Change the LAPACK definition in the corresponding make.inc-$SYSTEM.

3. Compiling the code by "make"

4. Modify the alias AWK in gap_init if necessary

5. Necessary or optional changes in WIEN2k source code

6. Now,

    3999 continue 
         WRITE(29,2003)  NQN(J),TITRE(J),DEN(j)*2.0d0+shift
         write(29,iform1) (dp(i),i=1,np-npxx)
         write(29,iform1) (dq(i),i=1,np-npxx)
   	  

    

Return to ToC

Using FHI-gap

gap_init -f casename -d gwdir -nkp number_of_kpoints_for_gw ......

gwdir is the working diretory for the gw calculation to be carried out. It is strongly recommended that gwdir be different from w2kdir. More detailed information about the options for gap_init can be obtained by "gap_init -h"

1. Run a full SCF calculation using WIEN2k
2. Prepare gw input files by running the shell script in the directory (denoated as w2kdir henceforth ) where the SCF calculation has been done
3. Modify gwdir/$file.ingw if necessary, and then run "gap" or "gap-mpi" (using mpirun) in the gwdir.
4. GW calculations are usually very heavy, and therefore parallelization is often necessary to obtain useful results within reasonable time. The way how the FHI-gap code is parallelized is explained in the CPC paper. Here a few guides are given to use the parallelization most efficiently.
   • Always use the parallelization over q-points first since it is almost perfectly parallelized, and this also the default.
   • When the number of q-points (same as the number of k-points n_k) is larger than the number of available processes (n_p), choosing n_p as an integer divider of n_k is helpful to obtain a balanced parallelization.
   • When n_k is smaller than n_p, then the parallization over unoccupied states and the frequencies is also used, for which the number of processes (n_row) is determined by n_p/n_k. Again it is desirable to choose n_row as an integer divider of the number of frequencies (n_w=16 by default).
   • It is possible to change the default parallelization setup by setting an environmental parameter, MPI_NPROC_COL (= n_row, meaning the number of processes per column or the number of rows) before starting gap-mpi (e.g. export MPI_NPROC_COL=...; mpirun -np ...)
5. Some useful information is given in the main output file, case.outgw. Further analysis can be done by using the C-shell script gap_analy (see post-processing)

Return to ToC

Input parameters

A sample of the master input file, case.ingw, reads as follows

Task = "gw"             # Option for task
Restart = F             # Option for whether retarting a previous calculation
nspin = 1               # 1 for spin-unpolarized and 2 for spin-polarized calculations 
ComplexVector = F       # T for for systems without inverse symmetry 
SymVector = F           # whether to use Kohn-Sham eigenvectors in the irreducible Brillouin zone 
barcevtol = 0.0         # tolerance used to reduce the bare Coulomb matrix eigenvectors as the basis set 
                        # for sp systems, barcevtol=0.6 is usually quite safe.
emingw = -2.0           # emingw and emaxgw ( unit --  Ry. ) are used to control the range of bands 
emaxgw =  2.0           # for which GW correction are calculated. Only states 
                        # falling between E_Fermi+emingw and E_Fermi+emaxgw are calculated 
 %SelfEnergy            # option for correlation self-energy 
 2  |  0   |   1        # <npol> |  <iopes>  |  <iopac> 
 %                      # Number of poles ( previous maxexp + 1, valid range: 2.. nomeg/2 ) 
                        # iopes: 0/1/2/3 - without or with itereration > 
                        # iopac:0/1     - Pade's approximation / multipole fitting 
 %BZConv                # BZ convolution options
 "tetra" | "imfreq"
 %
 %FreqGrid               # Frequency grid parameters
 3 | 16 | 0.42 |0.0 | 0  # iopfreq | nomeg | omegmax | omegmin | nomeg_blk
 %                       # iopfreq= 1 (equally spaced), 2 (Gauss-Laguerre) 
                         #  or 3 (double Gauss-Legendre) 
 %MixBasis               # Mixed basis parameters
 1.0 | 3 | 1.E-4         # Q, lmax_MB, lamda
 %
 ºreCoul               # Options to control the calculation of the bare Coulomb matrix 
 2.0  |  1.E-15
 %
 nvel = 8.0              # the number of valence electrons

A detailed description of all available input parameters can be found here.

Return to ToC

Post-processing

The main results of a G₀W₀ calculation by running the FHI-gap code are a set of quasi-particle energies on an equally spaced k-mesh. To calculate densities of states or plot the band structure diagram along high symmetry directions in the BZ, one usually needs QP energies on a fine k-mesh. Calculating QP energies at an arbitrary k is possible, but computationally expensive. Currently, the FHI-gap code uses the Fourier interpolation approach (Pickett et al. Phys. Rev. B 38, 2721 (1988)) to obtain QP energies on k-points that are not included in the original k-mesh. For that purpose, the C-shell scripts gap_analy and gap_gwnvf can be used.

The script gap_gwnvf generates two files, namely case.energy_gw and case.vector_gw, which contain interpolated quasiparticle energies, using data from case.energy and case.vector in the current directory. These QP interpolated energies and vectors can then be used for further analysis or other calculations.

The script gap_analy first call "x lapw1" to calculate KS energies and vectors on a new k-mesh, and then use gap_gwnvf to obtain QP corrected energies and vectors on the same k-mesh, from which either DOS ("x tetra") or band structure diagrams ("x spaghetti") can be obtained.

Note that both gap_analy and gap_gwnvf should be run in the w2kdir. More information about these two scripts can be obtained by running them with a command-line argument "-h".

To use gap_analy properly, some changes have to be made to x_lapw and SRC_spaghetti/spag.f in the WIEN2k package so that spaghetti reads band energies from case.energy instead of case.output1. The part related to spaghetti in x_lapw is changed as followings.

cat << theend > $def
5, '$file.insp',         'old',    'formatted',0
6, '$file.outputsp$updn',     'unknown','formatted',0
10,'$file.spaghetti${updn}_ene','unknown','formatted',0
11,'$file.spaghetti${updn}_ps', 'unknown','formatted',0
20,'$file.struct',       'old',    'formatted',0
30,'$file.irrep$so$updn',          'unknown','formatted',0
40,'$file.bands${updn}.agr',       'unknown','formatted',0
theend
echo "8, '$file.energy$so$updn',  'old','formatted',0" >>$def
if($?qtl) then
  echo "9, '$file.qtl$updn',          'unknown','formatted',0" >>$def
endif

The modified x_lapw and SRC_spaghetti/spag.f are included in the directory src_misc. Note: the included x_lapw is based on that of WIEN2k version 11.0; if you are using a different version of WIEN2k, you should modify x_lapw manually so that other functions of x_lapw are not affected.

Return to ToC

Citing FHI-gap

When publishing the work that has used the FHI-gap code, please cite the following paper

Hong Jiang, Ricardo I. Gomez-Abal, Xinzheng Li, Christian Meisenbichler, Claudia Ambrosch-Draxl, and Matthias Scheffler, FHI-gap: a GW code based on the All-electron augmented plane wave method, Computer Phys. Commun.,184, 348 (2013).

Return to ToC

Obtain FHI-gap

The FHI-gap code will be made available to the WIEN2k users free of charge. The stable version of the code will be available from the WIEN2k website. The latest version can be downloaded here.

Return to ToC

Publications using FHI-gap

1. Ferenc Karsai, Paul Tiwald, Robert Laskowski, Fabien Tran, David Koller, Stefanie Gr?fe, Joachim Burgd?rfer, Ludger Wirtz, and Peter Blaha, F center in lithium fluoride revisited: Comparison of solid-state physics and quantum-chemistry approaches, Phys. Rev. B 89, 125429 (2014).
2. Philipp Hansmann, Loig Vaugier, Hong Jiang and Silke Biermann, What about U on surfaces? Extended Hubbard models for adatom systems from first principles, J. Phys.: Condens. Matter 25 094005 (2013).
3. Hong Jiang, Ricardo I. Gomez-Abal, Xinzheng Li, Christian Meisenbichler, Claudia Ambrosch-Draxl, and Matthias Scheffler, FHI-gap: a GW code based on the All-electron augmented plane wave method, Computer Phys. Commun.,184, 348 (2013).
4. Loig Vaugier, Hong Jiang and Silke Biermann, Hubbard U and Hund's exchange J in transition metal oxides: screening vs localization tends from constrained random phase approximation, Phys. Rev. B 86, 165105(2012).
5. Hong Jiang, Patrick Rinke and Matthias Scheffler, Electronic properties of lanthanide oxides from the GW perspective, Phys. Rev. B, 86, 125115 (2012).
6. Hong Jiang, Electronic Band Structures of Molybdenum and Tungsten Dichalcogenides by the GW Approach J. Phys. Chem. C 116,7664-7671 (2012).
7. Xin-Zheng Li, Ricardo Gomez-Abal, Hong Jiang, Claudia Ambrosch-Draxl, and Matthias Scheffler, Impact of widely used approximations to the G0W0 method: an all-electron perspective, New J. Phys. 14, 023006 (2012).
8. Huihui Wang, Feng Wu and Hong Jiang, Electronic band structures of ATaO3 (A=Li, Na and K) from first-principles many-body theory, J. Phys. Chem. C 115, 16180 (2011).
9. Hong Jiang, Structural and electronic properties of ZrX2 and HfX2(X=S, and Se) from first principles calculations J. Chem. Phys. 134, 204705 (2011)
10. Hong Jiang , Ricardo I. Gomez-Abal, Patrick Rinke, and Matthias Scheffler, First-principles modeling of localized d states with the GW@LDA+U approach, Phys. Rev. B 82, 045108 (2010).
11. Hong Jiang, Ricardo I. Gomez-Abal, Patrick Rinke, and Matthias Scheffler, Quasiparticle Electronic Band Structure of Zirconia and Hafnia Polymorphs, Phys. Rev. B 81, 085119(2010).
12. Hong Jiang, The GW Method: Basic Principles, Latest Developments and Its Applications for d-and f-Electron Systems (Invited Review), Acta Phys.-Chim. Sin.,26, 1017(2010).
13. Hong Jiang , Ricardo I. Gomez-Abal, Patrick Rinke and Matthias Scheffler, Localized and itinerant states in lanthanide oxides united by GW@LDA+U , Phys. Rev. Lett. 102, 126403 (2009). pdf
14. R. Gómez-Abal, X. Li, M. Scheffler, and C. Ambrosch-Draxl,Influence of the core-valence interaction and of the pseudopotential approximation on the electron self-energy in semiconductors. Phys. Rev. Lett. 101, 106404 (2008).

Return to ToC