10.9. Electron-phonon interaction

The electron-phonon (eph) interaction is an important mechanism that affects a broad range of materials properties. For example, it is one the main scattering mechanisms in charge and heat transport, and it determines the lifetime of electrons and phonons. DFPT provides a way to obtain scattering rates due to electron-phonon interactions given by electron-phonon interaction (EPI) matrix, \(g_{\mathbf{k+q}m,\mathbf{k}n}^{\mathbf{q}\nu}\), which is the central quantity that characterizes eph interactions. Once EPI matrix is known, related physical properties can be calculated through various formalisms. In this section, we demonstrate simulating eph interactions in RESCU by looking at aluminum example and calculating electrical consuctivity.

EPI matrix for phonons with frequency \(\omega_{\mathbf{q}\nu}\) and polarization vector \(\eta_{\mathbf{q}\nu}\) and electronic bloch states of \(| \mathbf{k}+\mathbf{q}m \rangle\) and \(| \mathbf{k}n \rangle\) is given by [SS96]:

(10.9.1)\[g_{\mathbf{k}+\mathbf{q}m,\mathbf{k}n}^{\mathbf{q}\nu} = \sum_{R,\mu} \frac{\eta_{\mathbf{q}\nu}(R\mu)}{(M_{R} \omega_{\mathbf{q}\nu})^{1/2}} \langle \mathbf{k}+\mathbf{q}m| \frac{\delta V_{eff}}{\delta R_{\mu}} | \mathbf{k}n \rangle\]

where \(M_{R}\) is the nuclei mass. Therefore, following quantities are required for computing EPI matrix:

  • Ground state wavefunctions

  • Phonon modes

  • Perturbed self-consistent potentials for atomic displacements perturbations

In RESCU these quantities are obtained using the following workflow:

  • Performing a ground state calculation to obtain the ground state density and unperturbed wavefunctions

  • Performing a phonon (and dielectric in case of non-metallic systems) DFPT calculation to obtain dynamical matrix (dielectric constant and Born effective charges) and Perturbed self-consistent potentials

  • Performing a DFPT electron phonon calculation to obtain EPI matrix and calculate other eph related quantities


Electron-phonon calculations often require much denser k- and q-sampling compared to ground state and phonon (dielectric) DFPT calculations. RESCU uses a nonself-consistent calculation to compute quantities on dense k-mesh, and Fourier interpolation to compute dynamical matrix and perturbed self-consistent potentials on dense q-mesh, to speed up the calculations at eph level [B20].


k- and q-mesh used in this example are only for demonstration and does note generate accurate physical results

We begin by calculating the ground state of Al using the following input file:

info.calculationType = 'self-consistent'
info.savepath = './results/al_real_scf'
element(1).species = 'Al'
element(1).path = './Al_AtomicData.mat'
atom.element = [1]
atom.xyz = [0 0 0]
domain.latvec = 7.665*[0.5,0.5,0.0;0.0,0.5,0.5;0.5,0.0,0.5]
domain.lowres = 0.2
functional.libxc = true
functional.list = {'XC_LDA_X','XC_LDA_C_PW'}
kpoint.gridn = [3 3 3]
kpoint.sampling = 'gauss'   % required for dfpt for metals
kpoint.sigma = 0.025        % required for dfpt for metals
mixing.tol = 1e-8*[1 1]     % advised for dfpt
option.saveWavefunction = 1 % required for dfpt phonon (and dielectric)
diffop.method = 'fft'       % required for dfpt
domain.fourierInit = false  % required for dfpt

RESCU can be executed as

rescu -i al_real_scf

Next, we perform a DFPT phonon calculation as follows. Note that this time we also ask RESCU to save perturbed self-consistent potentials which are to be used later in electron-phonon step:

info.calculationType = 'dfpt-phonon'
info.savepath = './results/al_real_phonon'
rho.in = './results/al_real_scf.mat'
psi.in = './results/al_real_scf.h5'
mixing.tol = 1e-5*[1 1]       % advised for dfpt
option.saveWavefunction = 0   % No need to save perturbed wavefunctions
option.savePotential = 1      % MUST be saved to be used in eph calc
symmetry.timereversal   = 1
symmetry.spacesymmetry  = 0   % required for eph calc
symmetry.pointsymmetry  = 0   % required for eph calc


Perturbed self-consistent potentials for atomic displacements perturbation are required in electron-phonon calculations, hence they must be stored using:

option.savePotential = 1


RESCU currently supports timereversal symmetry for calculating perturbed self-consistent potentials. It must be ensured other symmetry settings are deactivated if phonon results are to be used later for an electron-phonon calculation:

symmetry.timereversal   = 1
symmetry.spacesymmetry  = 0
symmetry.pointsymmetry  = 0

RESCU is executed by

rescu -i al_real_phonon

The important outputs of this step are already explained in sections Dynamical matrix and Born effective charges and Metals. It should be noted that for system with an electronic band-gap we must perform a dielectric calculation before taking the next step. This provides additional quantities which are required to calculate LO-TO effect corresponding to phonons in long wave-length limit, e.g see Phonon band structure, [B20] and [GL97].

