Exciton-phonon coupling and luminescence hBN: Difference between revisions
| Line 248: | Line 248: | ||
in such a way results are stored in the <code>bse_Lbar</code> folder. | in such a way results are stored in the <code>bse_Lbar</code> folder. | ||
== Step 7: Obtain the electron-phonon matrix elements with LetzElPhC== | == Step 7: Obtain the electron-phonon matrix elements with yambopy and LetzElPhC== | ||
Go in the nscf folrder where it is present the SAVE directory and run | Go in the nscf folrder where it is present the SAVE directory and run | ||
Latest revision as of 13:17, 5 October 2026

This tutorial is alternative to the tutorial Exciton-phonon coupling and luminescence to calculated phonon-assisted luminescence.
(Also this tutorial is based on the following papers: the theory[1][2][3][4][5][6] and applications[7][8][9][10][11][12][13][14] of exciton-phonon physics.
This tutorial is specific for the Claudia cluster in Marseille
In this tutorial we will use hBN as example.
Step 0: Install YamboPy
First of all you have to install YamboPY. You can use the Miniconda module on Claudia
1) Install miniconda in your home directory
https://www.anaconda.com/docs/getting-started/miniconda/install#linux-2
reload bash environment
source ~/.bashrc
2) Download YamboPy from github
git clone https://github.com/yambo-code/yambopy.git
git checkout bug-fixes
3) Install YamboPY
cd yambopy
pip install ./
(if you want to do developer installation you can do pip install -e ./ )
Step 1: scf calculation
In order to study luminescence in a paradigmatic system, we switch to bulk hexagonal boron nitride and we repeat the workflow. As you can easily see, one can think about automatizing the execution of all these calculations via scripting or more advanced tools. However, in the case of very large simulations (memory-limited or disk-space limited) or for systems whose electronic and lattice properties are fragile with respect to tiny calculation details, one must be very careful and run many basic tests.
Fortunately, we are running a fast underconverged example. We use LDA pseudopotentials from the pseudo-dojo library.
You can download the psuedo-potentials here: B_LDA_dojo.UPF, N_LDA_dojo.UPF
Create the scf folder and and run the following job:
1. Input hbn.scf:
&control
calculation='scf',
prefix='hBN',
restart_mode='from_scratch'
pseudo_dir = '$PSEUDO_DIR',
outdir = './'
verbosity = 'high'
wf_collect=.true.
/
&system
ibrav = 4,
celldm(1) = 4.703675849
celldm(3) = 2.603711434
nat= 4,
ntyp= 2,
force_symmorphic=.true.
ecutwfc = 100,
/
&electrons
diago_david_ndim = 2
diago_full_acc=.true.
diago_thr_init=5.0e-6
mixing_mode = 'plain'
mixing_beta = 0.7
conv_thr = 1.0d-16
/
ATOMIC_SPECIES
B 10.81100 B_LDA_dojo.UPF
N 14.00674 N_LDA_dojo.UPF
ATOMIC_POSITIONS {crystal}
N 0.6666666670 0.3333333330 -0.250000000000
B 0.3333333330 0.6666666670 -0.250000000000
B 0.6666666670 0.3333333330 0.25000000000
N 0.3333333330 0.6666666670 0.25000000000
K_POINTS {automatic}
6 6 2 0 0 0
Step 2: nscf calculation for Yambo
Then we run the non self-consistent calculation to get valence and conduction bands that will be our basis for the calculations in Lumen.
Copy the SCF calculation in the nscf folder with the command:
cd nscf cp -r ../scf/hBN.save ./
then run the NSCF with QuantumEspresso:
&control
calculation='nscf',
prefix='hBN',
restart_mode='from_scratch'
pseudo_dir = '$PSEUDO_DIR'
outdir = './'
verbosity = 'high'
wf_collect=.true.
/
&system
ibrav = 4,
celldm(1) = 4.703675849
celldm(3) = 2.603711434
nat= 4,
ntyp= 2,
force_symmorphic=.true.
ecutwfc = 100,
nbnd = 120
/
&electrons
diago_david_ndim = 2
diago_full_acc=.true.
diago_thr_init=5.0e-6
mixing_mode = 'plain'
mixing_beta = 0.7
conv_thr = 1.0d-16
/
ATOMIC_SPECIES
B 10.81100 B_LDA_dojo.UPF
N 14.00674 N_LDA_dojo.UPF
ATOMIC_POSITIONS {crystal}
N 0.6666666670 0.3333333330 -0.250000000000
B 0.3333333330 0.6666666670 -0.250000000000
B 0.6666666670 0.3333333330 0.25000000000
N 0.3333333330 0.6666666670 0.25000000000
K_POINTS {automatic}
6 6 2 0 0 0
Step 3: dvscf phonon calculation
Create the ph folder and copy the scf results also here:
cd ph cp -r ../scf/hBN.save ./
Finally we run phonons on the same grid of the non self-consistent calculations with the intput hbn.ph :
hbn_dvscf &inputph tr2_ph=1.0d-12, verbosity='high' prefix='hBN', fildvscf = 'hBN-dvscf', electron_phonon = 'dvscf', fildyn='hBN.dyn', epsil=.false., ldisp=.true., recover=.true., nq1=6, nq2=6, nq3=2 /
Step 4: create Yambo SAVE directory
Go in the nscf folder and read the wave-functions with the command:
p2y
and then run the setup
yambo
Step 5: run a BSE calculation
Input bse.in (we include 2 valence and 2 conduction bands),
This input was generated with the command: yambo -X s -o b -k sex -y d -V resp
em1s # [R][Xs] Statically Screened Interaction optics # [R] Linear Response optical properties bss # [R] BSE solver bse # [R][BSE] Bethe Salpeter Equation. dipoles # [R] Oscillator strenghts (or dipoles) NLogCPUs=0 # [PARALLEL] Live-timing CPU`s (0 for all) PAR_def_mode= "balanced" # [PARALLEL] Default distribution mode DIP_CPU= "" # [PARALLEL] CPUs for each role DIP_ROLEs= "" # [PARALLEL] CPUs roles (k,c,v) DIP_Threads=0 # [OPENMP/X] Number of threads for dipoles X_and_IO_CPU= "" # [PARALLEL] CPUs for each role X_and_IO_ROLEs= "" # [PARALLEL] CPUs roles (q,g,k,c,v) X_and_IO_nCPU_LinAlg_INV=1 # [PARALLEL] CPUs for Linear Algebra (if -1 it is automatically set) X_Threads=0 # [OPENMP/X] Number of threads for response functions BS_CPU= "2 2 4" # [PARALLEL] CPUs for each role BS_ROLEs= "k eh t" # [PARALLEL] CPUs roles (k,eh,t) BS_nCPU_LinAlg_INV=1 # [PARALLEL] CPUs for Linear Algebra (if -1 it is automatically set) BS_nCPU_LinAlg_DIAGO=1 # [PARALLEL] CPUs for Linear Algebra (if -1 it is automatically set) K_Threads=0 # [OPENMP/BSK] Number of threads for response functions Chimod= "HARTREE" # [X] IP/Hartree/ALDA/LRC/PF/BSfxc BSKmod= "SEX" # [BSE] IP/Hartree/HF/ALDA/SEX/BSfxc BSEmod= "resonant" # [BSE] resonant/retarded/coupling BSSmod= "d" # [BSS] (h)aydock/(d)iagonalization/(s)lepc/(i)nversion/(t)ddft` BSENGexx= 31583 RL # [BSK] Exchange components BSENGBlk=-1 RL # [BSK] Screened interaction block size [if -1 uses all the G-vectors of W(q,G,Gp)] #WehCpl # [BSK] eh interaction included also in coupling Lkind= "full" # [BSE,X] bar(default)/full/tilde % KfnQP_E 1.259970 | 1.088160 | 1.126830 | # [EXTQP BSK BSS] E parameters (c/v) eV|adim|adim % EXCTemp= 0.000000 eV # Excitonic Temperature (for luminescence spectra and EXCPH) % BSEQptR 1 | 14 | # [BSK] Transferred momenta range % % BSEBands 7 | 10 | # [BSK] Bands range % % BEnRange 0.500000 | 8.000000 | eV # [BSS] Energy range % % BDmRange 0.050000 | 0.050000 | eV # [BSS] Damping range % BEnSteps= 1000 # [BSS] Energy steps % BLongDir 1.000000 | 1.000000 | 0.000000 | # [BSS] [cc] Electric Field versor % BSEprop= "abs" # [BSS] Can be any among abs/jdos/kerr/asymm/anHAll/magn/dich BSEdips= "none" # [BSS] Can be "trace/none" or "xy/xz/yz" to define off-diagonal WRbsWF # [BSS] Write to disk excitonic the WFs BSSldiaLib= "s" # [BSS] Library used with ldiago solver `(s)calapack/(e)lpa` BSSNEig=0 # [SLEPC/DIAGO] Number of eigenvalues to compute % BndsRnXs 1 | 120 | # [Xs] Polarization function bands % NGsBlkXs= 10000 mRy # [Xs] Response block size % LongDrXs 1.000000 | 1.000000 | 0.000000 | # [Xs] [cc] Electric Field %
run the BSE in such a way to save rusults in the bse_LFull folder.
nohup mpirun -np 8 yambo -F bse.in -J bse_Lfull -C bse_Lfull
Step 6: Run BSE for Lbar
Importantly, since we want to describe the phonon-assisted recombination process of an *optical* exciton (i.e., emitting a transverse photon), this time we also run an additional calculation at `Q=0` omitting the nonanalytic long-range Coulomb exchange. Make a second input bse_Lbar.in with the following changes:
Lkind= "Lbar" # [BSE] Lbar (default) / full
% BSEQptR
1 | 1 | # [BSK] Transferred momenta range
%
4b. So now we make a second BSE run in a different directory specified by -J. Here, we also pass to yambo the directory of the previous run as it includes the important screening databases ndb.em1s* that we do not want to recompute from scratch.
mpirun -np 8 yambo -F bse_Lbar.in -J bse_Lbar,bse_Lfull -C bse_Lbar
in such a way results are stored in the bse_Lbar folder.
Step 7: Obtain the electron-phonon matrix elements with yambopy and LetzElPhC
Go in the nscf folrder where it is present the SAVE directory and run
module load yambo/5.3.0
yambopy l2y -ph ../../ph/hbn.ph -b 7 10 -par 1 1 -lelphc /home/attaccalite/SOFTWARE/LetzElPhC/src/lelphc --no_gkkp
Important points:
- We added the flag --no_gkkp to not produce Yambo databases (this part still crashes for me)
- The bands in LetzElphC are between 7 and 10 because LetzElPhC use the bands convention of QE starting from 1, while YamboPY starts from 0.
If you have a big system you can run this part in parallel, for example if you want to run on 8 cores:
#!/bin/bash
#SBATCH -J ypy
#SBATCH -o %x.%j.out
#SBATCH -e %x.%j.status
#SBATCH --cpus-per-task=3
#SBATCH --ntasks 8
#SBATCH -p c6
# Module setup
module purge
module load gnu12/12.2.0
module load openmpi4/4.1.5
module load phdf5/1.14.0
module load netcdf/4.9.0
module load netcdf-fortran/4.6.0
module load scalapack/2.2.2
module load fftw/3.3.10
module load openblas/0.3.21
module load petsc/3.22.0_s
module load slepc/3.22.0_s
module load libxc/5.2.3
module load pmix/4.2.9
module load prun/2.2
# Required parameters
export OMP_NUM_THREADS=${SLURM_CPUS_PER_TASK:-1}
# Report on setup
echo "Starting at $(date)"
# echo " QE PATH: ${QE_ROOT}"
echo " Running on partition: ${SLURM_JOB_PARTITION}"
echo " Running on machines: ${SLURM_JOB_NODELIST}"
module list
ulimit -l unlimited
# Start actual computation
echo "Starting Code"
yambopy l2y -ph ../../ph/hbn.ph -b 7 10 -par 4 2 -lelphc /home/attaccalite/SOFTWARE/LetzElPhC/src/lelphc --no_gkkp
echo "Finished at $(date)"
do not forget that the number of tasks ntasks=nkpool x nqpool in our example 8 = 4 x 2
Step 8: Compute phonon-assisted luminescence
Now you can compute the luminescence with the script:
import numpy as np
import matplotlib.pyplot as plt
from yambopy.exciton_phonon.excph_luminescence import exc_ph_luminescence
from yambopy.exciton_phonon.excph_input_data import exc_ph_get_inputs
path = './'
# Path to BSE calculation (Lout--> response is Lfull)
bsepath = f'{path}/bse_Lfull' # ndb.BS_diago_Q* databases are needed
# Path to BSE calculation for optically active exciton (Lin --> response is Lbar)
bseBARpath = f'{path}/bse_Lbar' # ndb.BS_diago_Q1 database is needed
# Path to electron-phonon calculation
elphpath = path # ndb.elph is needed
# Path to unprojected dipoles matrix elements (optional)
dipolespath = bsepath # ndb.dipoles is needed (optional)
# Path to lattice and k-space info
savepath = f'{path}/SAVE' # ns.db1 database is needed
bands_range=[6,10] # 2 valence, 2 conduction bands
phonons_range=[0,12] # All phonons
nexc_out = 12 # 12 excitonic states at each momentum (Lout)
nexc_in = 12 # 12 excitonic states at Q=0 (Lin)
T_ph = 10 # Lattice temperature
T_exc = 10 # Effective excitonic temperature
emin=4.2 # Energy range and plot details (in eV)
emax=5.0
estep=0.0002
broad = 0.005 # Broadening parameter for peak width (in eV)
# We calculate and load all the inputs:
# * Exciton-phonon matrix elements
# * Excitonic dipole matrix elements
# * Exciton energies
# * Phonon energies
# * We specify bse_path2=bseBARpath meaning we use Lbar calculation for Q=0 excitons
input_data = exc_ph_get_inputs(savepath,elphpath,bsepath,\
bse_path2=bseBARpath,dipoles_path=dipolespath,\
nexc_in=12,nexc_out=12,\
bands_range=[6,10],phonons_range=[0,12])
ph_energies, exc_energies, exc_energies_in, G, exc_dipoles = input_data
# We calculate the luminescence spectrum including the input data from before
w,PL = exc_ph_luminescence(T_ph,ph_energies,exc_energies,exc_dipoles,G,\
exc_energies_in=exc_energies_in,exc_temp=T_exc,\
nexc_out=nexc_out,nexc_in=nexc_in,emin=emin,emax=emax,\
estep=estep,broad=broad)
# Save to file
data = np.column_stack((w, PL))
np.savetxt("hBN_luminescence_12x12x1.dat", data, fmt="%.8f")
and the result is displayed here:

