X-Ray Absorption Spectra Using BSE

for helium version of Exciting


Purpose: In this tutorial we consider the X-ray absorption near-edge structure (XANES) spectra of core-states by solving the Bethe-Salpeter equation (BSE). The general way to set up the BSE computations for XANES in exciting is shown, together with examples for the Li-K edge in lithium fluoride and N-K edge in cubic boron nitride. The 1s state corresponds to the K edge in spectroscopic notation.



Introduction

In order to solve the BSE, we first need eigenenergies and wave functions from a ground-state calculation as input. In the present methodology the core-edge of interest must be in the valence or otherwise set as a local orbital, and not core state.

If the corresponding core-edge state is not in the valence and not already set as a local orbital, it is possible to do so by hand. This procedure will be shown explicitly in the example for the N-K edge in cubic BN further down below.

A description of basic setup and parameters can be found in Excited States From TDDFT. Further, all the input parameters are described in Input Reference.

Observe that for production results the X-ray absorption spectra need to be converged regarding several parameters such as the k-mesh, see Excited States From TDDFT.

Notice that:

  • The EXCITINGROOT variable has to be replaced in the subsequent input files by its actual value (the path itself, e.g., /home/user/exciting).

The Li-K edge in lithium fluoride

Read through the text below in order to complete 'Exercise 1'

Let us start with the Li-K edge in LiF. The Li 1s state is quite shallow in energy and it is not put into the core by default. This can be directly seen by checking the Li atom species file EXCITINGROOT/species/Li.xml

...
    <atomicState n="1" l="0" kappa="1" occ="2.00000" core="false"/>
...
    <lorb l="0">
      <wf matchingOrder="0" trialEnergy="0.1500" searchE="false"/>
      <wf matchingOrder="1" trialEnergy="0.1500" searchE="false"/>
      <wf matchingOrder="0" trialEnergy="-1.8034" searchE="true"/>
    </lorb>
...

Note that core="false". This means that we don't need to include any new local orbital, and we can directly use the Li species file EXCITINGROOT/species/Li.xml.

Below is an example input.xml file which can be used for the LiF <groundstate> and the Li-K edge calculations <xs> in the exercise (do not forget to change speciespath)

<input>
 
  <title>Lithium fluoride</title>
 
  <structure speciespath="EXCITINGROOT/species">
    <crystal>
      <basevect>3.80402 3.80402 0.00000</basevect>
      <basevect>3.80402 0.00000 3.80402</basevect>
      <basevect>0.00000 3.80402 3.80402</basevect>
    </crystal>
    <species speciesfile="Li.xml">
      <atom coord="0.0000  0.0000  0.0000" bfcmt="0.0 0.0 0.0" />
    </species>
    <species speciesfile="F.xml">
      <atom coord="0.5000  0.5000  0.5000" bfcmt="0.0 0.0 0.0" />
    </species>
  </structure>
 
  <groundstate lradstep="2" lmaxvr="8" ngridk="4  4  4" vkloff="0.05 0.15 0.25" />
 
  <xs xstype="BSE" nosym="true" 
   ngridq="2 2 2" ngridk="2 2 2"
   vkloff="0.05 0.15 0.25" nempty="20" lmaxapwwf="3"
   lmaxemat="3" gqmax="2.0" broad="0.0073499"
   tevout="true" scissor="0.00000">
    <energywindow intv="0.0 5.0" points="2000" />
    <screening screentype="full" nempty="115" />
    <BSE bsetype="singlet" nstlbsemat="1 1 1 15" nstlbse="1 1 1 15" />
    <qpointset>
      <qpoint>0.0 0.0 0.0</qpoint>
    </qpointset>
  </xs>
 
</input>

nstlbsemat and nstlbse in <BSE> give the band ranges for the computation of the matrix elements and for diagonalizing the Hamiltionian matrix, respectively.

Here the band ranges are set to "1 1 1 15", where the first part, 1 1, refers to the core-edge state as shown, e.g., in the ground-state output EIGVAL.OUT and the second part, 1 15, gives the conduction states counted from the Fermi level.

