Atomic setups for PAW calculations

1 Introduction
2 What defines a setup?
3 Specification of the elements
- 3.1 The header
- 3.2 A comment
- 3.3 The atom element
- 3.4 Exchange-correlation
- 3.5 Generator
- 3.6 Energies
- 3.7 Valence states
- 3.8 Radial grids
- 3.9 Shape function for the compensation charge
- 3.10 Radial functions
- 3.11 Kinetic energy differences
4 The Kresse-Joubert formulation
5 How to use the setups
6 Plotting the radial functions
- 6.1 Usage
7 References

1 Introduction

This page contains information about the PAW-XML data format for the atomic setups necessary for doing projector-augmented wave calculations[4]. We use the term setup because the PAW method is not a pseudopotential method.

An example XML file for nitrogen PAW setup using LDA can be seen here: N.LDA.

This document as a pdf file: pawxml.pdf

Note

Hartree atomic units are used in the XML file ( $ℏ = m = e = 1$ ).

2 What defines a setup?

The following quantities defines a minimum PAW setup (the notation from Ref. [1] is used here):

Quantity Description

$Z$ atomic number

$E_{XC} [n]$ exchange-correlation functional

$E_{c}^{kin}$ kinetic energy of the core electrons

$g_{ℓ m} (r)$ shape function for compensation charge

$n_{c} (r)$ all-electron core density

${\tilde{n}}_{c} (r)$ pseudo electron core density

${\tilde{n}}_{v} (r)$ pseudo electron valence density

$\overline{v} (r)$ zero potential

$ϕ_{i} (r)$ all-electron partial waves

${\tilde{ϕ}}_{i} (r)$ pseudo partial waves

${\tilde{p}}_{i} (r)$ projector functions

$Δ E_{i j}^{kin}$ kinetic energy differences

Contents of a PAW setup (in order of appearance).

Quantity	Description
$Z$	atomic number
$E_{XC} [n]$	exchange-correlation functional
$E_{c}^{kin}$	kinetic energy of the core electrons
$g_{ℓ m} (r)$	shape function for compensation charge
$n_{c} (r)$	all-electron core density
${\tilde{n}}_{c} (r)$	pseudo electron core density
${\tilde{n}}_{v} (r)$	pseudo electron valence density
$\overline{v} (r)$	zero potential
$ϕ_{i} (r)$	all-electron partial waves
${\tilde{ϕ}}_{i} (r)$	pseudo partial waves
${\tilde{p}}_{i} (r)$	projector functions
$Δ E_{i j}^{kin}$	kinetic energy differences

3 Specification of the elements

An element looks like this:

<name> ... </name>

or for an empty element:

<name/>

Tip

An XML-tutorial can be found here

3.1 The header

The first two lines should look like this:

<?xml version="1.0"?>
<paw_setup version="0.6">

The first line must be present in all XML files. Everything else is put inside an element with name paw_setup, and this element has an attribute called version. We are currently at version 0.6.

3.2 A comment

It is recommended to put a comment giving the units and a link to this web page:

<!-- Nitrogen setup for the Projector Augmented Wave method. -->
<!-- Units: Hartree and Bohr radii.                          -->
<!-- http://www.where.org/paw_setup.html                     -->

3.3 The `atom` element

<atom symbol="N" Z="7" core="2" valence="5"/>

The atom element has attributes symbol, Z, core and valence (chemical symbol, atomic number, number of core electrons and number of valence electrons).

3.4 Exchange-correlation

<xc_functional type="LDA", name="PW"/>

The xc_functional element defines the exchange-correlation functional used for generating the setup, and we currently have the following possibilities for LDA type functionals:

name reference

W Wigner[12]

HL Hedin, Lundqvist[11]

GL Gunnarson, Lundqvist[15]

VWN Vosko, Wilk, Nusair[14]

PZ Perdew, Zunger[10]

PW Perdew, Wang[6]

`name`	reference
`W`	Wigner[12]
`HL`	Hedin, Lundqvist[11]
`GL`	Gunnarson, Lundqvist[15]
`VWN`	Vosko, Wilk, Nusair[14]
`PZ`	Perdew, Zunger[10]
`PW`	Perdew, Wang[6]

For type="GGA" we have these possibilities:

name reference

BLYP Becke[16] + Lee, Yang, Parr[17]

PW91 Perdew, Wang[13]

PBE Perdew, Burke, Ernzerhof[7]

revPBE Zhang, Yang[8]

RPBE Hammer, Hansen, Nørskov[9]

`name`	reference
`BLYP`	Becke[16] + Lee, Yang, Parr[17]
`PW91`	Perdew, Wang[13]
`PBE`	Perdew, Burke, Ernzerhof[7]
`revPBE`	Zhang, Yang[8]
`RPBE`	Hammer, Hansen, Nørskov[9]

