ElectronPhononCoupling

class ElectronPhononCoupling(configuration, hamiltonian_derivatives, dynamical_matrix, kpoints_fractional=None, kpoints_cartesian=None, qpoints_fractional=None, qpoints_cartesian=None, electron_bands=None, phonon_modes=None, calculate_velocities=None, max_interaction_range=None, rotate_to_pure_spin_states=None, energy_tolerance=None, initial_state_energy_range=None, store_dense_coupling_matrices=None, coupling_matrix_tolerance=None)

Class for calculating the electron-phonon coupling matrix for a BulkConfiguration.

Parameters:
  • configuration (BulkConfiguration) – The BulkConfiguration for which to calculate the electron-phonon coupling.
  • hamiltonian_derivatives (HamiltonianDerivatives) – The HamiltonianDerivatives for the configuration.
  • dynamical_matrix (DynamicalMatrix) – The DynamicalMatrix for the configuration.
  • kpoints_fractional (MonkhorstPackGrid | RegularKpointGrid | list of list of float) – List of fractional k-points. This option is mutually exclusive to kpoints_cartesian, i.e. [[ 0.0, 0.0, 0.0], [0.5,0.0,0.0]].
    Default: [[0.0, 0.0, 0.0]].
  • kpoints_cartesian (PhysicalQuantity of type inverse length) – List of Cartesian k-points, i.e. [[ 0.0, 0.0, 0.0], [0.5,0.0,0.0]] * Angstrom**-1. This option is mutually exclusive to kpoints_fractional.
  • qpoints_fractional (MonkhorstPackGrid | RegularKpointGrid | list of list of float) – List of fractional q-points, i.e. [[0.0, 0.0, 0.0], [0.5, 0.0, 0.0]]. This option is mutually exclusive to qpoints_cartesian.
    Default: [[0.0, 0.0, 0.0]].
  • qpoints_cartesian (PhysicalQuantity of type inverse length) – List of Cartesian q-points i.e. [[0.0, 0.0, 0.0], [0.5, 0.0, 0.0]] * Angstrom**-1. This option is mutually exclusive to qpoints_fractional.
  • electron_bands (list of int | All) – The band indices of the Bloch states to include.
    Default: All
  • phonon_modes (list of int | All) – Phonon modes to include.
    Default: All
  • calculate_velocities (bool) – Boolean to control if the velocities should be calculated or not.
    Default: True
  • maximum_interaction_range (PhysicalQuantity of type length) – Set the maximum range of the interactions.
    Default: No maximum interaction range.
  • rotate_to_pure_spin_states (bool) – Boolean to control if the eigenstates should be rotated to pure spin states or not. This is only relevant for non-collinear or spin-orbit calculations. As default for non-collinear calculations the eigenstates are rotated to pure spin states (True) and for spin-orbit no eigenstates are rotated (False).
  • energy_tolerance (PhysicalQuantity of type energy) – Energy tolerance for calculating matrix elements. The matrix elements are only calculated if \((E_{final} - E_{initial} \pm E_{phonon})\) is smaller in magnitude than the energy_tolerance, where \(E_{initial}\) is the initial state energy with wavevector k, \(E_{final}\) is the final state energy with wavevector k \(\pm\) q, and \(E_{phonon}\) is the phonon energy at phonon wavevector q.
    Default: 0.1 * eV
  • initial_state_energy_range (PhysicalQuantity of type energy) – Energy range for the initial state. Electron-phonon coupling elements are only calculated if the energy of the initial state is within the initial state energy range.
    Default: [-100, 100] * eV (all matrix elements are calculated included).
  • store_dense_coupling_matrices (bool) – Boolean to control if the electron-phonon coupling elements should be stored as dense matrices. This should only be set to True if a large energy_tolerance ( > 0.1*eV) is used.
    Default: False
  • coupling_matrix_tolerance (PhysicalQuantity of type energy) – Tolerance for storing the coupling matrix elements. The matrix element is only stored, if it is larger in magnitude than the coupling_matrix_tolerance. This parameter is only used, if the parameter store_dense_coupling_matrices is False.
    Default: 1e-9 * eV