bsetype="singlet" in <BSE> corresponds to including all terms in the effective Hamiltonian (diagonal, direct and exchange).
bsetype="RPA" gives the random phase approximation, which does not include the direct term and, therefore, does not give rise to bound core-exciton states.

Important and useful settings for input.xml

In order to not recalculate the ground-state

<groundstate do="skip" ...

To avoid recalculating the screening and matrix elements, set within <xs>

<plan> 
  <doonly task="bse" />
</plan>

This is useful for instance when computing a different bsetype afterwards, such as bsetype="RPA", or when changing the band ranges for the Hamiltonian diagonalization, e.g., nstlbse="1 1 1 25". Note that the present way of performing this setup might change in future versions of the code.

XANES spectra and core-exciton output from <xs>

INFOXS.OUT shows the general messages for the calculation procedure.

The absorption spectra is given in EPSILON_BSEsinglet_SCRfull_OC11.OUT, with a Lorentzian broadening broad defined in <energywindow> in <xs>.
The energy and oscillation strength of the excited states are listed in EXCITON_BSEsinglet_SCRfull_OC11.OUT, where negative energies correspond to bound core-excitons.

The XANES spectra can, e.g., be plotted using an xml visualization template, generating .agr files readable by xmgrace

$ xsltproc --stringparam filename dielectricfunction EXCITINGROOT/xml/visualizationtemplates/xmldielectricfunction2agr.xsl EPSILON_BSEsinglet_SCRfull_OC11.OUT.xml

The absorption spectra correspond to the imaginary part of the dielectric function, which can be plotted, e.g., by invoking xmgrace

$ xmgrace dielectricfunction_Im.agr &

Exercise 1

  • Groundstate calculation for LiF.
  • BSE calculations for the Li-K edge bsetype="singlet".
  • Plot the absorption spectra, e.g. using EPSILON_BSEsinglet_SCRfull_OC11.OUT.xml.
  • Check after bound core-exciton states in EXCITON_BSEsinglet_SCRfull_OC11.OUT.
  • RPA calculation (bsetype="RPA").
  • Plot RPA absorption spectra, e.g. using EPSILON_BSERPA_SCRfull_OC11.OUT.xml and compare with BSE.
  • Use a denser mesh for the BSE (singlet) calculations ngridq="4 4 4" ngridk="4 4 4".

The figure below shows the Li-K edge spectra in LiF using 15 conduction states for different choice of mesh and gqmax
LiF.png

The N-K edge in cubic boron nitride

'Exercise 2' is given at the end of the section

Here we consider the N-K edge in c-BN. Below is an example which can be used for the regular c-BN ground-state calculation (do not forget to change speciespath)

<input>
 
    <title>cubic boron nitride </title>
 
        <structure autormt="true" speciespath="EXCITINGROOT/species">
                <crystal scale="6.83136">
                        <basevect>0.5 0.0 0.5</basevect>
                        <basevect>0.0 0.5 0.5</basevect>
                        <basevect>0.5 0.5 0.0</basevect>
                </crystal>
                <species speciesfile="B_APWloLegacy.xml" >
                        <atom coord="0.0 0.0 0.0" bfcmt="0.0 0.0 0.0" />
                </species>
                <species speciesfile="N_APWloLegacy.xml">
                        <atom coord="0.25 0.25 0.25" bfcmt="0.0 0.0 0.0" />
                </species>
        </structure>
 
        <groundstate gmaxvr="14.0" ngridk="4 4 4"
            swidth="0.01d0" findlinentype="simple" mixer="lin" />
 
</input>

The N 1s state is relatively deep, around 1 keV, and it is by default a core-state in the species file EXCITINGROOT/species/N.xml

...
    <atomicState n="1" l="0" kappa="1" occ="2.00000" core="true"/>
...

In order to perform the BSE calculations, we must include N 1s as a local orbital.

Including a local orbital

First, we need to estimate the energy location of the state in order to make a good guess for the local orbital trial energy.

This is done by a regular ground-state calculation using the species file EXCITINGROOT/species/N_APWloLegacy.xml.

After successful convergence, check the energy of the state in EVALCORE.OUT, it may look like