3.5 Generator

<generator type="scalar-relativistic" name="MyGenerator-2.0">
  Frozen core: [He]
</generator>

This element contains character data describing in words how the setup was generated. The type attribute must be one of: non-relativistic, scalar-relativistic or relativistic.

3.6 Energies

<ae_energy kinetic="53.777460" xc="-6.127751"
           electrostatic="-101.690410" total="-54.040701"/>
<core_energy kinetic="43.529213"/>

The kinetic energy of the core electrons, $E_{c}^{kin}$ , is used in the PAW method. The other energies are convenient to have for testing purposes and can also be useful for checking the quality of the underlying atomic calculation.

3.7 Valence states

<valence_states>
  <state n="2" l="0" f="2"  rc="1.10" e="-0.6766" id="N-2s"/>
  <state n="2" l="1" f="3"  rc="1.10" e="-0.2660" id="N-2p"/>
  <state       l="0"        rc="1.10" e=" 0.3234" id="N-s1"/>
  <state       l="1"        rc="1.10" e=" 0.7340" id="N-p1"/>
  <state       l="2"        rc="1.10" e=" 0.0000" id="N-d1"/>
</valence_states>

The valence_states element contains several state elements. For this setup, the first two lines describe bound eigenstates with occupation numbers and principal quantum numbers. Notice, that the three additional unbound states should have no f and n attributes. In this way, we know that only the first two bound states (with f and n attributes) should be used for constructing an initial guess for the wave functions.

3.8 Radial grids

There can be one or more definitions of radial grids:

<radial_grid eq="r=a*i/(n-i)" a="0.400000" n="300" istart="0" iend="299" id="g1"/>

This defines one radial grid as:

\begin{matrix} r_{i} = \frac{a i}{n - i}, \end{matrix}

where $i$ runs from 0 to 299. All functions (densities, potentials, ...) that use this grid are given as 300 numbers defining the radial part of the function. The radial part of the function must be multiplied by a spherical harmonics: $f_{ℓ m} (r) = f_{ℓ} (r) Y_{ℓ m} (θ, ϕ)$ .

Each radial grid has a unique id:

<radial_grid eq="r=a*exp(d*i)" a="1.056e-4" d="0.05"
             istart="0" iend="249" id="log"/>
<radial_grid eq="r=d*i" d="0.01" istart="0" iend="99" id="lin"/>

and each numerical function must refer to one of these ids:

<function grid="lin">
  ... ... ...
</function>

In this example, the function element should contain 100 numbers ( $i = 0, . . ., 99$ ). Each number must be separated by a <newline> character or by one or more <tab>'s or <space>'s (no commas). For numbers with scientific nutation, use this format: 1.23456e-5 or 1.23456E-5 and not 1.23456D-5.

There are currently five types of radial grids:

eq parameters

r=a*exp(d*i)) a and d

r=a*i/(n-i) a and n

r=a*(exp(d*i)-1) a and d

r=d*i d

r=(i/n+a)^5/a-a^4 a and n

Types of radial grids. Note that the eq attribute is not meant to be parsed - only recognized.

`eq`	parameters
`r=aexp(di))`	`a` and `d`
`r=a*i/(n-i)`	`a` and `n`
`r=a(exp(di)-1)`	`a` and `d`
`r=d*i`	`d`
`r=(i/n+a)^5/a-a^4`	`a` and `n`

The istart and iend attributes indicating the range of $i$ should always be present.

3.9 Shape function for the compensation charge

The compensation charge for an atom is expanded using the multipole moments $Q_{ℓ m}$ :

\begin{matrix} \sum_{ℓ m} Q_{ℓ m} {\tilde{g}}_{ℓ} (r) Y_{ℓ m} (θ, ϕ), \end{matrix}

where $g_{ℓ} (r) \propto r^{ℓ} k (r)$ , and $k (r)$ is the shape function defined in the file like this:

<shape_function type="gauss" rc="3.478505426185e-01"/>

Here, Gaussians are used, $exp (- (r / r_{c})^{2})$ , as the basis for the compensation charges as described in Ref. [4].

Another choice would be type="bessel" rc="...", where a sum of two Bessel functions is used[2]. Other possibilities are

\begin{matrix} {(\frac{sin (π r / r_{c})}{π r / r_{c}})}^{2} \end{matrix}

as described in Ref. [3], and $exp (- (r / r_{c})^{λ})$ .

type parameters

gauss rc

bessel rc

sinc rc

exp rc and lamb

`type`	parameters
`gauss`	`rc`
`bessel`	`rc`
`sinc`	`rc`
`exp`	`rc` and `lamb`

