Exciton-phonon coupling and luminescence hBN: Difference between revisions

From Lumen wiki
Jump to navigationJump to search
Davide (talk | contribs)
 
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

  1. ↑ Toyozawa, Yutaka, and Chris Oxlade, Optical processes in solids, Cambridge University Press, (2003).
  2. ↑ G. Antonius, S. G. Louie, Theory of exciton-phonon coupling, Phys. Rev. B, 105, 085111 (2022); arXiv1705.04245 (2017)
  3. ↑ 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
  4. ↑ 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)
  5. ↑ 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
  6. ↑ P. Lechifflart, Exciton-phonon coupling and phonon-assisted luminescence in hexagonal Boron Nitride nanostructures, PhD Thesis, University of Marseille (2023); From the yambo website
  7. ↑ 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
  8. ↑ 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
  9. ↑ 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)
  10. ↑ M. Zanfrognini et al., Distinguishing different stackings in layered materials via luminescence spectroscopy, Phys. Rev. Lett. 131, 206902 (2023); arXiv 2305.17554
  11. ↑ 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
  12. ↑ 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)
  13. ↑ 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)
  14. ↑ M. Nalabothula, S. Reichardt, L. Wirtz, Origin of Interlayer Exciton–Phonon Coupling in 2D Heterostructures, Nano Lett., 25, 15 (2025), arXiv 2407.16111 (2025)
  15. ↑ 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)