couplingMatrixTolerance()

Query method for the coupling_matrix_tolerance.

Returns:The coupling_matrix_tolerance used in the calculation.
Return type:PhysicalQuantity of type energy.
couplingMatrixVsQ(spin_index=None, phonon_mode=None, k_index=None, initial_band=None, final_band=None)

Query method for the scaled coupling matrix vs. q at specified spin, phonon mode, k-index, initial band and final band.

Parameters:
  • spin_index (int) – The spin index.
  • phonon_mode (int) – The phonon mode.
  • k_index (int) – The k-point index.
  • initial_band (int) – The initial band.
  • final_band (int) – The final band.
Returns:

The electron-phonon coupling matrix vs. q. The shape of the array is (number_of_qpoints).

Return type:

PhysicalQuantity of type energy.

eigenvaluesK()

Query method for the electron Bloch state energies as a function of k.

Returns:The electron Bloch state energies as a function of k. The shape of the array is (number_of_spins, number_of_kpoints, number_of_electron_bands).
Return type:PhysicalQuantity of type energy.
eigenvaluesKMinusQ()

Query method for the electron Bloch state energies as a function of k-q.

Returns:The electron Bloch state energies as a function of k-q. The shape of the array is (number_of_spins, number_of_kpoints, number_of_qpoints, number_of_electron_bands).
Return type:PhysicalQuantity of type energy.
eigenvaluesKMinusQVsQ(spin_index=None, k_index=None, band=None)

Query method for the eigenvalues_k_minus_q vs. q at specified spin, k-index, band.

Parameters:
  • spin_index (int) – The spin index.
  • k_index (int) – The k-point index.
  • band (int) – The electron band.
Returns:

The electron Bloch state energies at k-q as a function of q. The shape for the returned array is (N, ), where N is the number of q-points.

Return type:

PhysicalQuantity of type energy.

eigenvaluesKPlusQ()

Query method for the electron Bloch state energies as a function of k+q.

Returns:The electron Bloch state energies as a function of k+q. The shape of the array is (number_of_spins, number_of_kpoints, number_of_qpoints, number_of_electron_bands).
Return type:PhysicalQuantity of type energy.
eigenvaluesKPlusQVsQ(spin_index=None, k_index=None, band=None)

Query method for the eigenvalues_k_plus_q vs. q at specified spin, k-index, band.

Parameters:
  • spin_index (int) – The spin index.
  • k_index (int) – The k-point index.
  • band (int) – The electron band.
Returns:

The electron Bloch state energies at k+q as a function of q. The shape for the returned array is (N, ), where N is the number of q-points.

Return type:

PhysicalQuantity of type energy.

electronBands()

Query method for the electron bands.

Returns:The list of electronic band indices used in the calculation.
Return type:list of int.
energyTolerance()

Query method for the energy_tolerance.

Returns:The energy_tolerance used in the calculation.
Return type:PhysicalQuantity of type energy.
evaluate(spin_index=None, phonon_mode=None, k_index=None, initial_band=None, final_band=None)

Return the electron-phonon coupling matrix vs. q at specified spin, phonon mode, k-index, initial band and final band.

Parameters:
  • spin_index (int) – The spin index.
  • phonon_mode (int) – The phonon mode.
  • k_index (int) – The k-point index.
  • initial_band (int) – The initial band.
  • final_band (int) – The final band.
Returns:

The electron-phonon coupling matrix vs. q. The shape of the array is (number_of_qpoints).

Return type:

PhysicalQuantity of type energy.

initialStateEnergyRange()

Query method for the initial state energy range.

Returns:The initial_state_energy_range used in the calculation. This is a list of two energies [Emin, Emax] defining the allowed energy range for the initial states.
Return type:PhysicalQuantity of type energy.
inverseCharacteristicLength()

Query method for inverse characteristic length.