There is also the possibility that the shape function is given in numerical form. In that case, the shape function defined inside the shape_function element:

<shape_function grid="...">
  ... ...
</shape_function>

The grid attribute refers to one of the radial grids defined earlier.

3.10 Radial functions

Continuing, we have now reached the all-electron core density:

<ae_core_density grid="g1">
   6.801207147443e+02 6.801207147443e+02 6.665042896724e+02
   ... ...
</ae_core_density>
<pseudo_core_density grid="g1">
   ...
</pseudo_core_density>
<pseudo_valence_density grid="g1">
   ...
</pseudo_valence_density>
<zero_potential grid="g1">
   ...
</zero_potential>

The numbers inside the ae_core_density element defines the radial part of $n_{c} (r)$ . The radial part must be multiplied by $Y_{00} = (4 π)^{- 1 / 2}$ to get the full all-electron core density (which should integrate to the number of core electrons). The pseudo core density, the pseudo valence density and the zero potential, $\overline{v}$ , are defined similarly.

The ae_partial_wave, pseudo_partial_wave and projector_function elements contain the radial parts of the $ϕ_{i} (r)$ , ${\tilde{ϕ}}_{i} (r)$ and ${\tilde{p}}_{i} (r)$ functions for the states listed in the valence_states element above (five states in the nitrogen example). All functions must have an attribute state="..." referring to one of the states listed in the valence_states element:

<ae_partial_wave state="N-2s" grid="g1">
  -8.178800366898029e+00 -8.178246914143839e+00 -8.177654917302689e+00
  ... ...
</ae_partial_wave>
<pseudo_partial_wave state="N-2s" grid="g1">
  ...
</pseudo_partial_wave>
<projector_function state="N-2s" grid="g1">
  ...
</projector_function>
<ae_partial_wave state="N-2p" grid="g1">
  ...
</ae_partial_wave>
...
...

3.11 Kinetic energy differences

  <kinetic_energy_differences>
     1.744042161013e+00 0.000000000000e+00 2.730637956456e+00
     ...
  <kinetic_energy_differences>
</paw_setup>

This element contains the symmetric $Δ E_{i j}^{kin}$ matrix:

\begin{matrix} Δ E_{i j}^{kin} = 〈 ϕ_{i} | \hat{T} | ϕ_{j} 〉 - 〈 {\tilde{ϕ}}_{i} | \hat{T} | {\tilde{ϕ}}_{j} 〉 \end{matrix}

where $\hat{T}$ is the kinetic energy operator used by the generator. With $n$ states, we have an $n \times n$ matrix listed as $n^{2}$ numbers.

4 The Kresse-Joubert formulation

The Kresse-Joubert formulation of the PAW method[2] is very similar to the original formulation of Blöchl[4]. However, the Kresse-Joubert formulation does not use $\overline{v}$ directly, but indirectly through the local ionic pseudopotential, $v_{H} [{\tilde{n}}_{Z c}]$ . Therefore, the following transformation is necessary:

\begin{matrix} v_{H} [{\tilde{n}}_{Z c}] = v_{H} [{\tilde{n}}_{c} + (N_{c} - Z - {\tilde{N}}_{c}) g_{00} Y_{00}] + \overline{v} + v_{x c} [{\tilde{n}}_{v} + {\tilde{n}}_{c}] - v_{x c} [{\tilde{n}}_{v} + {\tilde{n}}_{c} + (N_{v} - {\tilde{N}}_{v} - {\tilde{N}}_{c}) g_{00} Y_{00}] \end{matrix}

where $N_{c}$ is the number of core electrons, $N_{v}$ is the number of valence electrons, ${\tilde{N}}_{c}$ is the number of electrons contained in the pseudo core density and ${\tilde{N}}_{v}$ is the number of electrons contained in the pseudo valence density. The Hartree potential from the density $n$ is defined as:

\begin{matrix} v_{H} [n] (r_{1}) = 4 π \int_{0}^{\infty} r_{2}^{2} d r_{2} \frac{n (r_{2})}{r_{>}}, \end{matrix}

where $r_{>}$ is the larger of $r_{1}$ and $r_{2}$ .

Note

In the Kresse-Joubert formulation, the symbol $\tilde{n}$ is used for what we here call ${\tilde{n}}_{v}$ and in the Blöchl formulation, we have $\tilde{n} = {\tilde{n}}_{c} + {\tilde{n}}_{v}$ .

It is also possible to add an element kresse_joubert_local_ionic_pseudopotential that contains the $v_{H} [{\tilde{n}}_{Z c}] (r)$ function directly, so that no conversion is necessary:

