Raman spectroscopy
GPAW offers two ways of calculating Raman intensities. One can use the ASE Raman utility together with the GPAW LRTDDFT module as shown in the Resonant Raman tutorial (Resonant-)Raman spectra of gas-phase water.
GPAW also implements Raman spectroscopy for zone-center phonons of extended systems using the electron-phonon coupling (see Electron-Phonon Coupling Theory and Electron-phonon coupling) within 3rd order perturbation theory Taghizadeh et a. [1] , which is discussed here. This method is currently only implementated for the LCAO mode.
The Stokes Raman intensity can be written as
where \(\nu\) denotes phonon modes and \(\alpha\), \(\beta\) denote polarisations of the incoming and outgoing laser light. The Raman tensor \(R_{\alpha \beta}^\nu\) has six terms and is given by Ref. [1] Eq. (10)
The first term is considered to be the resonant term of the expression, the other terms represent different time orderings of the interaction in the Feynman diagrams.
To compute the Raman intensity we need these ingredients: The momentum matrix
elements \(p_{ij}^\alpha=\langle i \mathbf{k} | \hat p^\alpha| j \mathbf{k} \rangle\),
the electron-phonon matrix \(g_{ij}^\nu = \langle i \mathbf{k} \vert \partial_{\nu{q=0}} V^{KS} \vert j \mathbf{k} \rangle\)
in the optical limit \(\mathbf{q}=0\) and of course knowledge of the electronic
states and phonon modes throughout the Brillouin zone.
For these calculations we can employ in the get_momentum_transitions()
method, the GPAW electron-phonon module Electron-phonon coupling and the ASE phonon module, respectively.
It is required to calculate all properties for the full Brillouin zone using \(symmetry=off\). By default the routine saves a file called mom_skvnm.npy
containing the
momentum matrix. This can be deactivated using the savetofile
switch. The
matrix is always the return value of get_momentum_transitions()
.
Energy changes for phonons and potential changes for electron-phonon couplings are both computed using a finite displacement technique. Both quantities can be obtained simultaenously. In principle the phonon modes can be obtained in any fashion, which yields an ASE phonon object though. For small systems the finite displacement method has the disadvantage of leading to an interaction of a displaced atom with its periodic images. This can lead to large errors especially in the electron-phonon matrix. This can be avoided by using a sufficiently large supercell for the finite displacement simulations.
If phonon and effective potential are calculated simultaenously, results are saved in the same cache files with default name \(elph\).
Some more details are elaborated in the related tutorial Raman spectroscopy for extended systems.
References
Code
- class gpaw.elph.ResonantRamanCalculator(calc, wph_w, momname='mom_skvnm.npy', elphname='gsqklnn.npy', raman_name='Rlab')[source]
Resonant Raman Matrix calculator
- Parameters:
- calculate(w_in, d_i, d_o, gamma_l=0.1, limit_sum=False)[source]
Calculate resonant Raman matrix
- Parameters:
- class gpaw.elph.RamanData(raman_name='Rlab', gridspacing=1.0)[source]
Raman Spectroscopy data.
- Parameters:
- calculate_raman_intensity(d_i, d_o, T=300, sigma=5.0)[source]
Calculates Raman intensities from Raman tensor.
Returns bare \(|R|^2\) and and Bose occupation weighted values
- calculate_raman_spectrum(entries, T=300, sigma=5.0)[source]
Calculates Raman intensities from Raman tensor.
Returns Raman shift in eV, bare \(|R|^2\) and Bose occupation weighted \(|R|^2\) values
- gpaw.lcao.dipoletransition.get_momentum_transitions(wfs, savetofile=True)[source]
Finds the momentum matrix elements: <nk|p|mk> = k delta_nm - i <nk|nabla|mk>
- Parameters:
wfs – LCAO WaveFunctions object
savetofile (bool) – Determines whether matrix is written to the file mom_skvnm.npy (default=True)