Returns:The inverse characteristic length. The shape of the array is (number_of_degrees_of_freedom, number_of_qpoints)
Return type:PhysicalQuantity of type inverse length.
inverseLifeTime(phonon_modes=None, electron_bands=None, temperature=None, fermi_shift=None, integration_method=None, refinement=None)

Calculates the inverse life time due to electron-phonon coupling.

Parameters:
  • phonon_modes (List of indices or NLFlag All.) – Phonon modes to include.
    Default: All.
  • electron_bands (List of indices or NLFlag All.) – Electronic bands to include.
    Default: All (Include all the calculated electron bands).
  • temperature (PhysicalQuantity of type temperature) – The temperature for the thermal smearing of the Bose and Fermi distributions.
    Default: 300 * Kelvin
  • fermi_shift (PhysicalQuantity of energy unit.) – Shift of the Fermi level.
    Default: 0.0 * eV
  • integration_method (GaussianBroadening | TetrahedronMethod) – The method to use for calculating the q-integrale.
    Default: TetrahedronMethod
  • refinement (Positive int) – The number of times the q-grid is refined in each direction.
    Default: 1 (no interpolation)
Returns:

An averaged inverse life time, and phonon mode, k-point, and electron band resolved life time.

Return type:

Two PhysicalQuantity of type inverse time of shape (1), and (phonon modes, k-points, electron bands).

inverseSpinLifeTime(phonon_modes=None, electron_bands=None, temperature=None, fermi_shift=None, integration_method=None, refinement=None)

Calculates the inverse spin life time due to electron-phonon coupling. This can only be done for spin-orbit calculations.

Parameters:
  • phonon_modes (List of indices or NLFlag All.) – Phonon modes to include.
    Default: All.
  • electron_bands (List of indices or NLFlag All.) – Electronic bands to include.
    Default: All (Include all the calculated electron bands).
  • temperature (PhysicalQuantity of type temperature) – The temperature for the thermal smearing of the Bose and Fermi distributions.
    Default: 300 * Kelvin
  • fermi_shift (PhysicalQuantity of energy unit.) – Shift of the Fermi level.
    Default: 0.0 * eV
  • integration_method (GaussianBroadening | TetrahedronMethod) – The method to use for calculating the q-integrale.
    Default: TetrahedronMethod
  • refinement (Positive int) – The number of times the q-grid is refined in each direction.
    Default: 1 (no interpolation)
Returns:

An averaged inverse spin life time, and phonon mode, k-point, and electron band resolved spin life time.

Return type:

Two PhysicalQuantity of type inverse time of shape (1), and (phonon modes, k-points, electron bands).

kpoints()

Query method for the k-points in fractional coordinates.

Returns:The list of fractional k-points used in the calculation.
Return type:list
kpointsCartesian()

Query method for the k-points in Cartesian coordinates.

Returns:The list of Cartesian k-points used in the calculation.
Return type:PhysicalQuantity of type inverse length.
metatext()
Returns:The metatext of the object or None if no metatext is present.
Return type:str | unicode | None
nlprint(stream=None)

Print a string containing an ASCII table useful for plotting the AnalysisSpin object.

Parameters:stream (python stream) – The stream the table should be written to.
Default: NLPrintLogger()
phononEnergies()

Query method for the phonon energies.

Returns:The phonon energies as a function of q. The shape of the array is (number_of_degrees_of_freedom, number_of_qpoints), where number_of_degrees_of_freedom = 3*number_of_atoms.
Return type:PhysicalQuantity of type energy.
phononModes()

Query method for the phonon modes.

Returns:The list of phonon mode indices used in the calculation.
Return type:list of int.
qpoints()

Query method for the q-points in fractional coordinates.

Returns:The list of fractional q-points used in the calculation.
Return type:list
qpointsCartesian()

Query method for the q-points in Cartesian coordinates.

Returns:The list of Cartesian q-points used in the calculation.
Return type:PhysicalQuantity of type inverse length.
setMetatext(metatext)

Set a given metatext string on the object.

Parameters:metatext (str | unicode | None) – The metatext string that should be set. A value of “None” can be given to remove the current metatext.
storeDenseCouplingMatrices()