<kresse_joubert_local_ionic_pseudopotential grid="log">
   ...
</kresse_joubert_local_ionic_pseudopotential>

5 How to use the setups

Most likely, the radial functions will be needed on some other type of radial grid than the one used in the setup. The idea is that one should read in the radial functions and then transform them to the radial grids used by the specific implementation. After the transformation, some sort of normalization may be necessary.

6 Plotting the radial functions

The first 10-20 lines of the XML-setups, should be pretty much human readable, and should give an overview of what kind of setup it is and how it was generated. The remaining part of the files contain numerical data for all the radial functions. To get an overview of these functions, you can extract that data with the pawxml.py program and then pass it on to your favorite plotting tool.

Note

The pawxml.py program is very primitive and is only included in order to demonstrates how to parse XML using SAX from a Python program. Parsing XML from Fortran or C code with SAX should be similar.

6.1 Usage

It works like this:

$ pawxml.py [options] setup[.gz]

Options:

`--version`	Show program's version number and exit.
`-h, --help`	Show this help message and exit.
`-x <name>, --extract=<name>`	Function to extract.
`-s<channel>, --state=<channel>`	Select valence state.
`-l, --list`	List valence states

Examples:

[~]$ pawxml.py -x pseudo_core_density N.LDA | xmgrace -
[~]$ pawxml.py -x ae_partial_wave -s N2p N.LDA > N.ae.2p
[~]$ pawxml.py -x pseudo_partial_wave -s N2p N.LDA > N.ps.2p
[~]$ xmgrace N.??.2p

7 References

[1]	P. E. Blöchl, C. J. Forst and J. Schimpl, Projector augmented wave method: Ab initio molecular dynamics with full wave functions, Bulletin of Materials Science 26, 33-41 (2003)

[2]	(1, 2) G. Kresse and D. Joubert, Form ultrasoft pseudopotentials to the projector augmented-wave method, Phys. Rev. B 59, 1758-1775 (1999)

[3]	N. A. W. Holzwarth, A. R. Tackett, and G. E. Matthews, A Projector Augmented Wave (PAW) code for electronics structure calculations: Part I atompaw for generating atom-centered functions, Computer Physics Communications 135, 329-347 (2001)

[4]	(1, 2, 3) P. E. Blöchl, Projector augmented-wave method, Phys. Rev. B 50, 17953-19979 (1994)

[5]	D. Vanderbilt, Soft Self-Consistent Pseudopotentials in a Generalized Eigenvalue Formalism, Phys. Rev. B 41 (Rapid Communications), 7892 (1990)

[6]	J. P. Perdew and Y. Wang, Accurate and simple analytic representation of the electron-gas correlation energy Phys. Rev. B 45, 13244-13249 (1992)

[7]	J. P. Perdew, K. Burke, and M. Ernzerhof, Generalized Gradient Approximation Made Simple, Phys. Rev. Lett. 77, 3865 (1996)

[8]	Y. Zhang and W. Yang, Comment on "Generalized Gradient Approximation Made Simple", Phys. Rev. Lett. 80, 890 (1998)

[9]	B. Hammer, L. B. Hansen and J. K. Nørskov, Improved adsorption energetics within density-functional theory using revised Perdew-Burke-Ernzerhof functionals, Phys. Rev. B 59, 7413 (1999)

[10]	J. P. Perdew and A. Zunger, Self-interaction correction to density-functional approximations for many-electron systems, Phys. Rev. B 23, 5048 (1981)

[11]	L. Hedin and B. I. Lundqvist, Explicit local exchange-correlation potentials, J. Phys. C 4, 2064 (1971)

[12]	E. P. Wigner, Effects of the electron interaction on the energy levels of electrons in metals, Trans. Faraday Soc. 34, 678 (1938)

[13]	J. P. Perdew and Y.Wang, Atoms, molecules, solids, and surfaces: Applications of the generalized gradient approximation for exchange and correlation, Phys. Rev. B 46, 6671 (1992)

[14]	S. H. Vosko, L. Wilk and M. Nusair, Accurate spin-dependent electron liquid correlation energies for local spin density calculations: a critical analysis, Can. J. Phys. 58, 1200 (1980)

[15]	O. Gunnarson and B. I. Lundqvist, Exchange and correlation in atoms, molecules, and solids by the spin-density-functional formalism, Phys. Rev. B 13 ,4274 (1976)

[16]	A. D. Becke, Density-functional exchange-energy approximation with correct asymptotic behavior, Phys. Rev. A 38, 3098 (1988)

[17]	C. Lee, W. Yang and R. G. Parr, Development of the Colle-Salvetti correlation-energy formula into a functional of the electron density, Phys. Rev. B 37, 785 (1988)