# The GWinput file

This guide documents the structure of the *GWinput* file used by the *GW* package. The *GW* suite does not use the *ctrl* file, but reads tags from *GWinput*. You can, however, create a template *GWinput* from the *ctrl* file, as explained for example, in this tutorial. This is particularly convenient because *GWinput* is not easy to read, or to construct by hand.

*Note*: in future *GWinput* will be subsumed into the *ctrl* file.

This documentation denotes Ref. I as a reference to PRB76 165106 (2007).

### Table of Contents

- How input is organized
- The general category
- The PRODUCT_BASIS category
- The QPNT category
- The PBASMAX category
- The QforEPS and QforEPSL categories

### How input is organized

File *GWinput* is the main put file for the *GW* package. Its structure has the tag-value style roughly similar to the *ctrl* file, but it is much less sophisticated. All but the first category is delimited html-style, e.g.

```
<PRODUCT_BASIS>
...
</PRODUCT_BASIS>
```

Categories must be ordered as follows:

```
(general, no name or delimiter)
<PRODUCT_BASIS>
...
</PRODUCT_BASIS>
<QPNT>
...
</QPNT>
<PBASMAX>
...
</PBASMAX>
<QforEPS>
...
</QforEPS>
```

Some lines within a category are comments; some require a fixed order though the first (general) category do not. The user is strongly advised to make a *GWinput* template using `lmfgwd --job=-1`

, and edit it.

**<PBASMAX>** is optional; **<QforEPS>** is needed only to compute macroscopic the dielectric function or magnetic susceptibility at specified **k** points.

### The general category

The general category has the structure

```
tag data
```

Some tags are optional. `data`

may be a string, a real number, an integer or a logical, and sometimes a vector of numbers. Lines beginning with, or text following, a **!** are comments, and ignored by the reader. Below is a sample general category. The first column was added to label lines; they are described following the box.

```
!Verbose 41 ! 0-->default; 100-->debug
GWversion 12 ! 0-->Sep12 compatible.
!Q0P_Choice 1 ! 1(default):Near q->0 limit 2:1/q^2 average in Gamma region
!CoreOrth off ! off --> Do not orthogonalize core to valence (default)
!EXonly .15 ! for exchange-only calculations
!EIBZmode off ! turn off to suppress use of symmetrization in calculating polarization
!TimeReversal off ! when time-reversal symmetry is broken
KeepEigen on ! keep eigenfunctions in memory
KeepPPOVL on ! keep PPOVL in memory
!Chi_RegQbz on ! Offset Gamma point mesh for chi. on => no offset
BZmesh 1 ! Offset Gamma point mesh for Sigma=iGW
WgtQ0P 0.01 ! Weight used when BZmesh is 2
NormChk 0 ! 1,2 writes norm check files (diagonal or full)
n1n2n3 4 4 4 ! for GW BZ mesh
QpGcut_psi 2.7 ! |q+G| cutoff for eigenfunction
QpGcut_cou 2.2 ! |q+G| cutoff for coulomb int.
unit_2pioa off ! off --> units of 2 preceding Gcut are a.u.; on--> units are 2*pi/alat
alpha_OffG 1 ! (a.u.) parameter (auxiliary function, offset-Gamma method)
!nband_chi0 999 ! nband cutoff for chi0 (Optional)
!emax_chi0 999. ! emax cutoff for chi0, Ry (Optional)
!nband_sigm 999 ! nband cutoff for Sigma (Optional)
emax_sigm 2.5 ! Energy cutoff for Sigma, Ry (Optional)
dw 0.01 ! mesh spacing along Real axis (Hartree)
omg_c 0.1 ! mesh spacing increases linearly with w: becomes 2*dw at omg_c
iSigMode 3 ! QSGW mode switch (QSGW only)
niw 6 ! # freq. on Im axis; used for integration to make Sigma_c
delta -1e-4 ! delta-function broadening for calc. x0, a.u.. delta<0 => tetrahedron
deltaw 0.02 ! width in finite diff for sigma energy derivative, a.u.
esmr 3e-3 ! Broadening in the poles of G(LDA) (hsfp0)
GaussSmear on ! on => broadening of poles in G(LDA) by Gaussian
!mixbeta .25 ! mixing of input, output sigma for self-consistency
!mixpriorit 2
UseBSEW 0 ! include BSE (ladder diagram) contributions to W
NumCondStates 8 ! (BSE only) number of conduction states
Qvec_bse 1 1 1 ! (BSE only) q-vector for single-shot BSE dielectric function
```