Step 9: Exciton lifetimes
TO BE DONE!!!!
Using the exciton-phonon coupling it is possible also to calculate exciton life-time due to the scattering with phonon. For this kind of calculation input file can be generated with the command yambo -excph l:
excph # [R] Exction-phonon ExcGkkp # [R][EXCPH] Exciton-Phonon Matrix Elelements ExcPhLifeT # [R][EXCPH] Exciton-Phonon Life-Times BoseTemp= 0.000000 eV # Bosonic Temperature % ELPhExcStates 1 | 4 | # [EXCPH] Incoming (external) exciton states % % ELPhExcSum 1 | 8 | # [EXCPH] Outgoing (virtual) exciton states % LoutPath= "none" # [EXCPH] Path of the outgoing L FANdEtresh= 0.100000E-5 eV # [ELPH] Energy treshold for Fan denominator % ElPhModes 1 | 12 | # [ELPH] Phonon modes included % EXCTemp= 0.000000 eV # [EXCPH] Excitonic Temperature (for luminescence spectra) AlphaQ= 0.000000 # [EXCPH] Excitonic band structure 2D distortion #ExcPhOffDiago # [EXCPH] Exciton-Phonon off-diagonal self-energy #NoMatrxEl # [EXCPH] Set all exciton-phonon matrix elements to one (for testing purpose) LDamping= 0.1 meV # [EXCPH] Damping of exc-ph self-energy
notice that in this case LDamping is the broadening of the exciton-phonon self energy and it should be a very small value of the order of the phonon life-times around the meV.[15]. Running this input you will get:
<-> Exciton lifetimes: <-> Exciton [1] Lifetime: xxx meV <-> Exciton [2] Lifetime: xxx meV <-> Exciton [3] Lifetime: xxx meV <-> Exciton [4] Lifetime: xxx meV
the life-time of the first four excitons. Note that this result is not converged respect to the Q-points. In order to speedup convergence that double-grid[11] can be used also in the life-time calculations, see tutorial on electron-phonon coupling for tricks on the convergence of the double-grid.
References
- ↑ Toyozawa, Yutaka, and Chris Oxlade, Optical processes in solids, Cambridge University Press, (2003).
- ↑ G. Antonius, S. G. Louie, Theory of exciton-phonon coupling, Phys. Rev. B, 105, 085111 (2022); arXiv1705.04245 (2017)
- ↑ P. Cudazzo, First-principles description of the exciton-phonon interaction: A cumulant approach, Phys. Rev. B, 102, 045136 (2020); Open access pdf from Luxembourg University
- ↑ F. Paleari, First-principles approaches to the description of indirect absorption and luminescence spectroscopy: exciton-phonon coupling in hexagonal boron nitride, PhD thesis, University of Luxembourg (2019)
- ↑ F. Paleari, and A. Marini, Exciton-phonon interaction calls for a revision of the “exciton” concept, Phys. Rev. B, 106, 125403 (2022); arXiv 2205.02783
- ↑ P. Lechifflart, Exciton-phonon coupling and phonon-assisted luminescence in hexagonal Boron Nitride nanostructures, PhD Thesis, University of Marseille (2023); From the yambo website
- ↑ F. Paleari et al., Exciton-Phonon Coupling in the Ultraviolet Absorption and Emission Spectra of Bulk Hexagonal Boron Nitride, Phys. Rev. Lett. 122, 187401 (2019); arXiv1810.089776
- ↑ E. Cannuccia, B. Monserrat and C. Attaccalite, Theory of phonon-assisted luminescence in solids: Application to hexagonal boron nitride, Phys. Rev. B 99, 081109(R) (2019); arXiv1807.11797
- ↑ Exciton-Phonon Interaction and Relaxation Times from First Principles, Hsiao-Yi Chen, Davide Sangalli, and Marco Bernardi, Phys. Rev. Lett. 125, 107401 (2020); arXiv 2002.08913 (2020)
- ↑ M. Zanfrognini et al., Distinguishing different stackings in layered materials via luminescence spectroscopy, Phys. Rev. Lett. 131, 206902 (2023); arXiv 2305.17554
- ↑ 11.0 11.1 P. Lechifflart, F. Paleari, D. Sangalli, C. Attaccalite, First-principles study of luminescence in hexagonal boron nitride single layer: Exciton-phonon coupling and the role of substrate, Phys. Rev. M, 7 (2), 024006 (2023); arXiv2212.1047
- ↑ Y-h Chan, J. B. Haber, M. H. Naik, J. B. Neaton, D. Y. Qiu, F. H. da Jornada, S. G. Louie, Exciton Lifetime and Optical Line Width Profile via Exciton–Phonon Interactions: Theory and First-Principles Calculations for Monolayer MoS2, Nano Lett., 23, 9 (2023); arXiv 2212.08451 (2023)
- ↑ G. Marini, M. Calandra, P. Cudazzo, Optical absorption and photoluminescence of single layer boron nitride from a first principles cumulant approach, Nano Lett., 24, 20, 6017 (2024); arXiv 2402.03826 (2024)
- ↑ M. Nalabothula, S. Reichardt, L. Wirtz, Origin of Interlayer Exciton–Phonon Coupling in 2D Heterostructures, Nano Lett., 25, 15 (2025), arXiv 2407.16111 (2025)
- ↑ First-principles calculations of phonon frequencies, lifetimes, and spectral functions from weak to strong an-harmonicity: The example of palladium hydrides, L. Paulatto, et al. Phys. Rev. B 91, 054304 (2015)