A Portal Site of Materials Science Simulation

There are many possible ways to partition the total electron density into contributions from each atom. One popular method is Bader partitioning, where the electron density is divided at zero flux surfaces. A merit in the Bader partitioning method is that in comparison to Mulliken analysis, which is another popular method, the results depend little on the employed basis set. Here, we will walk through the process of performing Bader charge analysis on NaCl crystal using Quantum Espresso. The following procedure has been tested on MateriApps Live! 2.4.

Installation of bader

Quantum Espresso comes preinstalled in MateriApps Live!. However, we will need to download the bader program from Henkelman Group’s web site:

$ wget  http://theory.cm.utexas.edu/henkelman/code/bader/download/bader_lnx_64.tar.gz
$ tar xzvf bader_lnx_64.tar.gz
$ mv bader ~/bin

Self consistent calculation of electron density using Quantum Espresso

We will analyze results obtained using Quantum Espresso. Let’s create an input file for self-consistent calculation and name it “NaCl.scf.in”.

&control
calculation = 'scf'
prefix = 'NaCl'
pseudo_dir = './'
wf_collect = .true.
/
&system
ibrav = 2
celldm(1) = 10.66
nat = 2
ntyp = 2
ecutwfc = 60
ecutrho = 244
/
&electrons
mixing_mode = 'plain'
mixing_beta = 0.7
conv_thr = 1.0d-8
/
ATOMIC_SPECIES
Na 1 Na.pbesol-spn-kjpaw_psl.1.0.0.UPF
Cl 1 Cl.pbesol-n-kjpaw_psl.1.0.0.UPF
ATOMIC_POSITIONS
Na 0.00 0.00 0.00
Cl 0.5 0.5 0.5
K_POINTS {automatic}
8 8 8 0 0 0

We use the following PBEsol pseudopotentials from the Quantum Espresso pslibrary:

• Cl: https://www.quantum-espresso.org/upf_files/Cl.pbesol-n-kjpaw_psl.1.0.0.UPF
• Na: https://www.quantum-espresso.org/upf_files/Na.pbesol-spn-kjpaw_psl.1.0.0.UPF

They can be downloaded to the current directory as

$ wget https://www.quantum-espresso.org/upf_files/Cl.pbesol-n-kjpaw_psl.1.0.0.UPF
$ wget https://www.quantum-espresso.org/upf_files/Na.pbesol-spn-kjpaw_psl.1.0.0.UPF

Now run Quantum Espresso as

$ pw.x < NaCl.scf.in > NaCl.scf.out

Post-processing to obtain cube files containing density data on a 3D mesh

Once that’s successfully done, we need to convert Quantum Espresso output to electron density files that are readable by bader. This can be achieved by pp.x (post processing) utility provided with the Quantum Espresso package. To perform the Bader partitioning, we need to prepare both the all-electron density and the valence electron density (an explanation can be found here or here). Here are the input files for the valence density (let’s name this NaCl_val.pp.in):

&inputpp
    prefix  = 'NaCl'
    outdir = './'
    filplot = 'NaCl_valence'
    plot_num= 0
/
&plot
nfile = 1
iflag = 3
output_format = 6
fileout = 'NaCl_valence.cube'
/

and all-electron density (NaCl_all.pp.in):

&inputpp
    prefix  = 'NaCl'
    outdir = './'
    filplot = 'NaCl_allelec'
    plot_num= 21
/
&plot
nfile = 1
iflag = 3
output_format = 6
fileout = 'NaCl_allelec.cube'
/

Make sure you have an empty newline at the end of each file. Otherwise pp.x will be unable to parse the input and will terminate with error. The most important difference between the two input files in our case is the plot_num line which specifies the physical quantities (all-electron or valence densities) that are written to the fileout file (the manual describing all options can be found here). The pp.x can be run as follows:

$ pp.x < NaCl_val.pp.in > NaCl_val.pp.out
$ pp.x < NaCl_all.pp.in > NaCl_all.pp.out

Bader charge analysis

Now we can perform Bader partitioning on the .cube files created by pp.x:

$ bader NaCl_valence.cube -ref NaCl_allelec.cube

The Bader charges on each atom is contained in the ACF.dat file which should look similar to the following:

    #         X           Y           Z       CHARGE      MIN DIST   ATOMIC VOL
 --------------------------------------------------------------------------------
    1  -10.660000   10.660000   10.660000    8.140837     1.940149    64.992061
    2   -5.330000    5.330000    5.330000    7.859147     2.931500   237.846813
 --------------------------------------------------------------------------------
    VACUUM CHARGE:               0.0000
    VACUUM VOLUME:               0.0000
    NUMBER OF ELECTRONS:        16.0000

This shows that out of 16 valence electrons, 8.14 is assigned to Na and 7.86 to Cl. To convert this to ionic charge, we will need to check the number of electrons that are fixed in the core. This can be achieved by looking at these lines in the pseudopotential files for Na

    Valence configuration:
    nl pn  l   occ       Rcut    Rcut US       E pseu
    2S  1  0  2.00      1.000      1.250    -4.158089
    3S  2  0  1.00      1.000      1.250    -0.198813
    2P  2  1  6.00      0.900      1.300    -2.106125

and Cl

    Valence configuration:
    nl pn  l   occ       Rcut    Rcut US       E pseu
    3S  1  0  2.00      1.200      1.600    -1.515413
    3P  2  1  5.00      1.300      1.600    -0.629258

The Na pseudopotential indicates that two 1s electrons are fixed to the core, while the Cl pseudopotential indicates that two 1s, two 2s, and six 2p electrons are fixed to the core. From this, we can calculate the total electron number on Na as (2 + 8.14) = 10.14 and on Cl as (2 + 2 + 6 + 7.86) = 17.86. By subtracting these values from the atomic numbers for Na (11) and Cl (17), we arrive at Bader charge-based ionic charge states of Na^0.86+ and Cl^0.86-. This is slightly smaller than the nominal ionic valences of ± 1, and indicates that there is some covalent character even in a system like NaCl which is considered a prototypical ionic crystal.