Query method for the store_dense_coupling_matrices Boolean.

Returns:The Boolean determining if the electron-phonon coupling matrix is stored in a dense format (True) or in a sparse format (False).
Return type:bool.
unscaledCouplingMatrixVsQ(spin_index=None, phonon_mode=None, k_index=None, initial_band=None, final_band=None)

Query method for the unscaled coupling matrix vs. q at specified spin, phonon mode, k-index, initial band and final band.

Parameters:
  • spin_index (int) – The spin index.
  • phonon_mode (int) – The phonon mode.
  • k_index (int) – The k-point index.
  • initial_band (int) – The initial band.
  • final_band (int) – The final band.
Returns:

The unscaled electron-phonon coupling matrix vs. q. The shape of the array is (number_of_qpoints).

Return type:

PhysicalQuantity of type energy / length.

velocitiesK()

Query method for the velocities as function of k.

Returns:The electron velocities as a function of k. The shape of the array is (number_of_spins, number_of_kpoints, number_of_electron_bands, 3).
Return type:PhysicalQuantity of type Meter/Second.
velocitiesKMinusQ()

Query method for the velocities as function of k-q.

Returns:The electron velocities as a function of k-q. The shape of the array is (number_of_spins, number_of_kpoints, number_of_qpoints, number_of_electron_bands, 3).
Return type:PhysicalQuantity of type Meter/Second
velocitiesKMinusQVsQ(spin_index=None, k_index=None, band=None)

Query method for the velocities_k_minus_q vs. q at specified spin, k-index, band.

Parameters:
  • spin_index (int) – The spin index.
  • k_index (int) – The k-point index.
  • band (int) – The electron band.
Returns:

The electron Bloch state velocities at k-q as a function of q. The shape for the returned array is (N, 3), where N is the number of q-points.

Return type:

PhysicalQuantity of type velocity.

velocitiesKPlusQ()

Query method for the velocities as function of k+q.

Returns:The electron velocities as a function of k+q. The shape of the array is (number_of_spins, number_of_kpoints, number_of_qpoints, number_of_electron_bands, 3).
Return type:PhysicalQuantity of type Meter/Second
velocitiesKPlusQVsQ(spin_index=None, k_index=None, band=None)

Query method for the velocities_k_plus_q vs. q at specified spin, k-index, band.

Parameters:
  • spin_index (int) – The spin index.
  • k_index (int) – The k-point index.
  • band (int) – The electron band.
Returns:

The electron Bloch state velocities at k+q as a function of q. The shape for the returned array is (N, 3), where N is the number of q-points.

Return type:

PhysicalQuantity of type velocity.

Usage Examples

Calculate the electron-phonon coupling in graphene for k-point around the K-point, which is \((1/3, 1/3, 0)\) in fractional coordinates, and for phonon q-point around the \(\Gamma\) point.

# -------------------------------------------------------------
# Bulk configuration
# -------------------------------------------------------------
lattice = Hexagonal(2.4612*Angstrom, 6.709*Angstrom)

elements = [Carbon, Carbon]

fractional_coordinates = [[0.0           , 0.0           , 0.5],
                          [0.333333333333, 0.666666666667, 0.5]]

bulk_configuration = BulkConfiguration(
    bravais_lattice=lattice,
    elements=elements,
    fractional_coordinates=fractional_coordinates
    )

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
# Use only Gamma-point for HamiltonianDerivatives and DynamicalMatrix
numerical_accuracy_parameters = NumericalAccuracyParameters(
    k_point_sampling=(1, 1, 1),
    )

calculator = LCAOCalculator(
    numerical_accuracy_parameters=numerical_accuracy_parameters,
    )

bulk_configuration.setCalculator(calculator)
bulk_configuration.update()

# -------------------------------------------------------------
# Hamiltonian derivatives
# -------------------------------------------------------------
hamiltonian_derivatives = HamiltonianDerivatives(
    bulk_configuration,
    repetitions=(5, 5, 1),
    )