Tags by referenced by line:

^{1 } Specifying the number of divisions for a uniform mesh in Brillouin zone. Requires 3 integers: it is the analog of **BZ_NKABC** in the *ctrl* file.

^{2 } Controls changes in functionality as the *GW* code has evolved. **GWversion** consists of a compound of digits.

- 1s digit
- include local field correction when calculating W(q=0); same as LFC@Gamma, Sep12 version of aux function
- Same as 1, but Aug 14 auxiliary function

- 10s digit
- number of intermediate states are computed relative to , not (this is a bug fix)
- alternative analytic auxiliary functions to regularize singularity at q-> (
**q0irre.F**)

- 100s digit
- alternative estimate for average SE high lying states (
**hqpe.sc.m.F**) - rotate sigm in basis of LDA eigenfunctions (
**hqpe.sc.m.F**)

- alternative estimate for average SE high lying states (
- 1000s digit
- Use high accuracy but slow Bessel functions (not usually necessary)

**GWversion=12** is the recommended choice. Use **GWversion 0** for best compatibility with the Sep 2012 version of code.

^{3 } Controls how the offset gamma points are determined. **Q0P_Choice** should be one of:

- 6 points along the reciprocal lattice vectors ( for each vector). This the default.
- 6 points along Cartesian axes (± for each vector)
- 4 points in Qlat(1) and Qlat(2) (± for each vector)
- 2 points in Qlat(3) (± for each vector)

These points may be reduced by point group symmetry. Offset gamma points are generated by **qg4gw** and stored in file *Q0P*.

For slabs, **Q0Pchoice 2** may be better; but take care has this has not been fully analysed. Regardless, it is problematic to use unbalanced k points for anisotropic unit cells. see Computer Physics Comm. 176, 1 (2007).

One-dimensional chains seem to be well represented by **Q0P_Choice 3** together with **BZmesh 2**.

^{4 } If **on** cores are orthogonalized to valence partial waves and . Enforcing orthogonalization condition is needed to ensure the correct behavior in *ε* as **q**→0. However, orthogonalization may excessively deform shallow core orbitals, so usually it is safer not to use this switch. It is useful, however, to check how orthogonalization affects results. If it does, you can alternatively include these shallow cores in the secular matrix using local orbitals, which is very reliable and has no problems with orthogonalization.

^{5 } Use **EXonly** for exchange-only calculations. The value following **EXonly** tells the code how much Hartree Fock character to admix with LDA. Thus use **EXonly 1** for Hartree Fock.

^{6 } By default, symmetry operations are used in calculating the polarizability. We follow a procedure similar to that of Friedrich *et al* as described in Phys. Rev. B81, 125102. Use **EIBZmode off** to suppress use of symmetrization.

^{7 } Usually the polarizability assumes time-reversal symmetry Use **TimeReversal off** if it is broken.

^{8,9} To reduce disk access, eignvectors and the overlap matrix of the product basis are read and retained in memory by default. You can cause the *GW* code to read them from disk every time one is needed. This conserves memory but increases the (slow) disk access.

Use **KeepEigen off** to always read eigenfunctions from disk;

Use **KeepPPOVL off** to always read **PPOVL** from disk.

^{10} Use this switch to stagger the **k** mesh in calculating the dielectric function. It should be **on** (the default) or **off**. **off** offsets the mesh to avoid the point, analogous to what **BZ_BZJOB** does for **lmf** and related codes; see *ctrl* file input (but note that **BZJOB** can accept up to three switches, **Chi_RegQbz** only one).

At times the static limit of macroscopic dielectric function *ε _{M}*(

*ω*→0) can converge very slowly with the number of

*k*divisions

*n*, so that getting the converged result can be painful. It may happen (and often does) that

_{k}*ε*(

_{M}*ω*→0) converges with

*n*from below with the staggered mesh and from above with the usual one. In such a case the two meshes bound

_{k}*ε*(0) for finite

_{M}*n*, and should converge in the 1/

_{k}*n*→0 limit. By plotting

_{k}*ε*(0) as a function of 1/

_{M}*n*for both meshes at several values of

_{k}*n*, you can extrapolate each of them to 1/

_{k}*n*→0. If the two meshes extrapolate to the same point, it can be taken as a reliable value for

_{k}*ε*.

_{M}An instance of this is NH_{4}PbI_{3}. The Figure shows *ε _{M}*(0) plotted for the two meshes. Even while

*ε*(0) is not well converged for either mesh (especially the regular one), extrapolating to 1/

_{M}*n*=0 gives much confidence in the converged result.

_{k}See also the Brillouin integration page.

^{11} affects the **k** mesh in calculating the self-energy, analogous to **BZ_BZJOB** in the *ctrl* file. It should be **1** (the default) or **2**: **2** offsets the mesh to avoid the point. **1** is recommended.

^{12} The total weight for the offset gamma method when **BZmesh 2** (and is effective only for that case). More precisely it is fraction of the weight given to a regular mesh points **(n1 n2n3)^{−1}**. In principle, the final result should not depend on

**WgtQ0p**. It is usually OK to take 0.01 (default), but there is a little dependence and you can take smaller values. When

**WgtQ0p**is very small, e.g. 10

^{−5}, the QP levels become sensitive to the errors in normalization of the eigenfunction. You can improve the normalization with a larger

**QpGcut_psi**. Normalization is printed out in

*normchk.dia*.

^{13} Set to **1** or **2** causes the code to write norm check files (diagonal or full).

^{14} Number of *k* divisions along each reciprocal lattice vectors.

^{15} **G** cutoff (actually **q+G** cutoff) in the plane wave expansion of the eigenfunctions.

^{16} **G** cutoff (actually **q+G** cutoff) in the plane wave expansion of two particle objects such as and .

^{17} specifies whether **QpGcut_psi** and are in **QpGcut_cou** and are in atomic units (**off**), or in units of 2*π*/*a* (**on**).

^{18} corresponds to in Eq. 48 of Ref. I. **alpha_offG 1d0** corresponds to offset points close to zero, and is usually a good choice. You can check convergence in ***n1n2n3**, varying this number.

^{19} if present, limits the number of bands used in computing *ε*(*ω*).

^{20} if present, limits the energy in computing *ε*.

^{21} if present, limits the number of bands used in computing Σ.

^{22} if present, limits the energy in computing Σ.

^{23,24} **dw** and **omg_c** fix the frequency mesh for real-axis integration of the polarizability, Eq. 32 of Ref. I, in atomic units. For small the mesh points are linearly spaced, separated by **dw**. If **omg_c** is not present, mesh points are linearly spaced on the entire grid. Otherwise the spacing gradually increases with , becoming twice **dw** at . Points on the mesh are interpolated to determine . Point is given by

is approximately **dw**/2.

^{25} Defines kind of quasiparticle self-consistency (how static Σ is extracted from the dynamical one).

**iSigMode=3**: Use (Eq. 10, mode-A in Ref. I).

**iSigMode=1**: Use (Eq. 11, mode-B in Ref. I).

**iSigMode=5**: Use (Eigenvalue-only self-consistency).

There is no “right” choice but **3** has the best formal justification.

^{26} Number of integration points along the imaginary axis to get , as described in II-F in Ref.I. The integration points are , where is the set of points for a Legendre quadrature in the interval [0,1]. (There is also special treatment of the sharply peaked contribution near .) Convergence is rapid in **niw**. Tests show **niw=6** is sufficient to reach 0.01 eV accuracy in QP levels in Si.

^{27} Small imaginary part in Eq. 32 of Ref. I, in Hartree. turns on tetrahedron integration.

^{28} When is required, e.g. for the *Z* factor, it is computed numerically with a finite difference. **deltaw** is the energy spread used in the finite difference.

^{29,30} Broadening of the poles in the self-energy makers **hsfp0** and **hsfp0_sc**. Poles of the Green function are smeared out by a width **esmr**. If **GaussSmear** is **on**, the smearing has the shape of a Gaussian function; otherwise a rectangular smearing is used. For insulators the value should be smaller than the bandgap (you can use small numbers like 0.001 or even smaller). For metals, the number should be a little larger. Some experimenting is needed, as it depends on the case. For Fe, **esmr=1e-2** seems to be a good value.

^{31,32} In QS*GW*, self-energies are mixed with prior iterations to stabilize convergence. Provided a file *sigm* is available on disk (taken to be , i.e. the input ), **hqpe_sc** assembles results from the self-energy maker **hsfp0_sc** to make a new , and writes to *sigm* a new estimate for . It takes the linear combination

If **mixpriorit** is nonzero, **hqpe_sc** looks for pairs from prior iterations in file *mixsigma*. If they are available, Anderson mixing is used to estimate a new from a linear combination of and pairs from **mixpriorit** prior iterations (or as many as available, if fewer than **mixpriorit**). If **mixpriorit** is zero, is written. Anderson mixing used is the same as **lmf** uses to mix input and output charge densities, and some description can be found in the input file manual. Note in particular the discussion of parameters **tj** found the **hqpe_sc** output, file *lqpe*.

^{33,34,35} These are switches that cause the *GW* script to and add ladder diagram contributions to the screened coulomb interaction *W*, after *W* is calculated in the RPA. To add this correction, set **UseBSEW** to a nonzero value. The BSE branch also requires input for the number of conduction band states, **NumCondStates**.

### The PRODUCT_BASIS category

An instance of this category, taken from NiO, is shown below. (The first column was added to label lines; they are described following the box.) NiO has two Ni atoms and two O atoms. Thus the information for sites 1 and 2 should be identical, as and also for sites 3 and 4.

*Note:* when generating a template for *GWinput*, **lmfgwd** automatically generates equivalent product basis information for equivalent sites. Whenever you edit *GWinput* by hand it is a good idea to check for consistency. **lmfgwd** can perform a number of checks for you, including checking equivalency of equivalent sites. To run a check, do `lmfgwd --job=-2`

.

```
<PRODUCT_BASIS> ! Product basis block
tolerance = minimum eigenvalue in PB overlap to reduce linear dependency
3e-4 3e-4 1e-3
lcutmx(atom) = l-cutoff for the product basis
4 4 3 3
atom l nnvv nnc ! nnvv: num. radial functions (valence) for augmentation-waves. nnc = num. for core.
1 0 2 3
1 1 3 1
1 2 3 0
1 3 2 0
1 4 2 0
...
3 0 2 1
3 1 2 0
3 2 2 0
3 3 2 0
3 4 2 0
...
atom l n occ unocc :Valence(1=yes, 0=no)
1 0 1 1 1 ! 4S_p *
1 0 2 0 0 ! 4S_d
1 1 1 1 1 ! 4P_p
1 1 2 1 1 ! 4P_d
1 1 3 1 1 ! 3P_l
1 2 1 1 1 ! 3D_p
1 2 2 0 1 ! 3D_d
1 2 3 0 1 ! 4D_l
1 3 1 0 1 ! 4F_p
1 3 2 0 0 ! 4F_d
1 4 1 0 0 ! 5g_p
1 4 2 0 0 ! 5g_d
...
3 0 1 1 1 ! 2S_p *
3 0 2 0 0 ! 2S_d
3 1 1 1 1 ! 2P_p
3 1 2 0 0 ! 2P_d
3 2 1 1 1 ! 3D_p
3 2 2 0 0 ! 3D_d
3 3 1 0 1 ! 4F_p
3 3 2 0 0 ! 4F_d
3 4 1 0 0 ! 5g_p
3 4 2 0 0 ! 5g_d
...
atom l n occ unocc ForX0 ForSxc :CoreState(1=yes, 0=no)
1 0 1 0 0 0 0 ! 1S *
1 0 2 0 0 0 0 ! 2S
1 0 3 1 0 1 1 ! 3S
1 1 1 0 0 0 0 ! 2P
...
3 0 1 0 0 0 0 ! 1S *
...
```

^{2 } This line is used as a marker to identify the next line (keyword is **tolerance**), so it should be present.

^{3 } A list of floating-point numbers. The size of list should range between **1** and **2×lmax+1** numbers, where **lmax** is the largest integer in line 5. The code will populate any numbers omitted with the last number you supply.

These numbers refer to tolerances in the overlap matrix of product basis functions. Product basis functions of different *l* are automatically orthogonal, so each *l* can be treated independently of the others. The overlap matrix constructed from product basis for a particular site and *l* is diagonalized to form an orthonormal set. Orthogonalized functions with small eigenvalues add little to the hilbert space of functions. Those with eigenvalues below **tolerance** are discarded, both reducing the size of the basis and avoiding numerical issues stemming from linear dependencies.

For small *l* it is safer to use small tolerances; but for higher *l* tolerances can be coarser for the same accuracy. Computational time is more sensitive to the higher *l* tolerances because they have more functions of different *m* per radial function.

^{4 } This line is used as a marker to identify the next line (keyword is **lcutmx(atom)**), so it should be present.

^{5 } specifies the *l* cutoff **lmax** for each site. Product functions are expanded to *l* = 2×**lmax**. You must supply one **lmax** for each site. **lmax=3** is usually good enough for light *sp* elements, **lmax=4** for transition metals. But as the nucleus becomes heavier and larger, **lmax** needs to be made somewhat larger. We have not tested the convergence for *f* shell elements, but **lmax=6** is almost certainly safe. Smaller **lmax** reduces the computational cost, at the expense of accuracy.

By choosing rough tolerances considerable savings in time may be realized for systems with few atoms. But for larger systems, e.g. with 8 atoms or more, the time is limited by the plane wave part of the mixed basis, and there is little advantage to playing with this tolerance. The values autogenerated by **lmfgwd –job=-1** are usually pretty safe.

^{6 } This line is used as a marker to identify the next section (keyword is **atom**), so it should be present.

^{7-18} This table specifies the number of partial waves in the valence **nnvv**, and number of core levels **nnc**, one line for every site and *l*. **nnvv** will customarily be two for partial waves and , or three for channels that also have a local orbital . For the first (Ni) atom, there are three *s* core levels, but only a core because the is included in the valence as a local orbital (note **nnvv=3** for the *p*). O (site 3 or 4) has only one core level, the .

This table is made automatically by **lmfgwd –job=-1**, and it is tied to the basis. If you change the basis, be sure to check that the table is still valid (**lmfgwd –job=-2**) or regenerate the table.

^{19} This line is used as a marker to identify the next section (product basis for valence states, keyword **atom**), so it should be present.

^{12-18} This table specifies which partial waves should be used to construct the product basis. There is one line for every *l*, valence partial wave, and every site. The lines in this table must synchronize with the preceding one. Recall that two partial waves combine to make a single product basis function. From the table in this section, two lists of functions are assembled: one for the first function which enters into the bra (|…>) part of the matrix element and another for the second which enters into the ket (<…|) part of the matrix element.

Column

**l**refers to*l*quantum number of the partial waveColumn

**n**refers to the type of partial wave (1 for , 2 for , 3 for ).Column

**occ**if nonzero, is added to the “bra” listColumn

**unocc**if nonzero, is added to the “ket” list.

The remainder of the line is not used, but the template maker supplies a label so you can readily identify the kind of orbital it is. **_p** refers to principal wave (), **_d** refers to energy derivative () **_l** refers to local orbital ().

We have found that including the adds little to the effectiveness of product basis, so the template maker leaves them out. You can always include them by putting **1** in the **occ** and/or **unocc** columns.

In the sample section above, the default was altered. For Ni partial waves entering into the “bra” set is : , and the “ket” list is .

This basis is likely larger than needed, e.g. the inclusion of orbitals, and in the **unocc** column are probably not necessary.

^{44} This line is used as a marker to identify the next section (product basis for core states, keyword **atom**), so it should be present.

^{45-51} This section is the analog for core states of the preceding table. Columns **l**, **n**, **occ**, **unocc** serve the same function as for the preceding section. Column **ForX0** tells the *GW* code to include this state in the calculation of the polarizability. There is an orthogonalization issue, however, see line **4** in the general section for a discussion.

Column **ForSxc** tells the *GW* code to include this state in the calculation of the self-energy.

Return to Table of Contents

### The QPNT category

This category is designed for one-shot *GW*, to select which *k* points and QP levels for which to make the self-energy.

*Note*: if your system is metallic there are restrictions: you must calculate the self-energy for all *k* points and occupied states, in order to obtain the Fermi level.

A typical instance of this category is shown below. (The first column was added to label lines; they are described following the box.)

```
<QPNT> ! Specify particular k-points and bands (one-shot mode)
--- Specify qp and band indices at which to evaluate Sigma
*** Sigma at all q -->1; to specify q -->0. Second arg : up only -->1, otherwise 0
0 0
*** no. states and list of band indices to make Sigma and QP energies
4
4 5 6 7
*** q-points (must belong to mesh of points in BZ).
3
1 0.0000000000000000 0.0000000000000000 0.0000000000000000
3 -0.5000000000000000 0.5000000000000000 0.5000000000000000
7 0.0000000000000000 0.0000000000000000 1.0000000000000000
...
</QPNT>
```

The category as written above applies to 1-shot *GW* or QS*GW* calculations. (Extra information must be given when generating the dynamical self-energy; this will be described shortly.)

Lines preceding ‘*******’ or ‘**###**’ are not used (only line **2** in this case). Lines beginning with one of these are markers that flag the start of a subsection following the line.

^{3-4}. This is the marker for the first subsection and its contents. The two integers on line 4 should be 0 or 1. If the first is 1, QP levels are calculated for all **q** in the irreducible BZ. If it is zero, QP levels are computed only for a list following the third instance of ******* (lines **8-13** in this case). For a metal, you should always pick **1** so the Fermi level can be found.

The second integer is usually 0. Choosing 1 causes the code to calculate QP levels for the first spin only. (Useful for an antiferromagnetic system such as NiO).

^{5-7}. This is the marker for the second subsection and its contents. Its purpose is to select particular QP levels which calculate. (You can of course calculate them all but it is faster to make only a subset.) Line **6** tells the reader how many levels to calculate; line **7** lists the levels. It is an error for the list in line 7 to be smaller than number of levels you select.

^{8-13}. This is the marker for the third subsection and its contents. Its purpose is to select the **q** points. (If the first integer in line **4** is 1, this subsection is not used.) Line **9** tells the reader how many points **q** you want QP levels for; call this **nq**. Then follow **nq** lines specifying the **q** points.

- The first number is not used
- Each
**k**point must be a point on the mesh you specify; see**n1n2n3**, line**14**of the general category. The template maker generates a list of all the irreducible points (it was 8 points in this case). You can use a text editor to shuffle the rows (as was done here, to select the Γ, L and X points of a cubic material).

#### Extra input for the dynamical self-energy

**hsfp0 –job=4** is a special mode that makes the dynamical self-energy. In this mode the reader parses lines in the **QPNTS** category until one is found begining with ‘********’ (the line must be present). Two numbers must appear in the next line. The two lines should be inserted before the normal part of this category described above.

An example is:

```
***** ---Specify the q and band indices for which we evaluate the omega dependence of self-energy ---
0.01 2 (Ry) ! dwplot omegamaxin(optional) : dwplot is mesh for plotting.
```

**0.01** and **2** are the energy spacing and range over which is calculated. In this mode the usual **dw** and **omg_c** (lines **23** and **24** in the general section). The entire range is twice the second number: it extends both above and below the Fermi level.

### The PBASMAX category

Not documented yet.

### The QforEPS and QforEPSL categories

These categories apply to modes to calculate dielectric or magnetic response functions. Use one or the other.

Response functions are calculated at selected *k* points. If you add a line

```
QforEPSIBZ on
```

then the list becomes all the **q** in the irreducible BZ.

Otherwise, use a **QforEPS** or **QforEPSL** category. In the former, you specify a list of **q** points, e.g.

```
<QforEPS>
0d0 0d0 0.01d0
0d0 0d0 0.02d0
0d0 0d0 0.04d0
</QforEPS>
```

Alternatively specify **q** points along symmetry lines (useful for plotting a response function along particular symmetry lines). For example, for a cubic system, the following yields seven **q** points from X to Γ and another seven from Γ to L (the left end point is omitted).

```
<QforEPSL>
1d0 0d0 0d0 0d0 0d0 0d0 8
0d0 0d0 0d0 .5d0 .5d0 .5d0 8
</QforEPSL>
```