q-Dependent BSE Calculations¶
by Christian Vorwerk & Caterina Cocchi for exciting neon
(Jupyter notebook by Megha Arya and Mara Voiculescu)
Purpose: In this tutorial you will learn how to perform Bethe-Salpeter equation (BSE) calculations to obtain x-ray and electron scattering spectra for finite momentum transfer. As an example, the valence and core scattering spectra of LiF are calculated.
0. Before Starting¶
Read the following paragraphs before starting with the rest of this tutorial!
Before running any Jupyter tutorials, please refer to the 00_before_starting.md document on how to correctly set up the environment. This only needs to be done once. After which, the venv can be (re)activated from exciting's root directory:
source tools/excitingjupyter/venv/excitingvenv/bin/activate
1. Theoretical Background¶
$\Rightarrow$ Electron and X-Ray Scattering Spectroscopy
The inelastic, non-resonant scattering of electrons or x-ray photons with the electrons of a system are measured in angle-resolved electron energy-loss spectroscopy (AR-EELS) and inelastic X-ray scattering spectroscopy (IXS), respectively. Both spectroscopies probe the double differential cross section (DDCS) $\frac{\delta^{2}\sigma}{\delta\Omega\omega′}$, which describes the scattering probability of a particle with final energy $\omega′$ into the solid angle $\Omega$. For the scattering of x-ray photons, the DDCS is proportional to the dynamical structure factor S(q,$\Delta\omega$)
\begin{equation} \tag{1} \left.\frac{\delta^{2}\sigma}{\delta\Omega\omega′}\right\vert_{x-ray} \propto \mathit S(\mathbf q,\Delta\omega), \end{equation}where q = ki − kf is the momentum loss from the initial momentum ki to the final one kf, and $\Delta\omega=\omega_{i}−\omega_{f}$ is the corresponding energy loss.
For the scattering of electrons, the DDCS is proportional to the electron energy loss function $ \mathit L ( \mathbf q,\Delta\omega$):
\begin{equation} \tag{2} \left.\frac{\delta^{2}\sigma}{\delta\Omega\omega′}\right\vert_{electron} \propto \mathit L(\mathbf q,\Delta\omega), \end{equation}In linear response theory, both the dynamical structure factor and the loss function can be expressed in terms of the macroscopic dielectric function $\varepsilon_{M} (\mathbf q,\Delta\epsilon)$. We obtain
\begin{equation} \tag{3} \mathit S(\mathbf q, \Delta\omega) = -\frac{1}{\pi} \nu^{-1}(\mathbf q)â„‘ \Big[\frac{1}{\varepsilon_{M} (\mathbf q,\Delta\omega)}\Big], \end{equation}where $\nu(\mathbf q)$ is the Coulomb potential, and
\begin{equation} \tag{4} \mathit L(\mathbf q, \Delta\omega) = -â„‘ \Big[\frac{1}{\varepsilon_{M} (\mathbf q,\Delta\omega)}\Big]. \end{equation}In the theory part of the tutorial Excited states from BSE, the calculation of the macroscopic dielectric function $\varepsilon_{M} (\omega)$ is discussed in detail. Further in-depth discussion of this theory can be found here.
2. Preliminary Step: Ground-State Calculation¶
Important note: All input parameters are given in atomic units!
As a preliminary step to calculate excited-state properties from BSE, a ground-state calculation has to be performed. In this tutorial we consider as an example LiF. Create a directory named run_LiF_BSE_q as well as a sub-directory GS where the preliminary ground state calculation is performed.
%%bash
mkdir -p run_LiF_BSE_q && cd run_LiF_BSE_q
mkdir -p GS
Inside the GS sub-directory we create the input file for LiF. In the structure element we include the lattice parameter and basis vectors of LiF, which has a rock-salt cubic lattice, as well as the positions of the Li and F atoms. In the groundstate element, we include a 10×10×10 k-point mesh (ngridk) and a value of 14.0 for gmaxvr. This value, which is larger than the default, is needed in view of the excited-state calculation.
Create a new file input.xml that looks like the following:
<input>
<title>LiF-BSE: Ground-State Calculation</title>
<structure speciespath="$EXCITINGROOT/species/">
<crystal scale="7.608">
<basevect>0.5 0.5 0.0</basevect>
<basevect>0.5 0.0 0.5</basevect>
<basevect>0.0 0.5 0.5</basevect>
</crystal>
<species speciesfile="F.xml">
<atom coord="0.5000 0.5000 0.5000" />
</species>
<species speciesfile="Li.xml">
<atom coord="0.0000 0.0000 0.0000" />
</species>
</structure>
<groundstate
do="fromscratch"
ngridk="10 10 10"
xctype="GGA_PBE_SOL"
gmaxvr="14.0"/>
</input>
Replace the $EXCITINGROOT variable with the correct path by using the command:
%%bash
cd run_LiF_BSE_q/GS
python3 -m excitingscripts.setup.excitingroot
cd ../..
Now start the ground-state SCF calculation and check that it finishes gracefully.
%%bash
cd run_LiF_BSE_q/GS
time $EXCITINGROOT/bin/exciting_smp input.xml
cd ../..
In case of a successful run the files STATE.OUT and EFERMI.OUT should be present in the directory. These two files are needed as a starting point for the BSE calculation.
3. Optical BSE Calculation for Finite Momentum Values¶
To start the BSE calculations, move to the parent directory LiF-BSE-q and create a new folder with the name BSE:
%%bash
cd run_LiF_BSE_q
mkdir -p BSE
As anticipated, we need the files STATE.OUT and EFERMI.OUT obtained from the ground state calculation to be inside the new directory. Copy these two files from the GS folder into the BSE one.
%%bash
cd run_LiF_BSE_q/BSE
cp ../GS/{STATE.OUT,EFERMI.OUT,input.xml} ./
cd ../..
In this way, we can skip the ground-state calculation by setting do="skip" in the groundstate element of the input file. To perform an excited-state calculation we must include the xs element in the input file:
...
<xs
xstype="BSE"
ngridk="4 4 4"
vkloff="0.097 0.273 0.493"
ngridq="4 4 4"
nempty="30"
gqmax="3.0"
broad="0.007"
scissor="0.20947"
tevout="true">
<energywindow
intv="0.3 1.0"
points="1200"/>
<screening
screentype="full"
nempty="100"/>
<BSE
bsetype="singlet"
nstlbse="1 5 1 4"
iqmtrange="1 3"/>
<qpointset>
<qpoint>0.00 0.0 0.00</qpoint>
<qpoint>0.25 0.0 0.25</qpoint>
<qpoint>0.50 0.0 0.50</qpoint>
</qpointset>
</xs>
...
The structure of the xs element is nearly identical to the one described in the tutorial Excited states from BSE. The major changes are the following:
- The attribute
iqmtrangein theBSEelement defines the range of momentum transfer vectors q for which the calculation should be performed. The range has to lie within the boundaries of the the elementqpointset. - In the
qpointsetelement, moreqpointsubelements are added, which define momentum transfer vectors for the calculations. Note that the q-vectors are provided in units of the reciprocal lattice vectors!
We can now run the BSE calculation with the following commands:
%%bash
cd run_LiF_BSE_q/BSE
time $EXCITINGROOT/bin/exciting_smp input.xml
cd ../..
Once the calculation completes successfully (it will take a few minutes), a number of files and some folders will be generated in the working directory. Here, we are mostly interested in the calculated loss functions and dynamical structure factors, which can be found in the directory LOSS. The directory contains several files ending with _QMT00x.OUT, where x is the index of the momentum transfer vector. For q=0, the loss function is a 3×3 tensor, and the optical components Lxy (q=0, $\Delta\omega$) are written to the file ending with _OCxy.OUT.
To plot the scattering spectra from your command line, execute the following commands.
%%bash
cd run_LiF_BSE_q/BSE
cp LOSS/LOSS_BSE-singlet-TDA-BAR_SCR-full_OC11.OUT LOSS/QMT_1
cp LOSS/LOSS_BSE-singlet-TDA-BAR_SCR-full_QMT002.OUT LOSS/QMT_2
cp LOSS/LOSS_BSE-singlet-TDA-BAR_SCR-full_QMT003.OUT LOSS/QMT_3
python3 -m excitingscripts.plot.files -d LOSS -f QMT_1 QMT_2 QMT_3 -lx 'Energy [eV]' -ly 'Loss function' -x 8 27 -y 0 8 -lp 2
cd ../..
Notice that increasing the transferred momentum shifts the first peak to higher energies and reduces the oscillator strength. The pronounced maximum above 25 eV, which dominates the spectrum in the optical limit, is quenched at finite momentum transfer.
4. X-Ray BSE Calculation for Finite Momentum Values ¶
Equivalently to the calculation of the optical loss function above, we will now calculate the core loss function for excitations from the F 1s states, i.e., the F K edge of LiF. It is assumed that you are already familiar with this type of excited-state calculations, having performed the tutorial X-ray Absorption Spectra Using BSE. We move back to the parent directory LiF-BSE-q directory and create a new sub-directory NRIXS.
%%bash
cd run_LiF_BSE_q
mkdir -p NRIXS
As usual, we need to copy files from the GS folder into the NRIXS one, to avoid repeating the ground-state calculation.
%%bash
cd run_LiF_BSE_q/NRIXS
cp ../GS/{STATE.OUT,EFERMI.OUT,input.xml} ./
cd ../..
The ground-state calculation is skipped by setting do="skip" in the groundstate element of the input file. Then, we add the following xs subelement into the input file:
...
<xs
xstype="BSE"
ngridk="4 4 4"
vkloff="0.097 0.273 0.493"
ngridq="4 4 4"
nempty="30"
gqmax="3.0"
broad="0.007"
scissor="0.20947"
tevout="true">
<energywindow
intv="0.3 1.0"
points="1200"/>
<screening
screentype="full"
nempty="100"/>
<BSE
bsetype="singlet"
nstlbse="1 5 1 4"
iqmtrange="1 3"/>
<qpointset>
<qpoint>0.00 0.0 0.00</qpoint>
<qpoint>0.25 0.0 0.25</qpoint>
<qpoint>0.50 0.0 0.50</qpoint>
</qpointset>
</xs>
...
The input file is nearly identical to the one used for optical calculations above. We have simply added the attribute xas to trigger a core-level calculation, and we have specified which core states we want to excite using the attributes xasspecies, xasatom, and xasedge. The range of unoccupied states is now defined by the attribute nstlxas. Note that we also have to modify the attribute energywindow, since core excitations occur at much higher energies than the optical ones.
Run the BSE calculation with the following commands:
%%bash
cd run_LiF_BSE_q/NRIXS
time $EXCITINGROOT/bin/exciting_smp input.xml
cd ../..
After the run is successfully completed (it will take a few minutes), the usual output files and sub-folders will appear in the working directory.
Again, we are interested in the loss function, which is found inside the folder LOSS. The calculated loss function can be visualized from your command line with the same script used above:
%%bash
cd run_LiF_BSE_q/NRIXS
cp LOSS/LOSS_BSE-singlet-TDA-BAR_SCR-full_OC11.OUT LOSS/QMT_1
cp LOSS/LOSS_BSE-singlet-TDA-BAR_SCR-full_QMT002.OUT LOSS/QMT_2
cp LOSS/LOSS_BSE-singlet-TDA-BAR_SCR-full_QMT003.OUT LOSS/QMT_3
python3 -m excitingscripts.plot.files -d LOSS -f QMT_1 QMT_2 QMT_3 -lx 'Energy [eV]' -ly 'Loss function' -x 655 675
We observe less dependence of the loss function on momentum transfer than we saw in the previous case for the optical region. The most relevant changes in the spectrum with increasing momentum transfer are the emergence of a pre-peak that is dark at q = 0 and the increasing intensity of the first strong peak.
5. Exercise¶
- Modify the attribute
bsetypeto calculate the independent-particle (IP), random-phase approximation RPA and triplet loss function spectrum in the optical range for LiF. How is the spectrum modified? If you have already computed the singlet spectra, you only need to execute the last task of the BSE program flow, i.e., "bse", when switching thebsetype.