# -------------------------------------------------------------
# Dynamical matrix
# -------------------------------------------------------------
dynamical_matrix = DynamicalMatrix(
    bulk_configuration,
    repetitions=(5, 5, 1),
    max_interaction_range=10*Angstrom,
    )

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
# Use k-point sampling for the ElectronPhononCoupling object.
numerical_accuracy_parameters = NumericalAccuracyParameters(
    k_point_sampling=(5, 5, 1),
    )

calculator = LCAOCalculator(
    numerical_accuracy_parameters=numerical_accuracy_parameters,
    )

bulk_configuration.setCalculator(calculator)
bulk_configuration.update()

# -------------------------------------------------------------
# Electron Phonon Coupling
# -------------------------------------------------------------
k_a = numpy.linspace(0.32333, 0.34333, 20)
k_b = numpy.linspace(0.32333, 0.34333, 20)
k_c = numpy.linspace(0, 0, 1)
kpoints = [[a, b, c] for a in k_a for b in k_b for c in k_c]

q_a = numpy.linspace(-0.02, 0.02, 20)
q_b = numpy.linspace(-0.02, 0.02, 20)
q_c = numpy.linspace(0, 0, 1)
qpoints = [[a, b, c] for a in q_a for b in q_b for c in q_c]

electron_phonon_coupling = ElectronPhononCoupling(
    configuration=bulk_configuration,
    dynamical_matrix=dynamical_matrix,
    hamiltonian_derivatives=hamiltonian_derivatives,
    kpoints_fractional=kpoints,
    qpoints_fractional=qpoints,
    electron_bands=All,
    phonon_modes=All,
    energy_tolerance=0.01*eV,
    initial_state_energy_range=[-100,100]*eV,
    store_dense_coupling_matrices=False,
    )

nlsave('graphene.nc', electron_phonon_coupling)


By specifying the input parameter ‘energy_tolerance=0.01*eV’ as above, you will only be calculating a matrix element, if the initial- and final energies are separated by the phonon energy, within the specified energy tolerance. When using a small ‘energy_tolerance’, it will also be advantagerous to use a sparse data structure. This is achieved by setting ‘store_dense_coupling_matrices=False’, as above. These are the default options. If you will only use the ElectronPhononCoupling object as input to a Mobility calculation, it is advantageous and faster to only calculate the matric elements for those (k,q)-pairs which are relevant for the mobility calculation. In the equations used for the mobility, the electron-phonon coupling matrix element always occurs together with delta-functions, like \(\delta(E_i(k) - E_f(k+q) + \hbar\omega_\lambda(q))\) or \(\delta(E_i(k) - E_f(k+q) - \hbar\omega_\lambda(q))\), where \(E_i(k)\) and \(E_f(k+q)\) are the energies of the initial- and final Bloch state, respectively, and \(\hbar\omega_\lambda(q)\) is the phonon energy.

If you which to visually ananlyze the coupling matrix elements with the ElectronPhononCoupling analyzer, it is advantagerous to calculate the matrix elements for all (k,q)-points. This is achieved by specifying a large energy tolerance, ‘energy_tolerance=100*eV’, like in the script below. In that case it will also be advantagerous to use a dense matrix representation. This can be specified by setting ‘store_dense_coupling_matrices=True’.

electron_phonon_coupling = ElectronPhononCoupling(
    configuration=bulk_configuration,
    dynamical_matrix=dynamical_matrix,
    hamiltonian_derivatives=hamiltonian_derivatives,
    kpoints_fractional=kpoints,
    qpoints_fractional=qpoints,
    electron_bands=All,
    phonon_modes=All,
    energy_tolerance=100*eV,
    initial_state_energy_range=[-100,100]*eV,
    store_dense_coupling_matrices=True,
    )

nlsave('graphene.nc', electron_phonon_coupling)

Often it will only be relevant to consider initial states close to the Fermi energy. For materials with complicated Fermi surfaces, it is difficult to provide a list of k-point spanning only the Fermi surface. It is therefor possible to specify ‘initial_state_energy_range=[-0.1,0.1]*eV’, relative to the Fermi energy (see below). Only matrix elements, with inital states with energies \(E_i(k)\) in the specified range will be calculated.