Now we can proceed with the electron-phonon calculation using the following input file

info.calculationType = 'dfpt-elph'
info.savepath = './results/al_real_elph'
rho.in = './results/al_real_scf.mat'
dfpt.phononData = './results/al_real_phonon.mat' % path to DFPT phonon data
kpoint.gridn = [5,5,5]
dfpt.qpointGridn = [5,5,5]
dfpt.elphTemperature = linspace(5,300,6)    % required for elph
dfpt.elphZcut = 0.01                        % required for elph
units.energy = 'eV'

dfpt.saveEPImat = 0             % Optional, default value 0
dfpt.saveElectronVelocity = 0   % Optional, default value 0
dfpt.saveElectronSE = 0         % Optional, default value 0
dfpt.saveElectronTau   = 0      % Optional, default value 0

Description of important keywords follows:

  • info.calculationType determines the type of the calculation;;

  • info.savepath points to the file which stores the output data of the elph calculation;

  • rho.in this DFPT calculation loads the ground state electronic density at the beginning of the computation, and this keyword points to the file which contains the output of the ground state calculation, so that RESCU can load corresponding electronic density data;

  • dfpt.phononData points to the path of the file containing the dynamical matrix (and Born effective charges);

  • kpoint.grid determines the k-sampling of the reciprocal space used for nonself-consistent calculation of ground state wavefunctions and electronic structure;

  • dfpt.qpointGridn determines the q-sampling of the reciprocal space used in interplation of dynamical matrix for phonons and perturbed self-consistent potentials;


Current version of RESCU requires an identical, \(\Gamma\)-centered, and odd kpoint.gridn and dfpt.qpointGridn for DFPT calculations;

  • dfpt.elphZcut small imaginary parameter used in evaluating self-energy integrals. Here, it is combined with units.energy, so that it can be set in eV units;

  • dfpt.elphTemperature a row vector containing temperatures at which eph related variables are calculated;

Some of the large eph related outputs can be saved using following keywords:

  • dfpt.saveEPImat optional keyword to save EPI matrix to ./results/al_real_elph.h file. Default value is set to 0. Results can be found in dfpt/EPImat field;

  • dfpt.saveElectronVelocity optional keyword to save electron velocities to ./results/al_real_elph.h file. Default value is set to 0. Results can be found in dfpt/electronVelocitiy field;

  • dfpt.saveElectronSE optional keyword to save electron self-energy to ./results/al_real_elph.h file. Default value is set to 0. Results can be found in dfpt/electronSE field;

  • dfpt.saveElectronTau optional keyword to save electron lifetime to ./results/al_real_elph.h file. Default value is set to 0. Results can be found in dfpt/electronTau field.


Saving elph variables in h5 file may require a large disk space.

RESCU is executed as follows

rescu -i al_real_elph

Important output variables can be found in dfpt structure of ./results/bn_real_phonon.mat file:

>> load results/al_real_elph.mat
>> dfpt

      dfpt =

        struct with fields:

                          BornEC: [3x3 double]
                             IFC: [3x3x27 double]
                             asr: [1 1 1]
                            desc: [1x1 struct]
                dielectricTensor: [3x3 double]
                          dynMat: [3x3x125 double]
          electricalConductivity: [3x3x6 double]
                 elphTemperature: [5 64 123 182 241 300]
                  elphTransCoeff: [3x3x6 double]
                      phononFreq: [3x125 double]
                      phononMode: {125x1 cell}
                    qpointDirect: [125x3 double]
                     qpointGridn: [5 5 5]
                            type: 'elph'

RESCU computes generalized transport coefficients using relaxation time approximation (RTA) of linearized Boltzmann transport equations (LBTE), where electrons lifetimes are evaluated from electrons self-energy due to eph interactions [B20]:

(10.9.2)\[\begin{split}\begin{split} \Sigma_{n\mathbf{k}}^{FM}(\omega, \varepsilon_{F}, T) = & \sum_{m, \nu} \int_{BZ} \frac{d\mathbf{q}}{\Omega_{BZ}} |g_{mn\nu}(\mathbf{k}, \mathbf{q})|^{2} \\ & \times \Big[ \frac{n_{\mathbf{q}\nu}(T) + f_{m\mathbf{k+q}}(\varepsilon_{F},T)}{\omega-\varepsilon_{m\mathbf{k+q}}+\omega_{\mathbf{q}\nu} + i\eta} \\ & + \frac{n_{\mathbf{q}\nu}(T) + 1 - f_{m\mathbf{k+q}}(\varepsilon_{F},T)}{\omega-\varepsilon_{m\mathbf{k+q}} - \omega_{\mathbf{q}\nu} + i\eta}\Big] \end{split}\end{split}\]

where \(f_{m\mathbf{k+q}}(\varepsilon_{F},T)\) and \(n_{\mathbf{q}\nu}(T)\) correspond to the Fermi-Dirac and Bose-Einstein distributions. In the \(\eta \rightarrow 0^{+}\) limit, the imaginary part evaluated at the Kohn-Sham energy can be obtained by

