Van-der-Waals Energy Functional

for helium version of Exciting

Purpose: In this tutorial we show how to set up and execute van-der-Waals density-functional (vdW-DF) calculations using the exciting code. In particular, we consider the case of graphite, where the vdW forces are known to play a role. Structurally, graphite can be viewed as a stack of graphene (a plane of hexagonally packed carbon atoms) layers which are bound by weak vdW forces. Our aim is to find the equilibrium interlayer separation and the corresponding binding energy.

0. Define relevant shell variables and download scripts

Very important: Before starting, the following shell variable must be set by the user:

EXCITINGROOT = exciting installation directory, e.g.: /home/user/exciting .

To do it in a bash shell, one needs to type

$ export EXCITINGROOT=/local_path_to_exciting_root

The directory name /local_path_to_exciting_root is a dummy name which must be explicitly changed by the user to the appropriate ones.

As a second step, it will be necessary to download input files and scripts to start calculations. To execute the scripts one can use the command

$ bash <the corresponding script>

$ chmod u+x <the corresponding script>
$ ./<the corresponding script>

Requirements: Bash shell and awk.

1. Graphite

In this example, we examine the equilibrium interlayer distance in graphite and determine the corresponding binding energy.

i) Preparation of the input file

The first step is to create a directory for each system that you want to investigate. Therefore, one should create a directory, e.g., graphite. Inside this directory, we create (e.g., just copying from below or download from here) the file input_graphite.xml corresponding to a template of an exciting input file. This file contains all important parameters to perform self-consistent calculations. The only parameters which are not specified (denoted by 'XXXXX' and 'YYYYY') will be assigned during the script execution.

<?xml version="1.0" encoding="UTF-8"?>
<input xsi:noNamespaceSchemaLocation="$EXCITINGROOT/xml/excitinginput.xsd"
  xsltpath="$EXCITINGROOT/xml" scratchpath="/tmp">
  <structure speciespath="$EXCITINGROOT/species" autormt="true">
    <crystal scale="1.88973" stretch="2.461 2.461 XXXXX">
      <basevect>0.5000000000 0.8660254038 0.0000000000</basevect>
      <basevect>1.0000000000 0.0000000000 0.0000000000</basevect>
      <basevect>0.0000000000 0.0000000000 1.0000000000</basevect>
    <species speciesfile="C.xml">
      <atom coord="0.00000000 0.00000000 0.00000000"></atom>
      <atom coord="0.00000000 0.00000000 0.50000000"></atom>
      <atom coord="0.66666667 0.66666667 0.00000000"></atom>
      <atom coord="0.33333333 0.33333333 0.50000000"></atom>
    ngridk="4  4  2" 
    vkloff="0.5  0.5  0.5"
              <box grid="20 20 YYYYY">
            <origin coord="0 0 0" />
                    <point coord="1 0 0"/>
                    <point coord="0 1 0"/>
                    <point coord="0 0 1"/>

Please, note also that the input file for the exciting calculation (which is always called input.xml) will be created by the script Furthermore, for some parameters (rgkmax, ngridk, and box grid) rather low values are used as it would be required by convergence tests and used here only to speed up the calculations.

ii) Execute the calculations

The total energy and charge-density distribution

In this tutorial, the vdW-DF formalism is applied as a post-scf procedure (following the original paper by Dion and co-workers (PRL, 2004)). In order to do this, we perform first a self-consistent run using the PBE-GGA xc functional and generate the charge density, which is then used to evaluate the vdW-DF total energy given by expression:

\begin{align} E_{vdW-DF} = E^{PBE-GGA} - E_{xc}^{PBE-GGA} + \left( E_x^{revPBE-GGA} + E_c^{LDA}\right) + E_{c}^{NL}. \end{align}

To run calculations, you need to download the script (download it here).

This script performs a series of calculations for a set of graphite interlayer distances (specified inside the script) in order to obtain the first four energy contributions appearing in the r.h.s. of Eq. (1). It also generates a real-space 3D charge-density distribution in a format suitable for visualization with the XCrySDen program and stores results in a separated directory called output. The script should be executed in your work directory:

$ ./
Calculation of the non-local correlation energy correction

The charge-density distribution is the only input for the vdW-DF program, which evaluate the last term, EcNL, of Eq. (1). The non-local correlation correction EcNL is given by the six-dimensional integral:

\begin{align} E_c^{NL} = \frac{1}{2} \int \int n(\textbf{r}) \phi(\textbf{r},\textbf{r}') n(\textbf{r}') d\textbf{r} d\textbf{r}', \end{align}

where the integration is over the crystal volume. The routine vdWDF.x (to be found in $EXCITINGROOT/src/vdwdf/) performs the integration by using the Monte-Carlo (MC) method implemented in the external Cuba library (see here). For practical application, the user must create an input file, as shown below (this file is included in the script, to be downloaded here):

# Integration parameters ( FORMATTED INPUT (line order is important) )
  vdW-DF_Phi_read.dat    ! the tabulated vdW-DF kernel (\phi(r,r'))
  vdW-DF                 ! vdW-DF type (vdW-DF, vdW-DF2, VV09)
  8 8 8                  ! number of unit cells in a/b/c-direction
  1.d-16                 ! relative integration accuracy
  1.d-03                 ! absolute integration accuracy

# Density data file           

# Additional integration parameters (see Cuba-Library manual)
  0                      ! flags
  96                     ! mineval: minimal number of int. evaluations
  1000000000             ! maxeval: maximum number of int. evaluations
  100000                 ! key1: MC sampling number
  1                      ! key2: integration in the final int. phase
  1                      ! key3: refinement phase
  5                      ! maxpass
  0.00d0                 ! border
  1.00d0                 ! maxchisq
  0.25d0                 ! mindeviation

Important parameters, which must be specified, are

  1. the name and path to the tabulated vdW-DF kernel (stored in the source directory);
  2. the version of vdW-DF (experimental option);
  3. three integer numbers which determine the integration super cell (obtained by the unit cell translations along (+/-) a, b, and c directions). Another interpretation of these numbers is a cutoff radius of vdW interactions;
  4. the relative and (or) absolute accuracy of the MC integration estimator (the specifics of MC integration). Note that the super cell size and the integration accuracy are the quantities are system-dependent parameters and must be chosen carefully based on convergence tests.

To obtain EcNL in this example, one just needs to move to the results directory output and download the script Edit this script, choose 8 8 8 for the super cell size and 1d-3 for the absolute integration accuracy, and execute it:

$ ./
Visualization of the results

As a final step, we visualize the results of calculations. The script (download it here) processes all output files and prepare data for plotting:

$ ./ > plot.dat

Now one can easily plot the data with the Xmgrace program

$ xmgrace -nxy plot.dat &

To obtain nice plots one can download here the xmgr parameters file and load it in xmgrace (Plot / Load parameters …). The resulting plot should look like


One should note that used in this example integration accuracy $\bf{\epsilon=1e-3}$ is not sufficient to obtain the smooth EcNL curve and used here only for the execution time reason.

2. Additional exercises

  • Perform a convergence tests to find the optimal parameters (rgkmax, ngridk, vdW-DF integration volume, and integration accuracy).
  • Test the performance of the different exchange-correlation potentials (LDA, GGAs).
  • Compare the obtained results with the known literature: Adsorption of benzene and naphthalene on graphite, PRL 2006.
Unless otherwise stated, the content of this page is licensed under Creative Commons Attribution-ShareAlike 3.0 License