electron_phonon_coupling = ElectronPhononCoupling(
    configuration=bulk_configuration,
    dynamical_matrix=dynamical_matrix,
    hamiltonian_derivatives=hamiltonian_derivatives,
    kpoints_fractional=kpoints,
    qpoints_fractional=qpoints,
    electron_bands=All,
    phonon_modes=All,
    energy_tolerance=0.01*eV,
    initial_state_energy_range=[-0.1,0.1]*eV,
    store_dense_coupling_matrices=False,
    )

nlsave('graphene.nc', electron_phonon_coupling)


Theoretical notes

The theoretical details of the electron-phonon coupling in QuantumATK is described in [GMSB16].

In order to calculate the electron-phonon coupling three main ingredients are needed:

  • An updated configuration in order to get the Hamiltonian of the system and the electronic Bloch states.
  • A DynamicalMatrix object for calculating the phonon modes.
  • A HamiltonianDerivatives object for calculating the electron-phonon coupling matrix.

The configuration used as input to the ElectronPhononCoupling constructor must be exactly the same as used for the HamiltonianDerivatives and DynamicalMatrix objects.

From the dynamical matrix \(D\) we can obtain the phonon frequencies \(\omega_{\lambda,\mathbf{q}}\) and eigenvectors \(\mathbf{u}_{\lambda,\mathbf{q}}\) for a particular phonon band index \(\lambda\) and phonon wave number \(\mathbf{q}\) :

\[D(\mathbf{q}) \mathbf{u}_{\lambda,\mathbf{q}} = \omega_{\lambda,\mathbf{q}}^2 \mathbf{u}_{\lambda,\mathbf{q}},\]

where the phonon mode vectors are normalized, i.e. \(\mathbf{u}^T\cdot\mathbf{u}=1\). We further define mass-scaled phonon mode vectors according to

\[\mathbf{v}_{\lambda,\mathbf{q}}(\alpha) = \mathbf{u}_{\lambda,\mathbf{q}}(\alpha) \sqrt{ \frac{\hbar}{2M_\alpha \omega_{\lambda, \mathbf{q}}}},\]

where \(M_\alpha\) is the mass of the atom with index \(\alpha\). Let now \(j=(\alpha,\mu)\) denote a combined index for atomic site ( \(\alpha\) ) and orbital index ( \(\mu\) ), and let an electronic Bloch state be expanded in the LCAO basis as:

\[|\mathbf{k} \rangle = \sum_j c_j^\mathbf{k} | j \mathbf{k} \rangle,\]

where

\[|j\mathbf{k} \rangle = \frac{1}{N}\sum_n e^{2\pi i \mathbf{k}\cdot\mathbf{R}_n} | j \mathbf{R}_n \rangle,\]

and \(| j \mathbf{R}_n \rangle\) is an LCAO basis orbital in the n’th cell displaced by lattice vector \(\mathbf{R}_n\) from the central cell with index ‘0’. The electron-phonon coupling matrix for a given phonon mode, \(\{\mathbf{q},\lambda\}\) is then obtained as [KTJ12], [GMSB16]:

\[g_{ij}^\lambda(\mathbf{k},\mathbf{q}) = \sum_{ij} (c_i^{\mathbf{k}+\mathbf{q}})^*c_j^\mathbf{k} \sum_{mn}e^{2\pi (i\mathbf{k}\cdot(\mathbf{R}_n-\mathbf{R}_m) - i\mathbf{q}\cdot\mathbf{R}_m)}\langle i \mathbf{R}_m | \mathbf{v}_{\mathbf{q},\lambda}\cdot\nabla H_0(r)|j\mathbf{R}_n\rangle,\]

where the (ij)-sum runs over combined atom site and orbital indices, and the (mn)-sum runs over repeated unit cells in the HamiltonianDerivatives calculation. Specifically the derivative of the Hamiltonian \(\nabla H_0(r)\) is being calculated in the HamiltonianDerivatives object, where the subscript 0 indicates that the derivatives are only calculated for atoms in the unit cell with index ‘0’.