(10.9.3)\[\begin{split}\begin{split} \lim_{\eta \rightarrow 0^{+}} \Im\{ \Sigma_{n\mathbf{k}}^{FM}(\varepsilon_{F})\} = & \pi \sum_{m, \nu} \int_{BZ} \frac{d\mathbf{q}}{\Omega_{BZ}} |g_{mn\nu}(\mathbf{k}, \mathbf{q})|^{2} \\ & \times \Big[ ( n_{\mathbf{q}\nu} + f_{m\mathbf{k+q}}) \delta(\varepsilon_{n\mathbf{k}} - \varepsilon_{m\mathbf{k+q}} + \omega_{\mathbf{q}\nu} ) \\ & + ( n_{\mathbf{q}\nu} + 1 - f_{m\mathbf{k+q}} ) \delta (\varepsilon_{n\mathbf{k}}-\varepsilon_{m\mathbf{k+q}}-\omega_{\mathbf{q}\nu} ) \Big] \end{split}\end{split}\]

from which electrons lifetimes are calculated:

(10.9.4)\[\begin{aligned} \frac{1}{\tau_{n\mathbf{k}}} = 2 \lim_{\eta \rightarrow 0^{+}} \Im\{ \Sigma_{n\mathbf{k}}^{FM}(\varepsilon_{F})\}.\end{aligned}\]

These coefficients are stored in elphTransCoeff as a \(3\times 3\times N_{temperature}\) array in atomic units:

>> dfpt.elphTransCoeff(:,:,1)

      ans =

         1.0e+04 *

          1.4984   -0.0020   -0.0020
         -0.0020    1.4984   -0.0020
         -0.0020   -0.0020    1.4984

Also, electrical conductivity is stored as a \(3\times 3\times N_{temperature}\) array in (Siemens cm-1) units in electricalConductivity:

>> dfpt.electricalConductivity(:,:,1)

      ans =

         1.0e+07 *

          1.2244   -0.0016   -0.0016
         -0.0016    1.2244   -0.0016
         -0.0016   -0.0016    1.2244

For small data sizes, RESCU also outputs the following quantities in text file format:

  • Electrons lifetimes in results/Lifetime_electron

# Electrons lifetime (femto-seconds) with SERTA approximation
# nbands = 10
# T =    5.0 (K)      Efermi = -4.338795e-02 (Ha)
# spin =  1
# kpoint = [0.0000000   0.0000000   0.0000000]
+3.830255e-27 +6.056664e-28 +4.693947e-28 +7.562906e-28 +1.244860e-30 +6.569030e-31 +4.699217e-31 +1.647781e-28 +1.298661e-28 +2.911374e-31
# kpoint = [0.2000000   0.0000000   0.0000000]
  • Electrical conductivity in results/Conductivity_electron

# Electrical conductivity (Siemens/cm^-1) with SERTA approximation
# T =    5.0 (K)      Efermi = -4.338795e-02 (Ha)
+1.224412e+07   -1.645853e+04   -1.646234e+04
-1.645853e+04   +1.224413e+07   -1.647936e+04
-1.646234e+04   -1.647936e+04   +1.224416e+07
  • Electrons self-energy in results/Selfenergy_electron

# Electrons self-energy (Hartree) with SERTA approximation
# nbands = 10
# T =    5.0 (K)      Efermi = -4.338795e-02 (Ha)
# spin =  1
# kpoint = [0.0000000   0.0000000   0.0000000]
-2.845214e-35 -7.637874e-38 +2.644990e-36 -4.830218e-37 -3.199617e-35 -6.232497e-37 -1.452827e-35 -3.868223e-37 +2.390173e-34 -2.350065e-34
-6.871185e-35 -4.453474e-34 -6.396189e-34 -6.225507e-34 -2.519247e-35 -1.775418e-36 -6.435553e-35 -2.252705e-36 +8.060548e-34 -1.004852e-33
# kpoint = [0.2000000   0.0000000   0.0000000]

As mentioned above, these quantities as well as EPI matrix can also be saved in h5 format. For example EPI matrix information can be accessed using:

% (nbnad,nband) EPI matrix for [i_mode,i_qpt,i_kpt_ispin] = ind2sub([nmode,nqpt,nkpt,nspin],10)
>> h5info('results/al_real_elph.h5','/dfpt/EPImat10')
      ans =

        struct with fields:

            Filename: '/home/saeedb/Projects/nano/rescumat/doc-examples/al-eph/results/al_real_elph.h5'
                Name: 'EPImat10'
            Datatype: [1x1 struct]
           Dataspace: [1x1 struct]
           ChunkSize: []
           FillValue: 0
             Filters: []
          Attributes: [1x1 struct]

% real and imaginary parts
>> data = h5read('results/al_real_elph.h5','/dfpt/EPImat10');
>> size(data)

      ans =

          10    20
% (nband,nband) transition rates
>> g = data(:,1:10) + 1j*data(:,11:20)
>> size(g)

   ans =

       10    10