Species :    1 (B), atom :    1
 n =  1, l =  0, k =  1 :   -5.918319395    

Species :    2 (N), atom :    1
 n =  1, l =  0, k =  1 :   -13.31776506

Now we will make a separate ground-state calculation, including the local orbital.

Copy the input.xml file and the EXCITINGROOT/species/N_APWloLegacy.xml to a new working directory. It's useful to change the name of the species file to, e.g., N_APWloLegacy_XANES.xml to not confuse it with others.

In the new species file N_APWloLegacy_XANES.xml, remove the core-edge state from the core by setting core="true" to "false"

...
<atomicState n="1" l="0" kappa="1" occ="2.00000" core="false"/>
...

Insert a local orbital by adding

...
<!-- Here below a local orbital is added for 1s -->
  <lorb l="0">
    <wf matchingOrder="0" trialEnergy="0.1500" searchE="false"/>
    <wf matchingOrder="1" trialEnergy="0.1500" searchE="false"/>
    <wf matchingOrder="0" trialEnergy="-13.5000" searchE="true"/>
  </lorb>
...

The most important is to set a "good enough" trialEnergy in the last row, it is usually a bit lower than the result for the core-state. Here -13.5 Hartree is selected (compare with -13.31776506 Hartree from above).

Do not forget to edit speciespath in the new input.xml.

If calculations with a new local orbital does not converge, it might be necessary to set a different trialEnergy. Also, it might be necessary to explicitly set the parameters deband (increase default value) and epsband (decrease, if needed). See Input Reference for details.

After a successful ground-state calculation including the new local orbital, it is time for considering the XANES computation. An example is given below for bsetype="singlet"

    <xs xstype="BSE" nosym="true" ngridk="2 2 2" ngridq="2 2 2"
        vkloff="0.05 0.15 0.25"
        nempty="60" lmaxapwwf="5" lmaxemat="5" gqmax="2.0" broad="0.0183747"
        tevout="true">
        <energywindow intv="13.4 16.5" points="1550" />
               <screening screentype="full" nempty="115" />
               <BSE bsetype="singlet" nstlbsemat="1 1 1 15" nstlbse="1 1 1 15" />
        <qpointset>
            <qpoint>0.0 0.0 0.0</qpoint>
        </qpointset>
    </xs>

Exercise 2

  • A regular ground-state calculation for c-BN.
  • Include the N 1s local orbital in a new ground-state calculation, compare EIGVAL.OUT.
  • BSE calculation for the N-K edge bsetype="singlet".
  • Plot the absorption spectra, e.g. using EPSILON_BSEsinglet_SCRfull_OC11.OUT.xml.
  • Check after bound core-exciton states in EXCITON_BSEsinglet_SCRfull_OC11.OUT.
  • RPA calculation (bsetype="RPA").
  • Plot RPA absorption spectra, e.g. using EPSILON_BSERPA_SCRfull_OC11.OUT.xml and compare with BSE.
  • Use a denser mesh for the BSE (singlet) calculations ngridq="4 4 4" ngridk="4 4 4".

The figure below shows the B-K edge in c-BN using 15 conduction states and for different choice of mesh and gqmax
cBN.png

Exercise 3

  • Check the sensitivity of the LiF and c-BN spectra to different parameters, e.g., conduction states (nstlbsemat and nstlbse).
  • Calculate the boron-K edge in c-BN.

Some tips and hints for production calculations

  • Copy atom species files ending with _APWloLegacy.xml for including new local orbitals.
  • If convergence problems, check trialEnergy in species file, deband and epsband in input.xml.
  • Investigate the sensitivity of the spectra on the number of conduction states.
  • For studying the very near-edge spectra, fewer conduction states are typically needed.
  • For comparison with a low resolution experiment, a less dense mesh in <xs> could be sufficient.
  • In general, BSE computations may be impractical for unit cells containing many atoms.
  • Core-exciton binding energies are more difficult to converge than the general spectral shape.
  • L2,3-edge (2p) spectra are more costly than K-edge (more core states included, larger matrix to diagonalize).
Unless otherwise stated, the content of this page is licensed under Creative Commons Attribution-ShareAlike 3.0 License