Spin life time calculations

The spin life time is an important parameter for spintronics devices. At technologically relevant temperatures (> 100 K) the spin life time will, in many materials, be limited by electron-phonon interactions. These are mediated by spin-orbit coupling through the so-called Elliot-Yafet mechanism. The implementation in QuantumATK follows the theory in Ref. [RW12]. The phonon-limited spin life time can be calculated from an ElectronPhononCoupling object, provided that the calculator attached to the configuration, and used for the HamiltonianDerivatives object, has spin set to ‘Noncollinear Spin Orbit’. The inverse spin life time is calculated as:

# Calculate inverse spin life time
average_inv_time, all_inv_times = electron_phonon_coupling.inverseSpinLifeTime(
    phonon_modes=All,
    electron_bands=All,
    temperature=300*K,
    fermi_shift=0.0*eV,
    integration_method=TetrahedronMethod)

# Spin life time
T_1 = 1/average_inv_time

The function returns two physical quantities of unit inverse time. The first is a band-averaged inverse relaxation time given by

\[\langle \frac{1}{T_1} \rangle = \frac{\sum_n \int_0^{k_F} [1/T_1(n, k) \frac{\partial f}{\partial \epsilon}|_{\epsilon=\epsilon_F}]d^3k}{\sum_n \int_0^{k_F} [\frac{\partial f}{\partial \epsilon}|_{\epsilon=\epsilon_F}]d^3k}\]

where \(f=f(\epsilon_k)\) is the Fermi-Dirac distribution function, and the band- and \(k\)-resolved inverse spin life times are calculated as

\[\begin{split}\frac{1}{T_1(n,k)} &= \frac{2\pi}{\hbar}\sum_{q\lambda n'}|g_{k+qn'\uparrow, kn\downarrow}^{q \lambda}|^2 \{[f(\epsilon_{k+q,n'}) + n_{q\lambda} ]\delta(\epsilon_{k,n} - \epsilon_{k+q,n'} + \hbar\omega_{q\lambda}) \\ &\quad + [1 + n_{q\lambda} - f(\epsilon_{k+q,n'})]\delta(\epsilon_{k,n} - \epsilon_{k-q,n'} - \hbar\omega_{q\lambda}) \}\end{split}\]

where \(\hbar\omega_{q\lambda}\) are the phonon energies for phonon mode \(\lambda\) and momentum \(q\) with Bose-Einstein occupation factors \(n_{q\lambda}\). Notice that only matrix elements between spin up states (\(\uparrow\)) and spin-down states (\(\downarrow\)) are included.

The sum (integral) over \(q\) can be performed in two different ways, determined by the input variable integration_method which can either be TetrahedronMethod for a tetrahedron integration or GaussianBroadening(0.01*eV) for a gaussian integration, where the \(\delta\)-functions are replaced by gaussians with a certain broadening (in this case 0.01 eV). For three-dimensional sampling the TetrahedronMethod converges significantly faster, thus reducing the required number of q-points. For one-dimensional sampling, however, the GaussianBroadening is recomended.

[GMSB16](1, 2) T. Gunst, T. Markussen, K. Stokbro, and M. Brandbyge. First-principles method for electron-phonon coupling and electron mobility: Applications to two-dimensional materials. Phys. Rev. B, 93:035414, Jan 2016. doi:10.1103/PhysRevB.93.035414.
[KTJ12]K. Kaasbjerg, K. S. Thygesen, and K. W. Jacobsen. Phonon-limited mobility in n-type single-layer mos2 from first principles. Phys. Rev. B, 85:115317, Mar 2012. doi:10.1103/PhysRevB.85.115317.
[RW12]O. D. Restrepo and W. Windl. Full first-principles theory of spin relaxation in group-iv materials. Phys. Rev. Lett., 109:166604, 2012. doi:10.1103/PhysRevLett.109.166604.