Documentation

The complete PyPL reference is reported below.

Huang-Rhys solver

class pypl.hr_solver.hr_solver[source]

Bases: object

High-level solver for computing Huang-Rhys factors (HRFs), spectral density, and optical lineshapes from phonon and structural data.

compute_hrf_forces(phonon_freqs_in, phonon_modes_in, atomic_symbols, forces_in, mass_list=None)[source]

Compute Huang-Rhys factors using atomic forces.

Parameters:
  • phonon_freqs_in (ndarray of shape (M,)) – Phonon frequencies in THz.

  • phonon_modes_in (ndarray of shape (M, 3N)) – Phonon eigenmodes in Cartesian representation.

  • atomic_symbols (list of str) – Atomic symbols for all atoms.

  • forces_in (ndarray of shape (N, 3)) – Atomic forces in eV/Å.

  • mass_list (dict, optional) – Dictionary mapping element symbols to atomic masses (in a.u.).

Returns:

Dictionary with keys: - 'freqs' : ndarray of shape (M,), phonon frequencies in rad/s. - 'hr_factors' : ndarray of shape (M,), Huang–Rhys factors \(S_k\).

Return type:

dict

compute_hrf_dis(phonon_freqs_in, phonon_modes_in, atomic_symbols, gs_coord_in, es_coord_in, cell_parameters_in, mass_list=None)[source]

Compute Huang-Rhys factors using atomic displacements between ground-state (GS) and excited-state (ES) geometries.

Parameters:
  • phonon_freqs_in (ndarray of shape (M,)) – Phonon frequencies in THz.

  • phonon_modes_in (ndarray of shape (M, 3N)) – Phonon eigenmodes in Cartesian representation.

  • atomic_symbols (list of str) – Atomic symbols for all atoms.

  • gs_coord_in (ndarray of shape (N, 3)) – Ground-state atomic coordinates in Å.

  • es_coord_in (ndarray of shape (N, 3)) – Excited-state atomic coordinates in Å.

  • cell_parameters_in (ndarray of shape (3, 3)) – Lattice vectors in Å.

  • mass_list (dict, optional) – Dictionary mapping element symbols to atomic masses (in a.u.).

Returns:

Dictionary with keys: - 'freqs' : ndarray of shape (M,), phonon frequencies in rad/s. - 'hr_factors' : ndarray of shape (M,), Huang–Rhys factors \(S_k\).

Return type:

dict

static gaussian(x, mu, sigma)[source]

Compute the value of a Gaussian function.

\[G(x) = \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left(-\frac{(x - \mu)^2}{2 \sigma^2}\right)\]
Parameters:
  • x (float) – Evaluation point.

  • mu (float) – Mean of the Gaussian.

  • sigma (float) – Standard deviation.

Returns:

Value of the Gaussian at x.

Return type:

float

static f_sigma(freqs, sigma)[source]

Compute a frequency-dependent linewidth function by linear interpolation.

\[\sigma(\omega_k) = \sigma_0 - \frac{\sigma_0 - \sigma_1}{\max(\omega_k) - \min(\omega_k)} (\omega_k - \min(\omega_k))\]
Parameters:
  • freqs (ndarray of shape (M,)) – Phonon frequencies in rad/s.

  • sigma (list of float) – Two-element list [sigma_0, sigma_1] in meV, giving the linewidth at the minimum and maximum frequency.

Returns:

collect_sigma – Interpolated linewidths \(\sigma(\omega_k)\) in meV.

Return type:

ndarray of shape (M,)

compute_spectral_density(hrf_dict, energy_axis=None, sigma=[6.0, 1.5])[source]

Compute the phonon spectral density with Gaussian broadening.

\[S(\hbar\omega) = \sum_k S_k \, G(\hbar\omega; \hbar\omega_k, \sigma(\omega_k))\]

where \(S_k\) is the partial Huang–Rhys factor and \(G\) is a normalized Gaussian.

Parameters:
  • hrf_dict (dict) – Dictionary containing: - 'freqs' : ndarray of phonon frequencies (rad/s). - 'hr_factors' : ndarray of partial Huang–Rhys factors.

  • energy_axis (ndarray of shape (N,), optional) – Energy axis in meV over which to evaluate the spectral density. If None, defaults to np.linspace(0, 200, 2001).

  • sigma (list of float, optional) – Linewidth interpolation values [sigma_0, sigma_1] in meV.

Returns:

  • energy_axis (ndarray of shape (N,)) – Energy axis in meV.

  • spectral_density (ndarray of shape (N,)) – Spectral density \(S(\hbar\omega)\) in 1/meV.

compute_lineshape_numerical_integration(hrf_dict, temp=4, sigma=[6, 1.5], zpl_broadening=0.3, time_range=[0, 20000], time_resolution=200001, lineshape_energy_range=[-150, 550], lineshape_energy_resolution=701)[source]

Compute the temperature-dependent optical lineshape via direct numerical integration in the time domain.

\[A(\hbar\omega, T) = \int dt \, \exp\!\left(i \omega t - \frac{|t|\gamma_\mathrm{ZPL}}{\hbar}\right) G(t, T)\]

where \(G(t, T)\) is the generating function with phonon broadening, and \(\gamma_\mathrm{ZPL}\) is the zero-phonon line (ZPL) broadening.

Parameters:
  • hrf_dict (dict) – Dictionary containing phonon frequencies and HR factors.

  • temp (float, optional) – Temperature in Kelvin. Default is 4.

  • sigma (list of float, optional) – Linewidth interpolation values [sigma_min, sigma_max] in meV.

  • zpl_broadening (float, optional) – Additional Lorentzian broadening of the ZPL in meV. Default is 0.3.

  • time_range (list of float, optional) – Time domain [start, end] in fs. Default is [0, 20000].

  • time_resolution (int, optional) – Number of time points. Default is 200001.

  • lineshape_energy_range (list of float, optional) – Energy range [start, end] in meV. Default is [-150, 550].

  • lineshape_energy_resolution (int, optional) – Number of energy points. Default is 701.

Returns:

(energy_axis, A(E)), where: - energy_axis : ndarray of energies in meV. - A(E) : ndarray of the computed lineshape.

Return type:

tuple

compute_lineshape_fft(hrf_dict, temp=4, sigma=[6, 1.5], zpl_broadening=0.3, lineshape_energy_range=[-1000, 1000], lineshape_energy_resolution=2001)[source]

Compute the temperature-dependent optical lineshape using FFT of the time-domain correlation function.

\[A(\hbar\omega, T) = \frac{1}{2\pi\hbar} \int dt \, \exp\!\left(i \omega t - \frac{|t|\gamma_\mathrm{ZPL}}{\hbar}\right) G(t, T)\]
Parameters:
  • hrf_dict (dict) – Dictionary containing phonon frequencies and HR factors.

  • temp (float, optional) – Temperature in Kelvin. Default is 4.

  • sigma (list of float, optional) – Linewidth interpolation values [sigma_min, sigma_max] in meV.

  • zpl_broadening (float, optional) – Lorentzian broadening of the ZPL in meV. Default is 0.3.

  • lineshape_energy_range (list of float, optional) – Symmetric energy range [start, end] in meV (must satisfy start = -end). Default is [-1000, 1000].

  • lineshape_energy_resolution (int, optional) – Number of energy points. Default is 2001.

Returns:

(energy_axis, A(E)), where: - energy_axis : ndarray of energies in meV. - A(E) : ndarray of the computed lineshape.

Return type:

tuple

Raises:

ValueError – If lineshape_energy_range is not symmetric about zero.

compute_spectrum(ezpl=1.0, spectrum_type='PL', lineshape=None)[source]

Compute the normalized photoluminescence (PL) or absorption spectrum from a precomputed lineshape.

Parameters:
  • ezpl (float, optional) – Zero-phonon line (ZPL) energy in meV. Default is 1.0.

  • spectrum_type ({'PL', 'Abs'}, optional) – Type of spectrum to compute: - 'PL' : photoluminescence spectrum - 'Abs' : absorption spectrum

  • lineshape (tuple, optional) – Precomputed lineshape (energy_axis, A(E)). If None, uses self.lineshape.

Returns:

  • energy_axis_out (ndarray) – Energy axis in meV, shifted relative to the ZPL.

  • spectrum (ndarray) – Normalized spectral intensity.

Lineshape

class pypl.lineshape.lineshape(hrf)[source]

Bases: object

Class to compute optical lineshapes from Huang-Rhys factors.

Parameters:

hrf (dictionary) – hrf[‘freqs’] for phonon frequency. hrf[‘hr_factors’] for partial Huang-Rhys factors.

f_sigma()[source]

Compute phonon-dependent broadening parameters.

The broadening is interpolated linearly between two given values across the range of phonon frequencies:

\[\sigma(\omega_k) = \sigma_0 - \frac{\sigma_0 - \sigma_1}{\max(\omega_k) - \min(\omega_k)} \, (\omega_k - \min(\omega_k))\]

Notes

The computed broadening parameters are stored in the attribute self.all_sigma (ndarray), which contains one value per phonon frequency.

compute_lineshape_numerical_integration(temp, sigma, zpl_broadening, time_axis, ene_axis)[source]

Compute the lineshape function by direct numerical integration.

The lineshape is defined as

\[A(\hbar\omega, T) = \int dt \, G(t, T) \, e^{i \omega t},\]

where the generating function \(G(t, T)\) is constructed from Huang-Rhys factors, phonon occupations, and broadening terms.

Parameters:
  • temp (float) – Temperature in Kelvin.

  • sigma (list of float) – Two-element list [sigma0, sigma1] giving phonon broadening parameters in meV.

  • zpl_broadening (float) – Zero-phonon line Lorentzian broadening (HWHM) in meV.

  • time_axis (ndarray) – Time axis in femtoseconds.

  • ene_axis (ndarray) – Energy axis in meV.

Notes

The computed lineshape is stored in the attribute self.lineshape, a tuple (energy_axis, A(E)) where energy_axis is in meV and A(E) is the computed lineshape.

compute_lineshape_fft(temp, sigma, zpl_broadening, energy_axis)[source]

Compute the lineshape function by FFT of the time-domain correlation.

The time-domain generating function is

\[G(t) = \exp \Big[ S(t) - S(0) + 2C(t) - 2C(0) \Big],\]

where

\[\begin{split}S(t) &= \sum_k S_k \exp\!\left(-\frac{\sigma_k^2t^2}{2\hbar^2}\right) e^{i \omega_k t}, \\\\ C(t) &= \sum_k n_k S_k \exp\!\left(-\frac{\sigma_k^2t^2}{2\hbar^2}\right) \cos(\omega_k t).\end{split}\]

The frequency-domain lineshape is then obtained as

\[A(\hbar\omega) = \frac{1}{2\pi\hbar} \!\int dt \exp \left(i \omega t - \frac{|t| \gamma_\mathrm{ZPL}}{\hbar} \right) G(t).\]
Parameters:
  • temp (float) – Temperature in Kelvin.

  • sigma (list of float) – Two-element list [sigma0, sigma1] for phonon broadening parameters (in meV).

  • zpl_broadening (float) – Zero-phonon line Lorentzian broadening (HWHM) in meV.

  • energy_axis (ndarray) – Target energy axis in meV.

Notes

The computed lineshape is stored in the attribute self.lineshape, a tuple (energy_axis, A(E)) where energy_axis is in meV and A(E) is the computed lineshape.

Huang-Rhys factors

class pypl.hr_factors.hr_factors(freqs, modes, atomic_symbols, mass_list=None)[source]

Bases: object

Compute Huang-Rhys factors and spectral density from phonon frequencies, eigenmodes, and either forces or displacements. All inputs, except atomic masses, should be in SI units.

Parameters:
  • freqs (ndarray of shape (3N,)) – Phonon frequencies in rad/s.

  • modes (ndarray of shape (3N, N, 3)) – Phonon eigenmodes.

  • atomic_symbols (list of str) – Atomic chemical symbols.

  • mass_list (dict, optional) – Dictionary mapping element symbols to atomic masses in unified atomic mass units (u). If not provided, IUPAC-recommended values are used.

set_masses(mass_list=None)[source]

Set atomic masses in kilograms for all atoms in the system.

Parameters:

mass_list (dict, optional) – Dictionary mapping element symbols to atomic masses in unified atomic mass units (u). If not provided, IUPAC-recommended values are used.

Notes

The converted atomic masses are stored in self.masses in kilograms.

compute_hrf_forces(forces)[source]

Compute Huang-Rhys factors using atomic forces.

\[S_k = \frac{\omega_k \Delta Q_k^2}{2 \hbar}\]

where

\[\Delta Q_k = \frac{1}{\omega_k^2} \sum_{\alpha=1}^N \sum_{i=x,y,z} \frac{F_{\alpha i}}{\sqrt{M_\alpha}} e_{k, \alpha i}\]
Parameters:

forces (ndarray of shape (N, 3)) – Forces on atoms in J/m (SI units).

Notes

The computed Huang-Rhys factors are stored in self.hrf.

compute_hrf_dis(gs_coord, es_coord, cell_parameters)[source]

Compute Huang-Rhys factors using atomic displacements.

\[S_k = \frac{\omega_k \Delta Q_k^2}{2 \hbar}\]

with

\[\Delta Q_k = \sum_{\alpha=1}^N \sum_{i=x,y,z} \sqrt{M_\alpha} \Delta R_{\alpha i} e_{k, \alpha i}\]
Parameters:
  • gs_coord (ndarray of shape (N, 3)) – Ground-state atomic coordinates (in meters).

  • es_coord (ndarray of shape (N, 3)) – Excited-state atomic coordinates (in meters).

  • cell_parameters (ndarray of shape (3, 3)) – Lattice vectors of the simulation cell (in meters).

Notes

Displacements are wrapped back into the unit cell. The computed Huang-Rhys factors are stored in self.hrf.

1D configurational coordinate solver

class pypl.config_coord_1d_solver.config_coord_1d_solver(omega_i, omega_f, delta_q, prec=100)[source]

Bases: object

One-dimensional configuration coordinate diagram solver.

Parameters:
  • omega_i (float) – Initial vibrational frequency (in meV).

  • omega_f (float) – Final vibrational frequency (in meV).

  • delta_q (float) – Displacement along the configuration coordinate (in Å·amu^{1/2}).

  • prec (int, optional) – Precision for decimal operations, default is 100.

compute_franck_condon_integrals(ni, nf)[source]

Compute Franck-Condon (FC) overlap integrals between vibrational states in displaced harmonic oscillators.

The integrals are defined as

\[F_{ij} = \langle \chi_i(Q + \Delta Q) \mid \chi_j(Q) \rangle = \langle \chi_i(Q) \mid \chi_j(Q - \Delta Q) \rangle\]

with recurrence relations given in [Peder Thusgaard Ruhoff, Chem. Phys. 186, 355-374 (1994)].

Parameters:
  • ni (int) – Number of vibrational states in the initial potential.

  • nf (int) – Number of vibrational states in the final potential.

Notes

self.fc_ints, a ndarray of shape (ni, nf), stores the Franck-Condon overlap integrals.

static Gaussian(x, mu, sigma_r, sigma_factor=160)[source]

Gaussian line shape function.

\[G(x; \mu, \sigma) = \frac{1}{\sqrt{2\pi\sigma^2}} \exp\!\left( - \frac{(x - \mu)^2}{2\sigma^2} \right)\]

where

\[\sigma = ( \sigma_r[1] - \sigma_r[0] ) \cdot \frac{|\mu|}{\text{sigma\_factor}} + \sigma_r[0]\]
Parameters:
  • x (ndarray) – Energy axis.

  • mu (float or ndarray) – Center position.

  • sigma_r (tuple of float) – Range for Gaussian broadening (min, max).

  • sigma_factor (float, optional) – Scaling factor for energy dependence of sigma, default is 160.

Returns:

Gaussian line shape values.

Return type:

ndarray

static Lorentzian(x, mu, gamma)[source]

Lorentzian line shape function.

\[L(x; \mu, \gamma) = \frac{1}{\pi \gamma} \frac{\gamma^2}{(x - \mu)^2 + \gamma^2}\]
Parameters:
  • x (ndarray) – Energy axis.

  • mu (float) – Center position.

  • gamma (float) – Broadening parameter.

Returns:

Lorentzian line shape values.

Return type:

ndarray

static BZ(ene, T)[source]

Boltzmann factor.

\[BZ(E, T) = \exp\!\left( -\frac{E}{k_B T} \right)\]
Parameters:
  • ene (float or ndarray) – Vibrational energy in meV.

  • T (float) – Temperature in Kelvin.

Returns:

Boltzmann factor.

Return type:

float or ndarray

bulid_fc_lsp(eneaxis=None, temp=4, sigma=10, zpl_lorentzian=True, gamma=1)[source]

Build the Franck-Condon (FC) lineshape function.

The spectrum is given by

\[I(E) = \sum_{i, j} p_i \, |F_{ij}|^2 \, G(E; E_{ij}, \sigma)\]

where \(p_i\) is the Boltzmann prefactor, \(F_{ij}\) are Franck-Condon integrals, and \(G\) is a Gaussian (or Lorentzian for the ZPL).

Parameters:
  • eneaxis (ndarray, optional) – Energy axis in meV. Default is np.linspace(-1000, 1000, 2001).

  • temp (float, optional) – Temperature in Kelvin. Default is 4 K.

  • sigma (float, optional) – Gaussian broadening width. Default is 10 meV.

  • zpl_lorentzian (bool, optional) – Replace zero-phonon line (ZPL) Gaussian with Lorentzian. Default is True.

  • gamma (float, optional) – Lorentzian broadening parameter for ZPL. Default is 1 meV.

Returns:

  • eneaxis (ndarray) – Energy axis (meV).

  • fc_lineshape (ndarray) – Computed FC lineshape intensity.

compute_spectrum(eneaxis=None, linshape=None, tdm=1.0, zpl=0.0, spectrum_type='PL')[source]

Compute the optical spectrum from the FC lineshape.

For photoluminescence (PL):

\[I_\text{PL}(E) \propto \text{linshape}(E) \cdot \mu^2 \cdot (E_\text{ZPL} - E)^3\]

For absorption (ABS):

\[I_\text{ABS}(E) \propto \text{linshape}(E) \cdot \mu^2 \cdot (E_\text{ZPL} + E)\]
Parameters:
  • eneaxis (ndarray, optional) – Energy axis (meV). Default is self.eneaxis.

  • linshape (ndarray, optional) – Franck-Condon lineshape. Default is self.lineshape.

  • tdm (float, optional) – Transition dipole moment. Default is 1.0.

  • zpl (float, optional) – Zero-phonon line (ZPL) energy. Default is 0.0 meV.

  • spectrum_type ({'PL', 'Abs'}, optional) – Type of spectrum to compute: ‘PL’ (photoluminescence) or ‘Abs’ (absorption). Default is ‘PL’.

Returns:

  • eneaxis_out (ndarray) – Shifted energy axis (meV).

  • spectrum (ndarray) – Spectrum intensity.

Utils

pypl.utils.parse_atoms_qexml(fileName)[source]

Parse atomic structure information from a Quantum ESPRESSO XML output file.

Parameters:

fileName (str) – Path to the QE XML file.

Returns:

  • atomic_symbols (list of str) – Atomic symbols for all atoms in the structure.

  • atomic_positions (ndarray of shape (N, 3)) – Cartesian atomic coordinates in Ångström.

  • cell_parameters (ndarray of shape (3, 3)) – Lattice vectors (cell parameters) in Ångström.

pypl.utils.parse_forces_qexml(fileName)[source]

Parse atomic forces from a Quantum ESPRESSO XML output file.

Parameters:

fileName (str) – Path to the QE XML file.

Returns:

  • atomic_symbols (list of str) – Atomic symbols for all atoms in the structure.

  • forces (ndarray of shape (N, 3)) – Forces on atoms in eV/Å.

Raises:

ValueError – If the <forces> tag cannot be found in the XML file.

pypl.utils.parse_total_energy_qexml(fileName)[source]

Parse the final total energy from a Quantum ESPRESSO XML output file.

Parameters:

fileName (str) – Path to the QE XML file.

Returns:

total_energy – Total energy in eV.

Return type:

float

Raises:

ValueError – If the <total_energy/etot> tag cannot be found in the XML file.

pypl.utils.parse_phonopy_h5(fileName, real_tol=1e-10)[source]

Parse phonon frequencies and eigenmodes from a Phonopy HDF5 file.

Parameters:
  • fileName (str) – Path to the Phonopy HDF5 file (e.g., phonon.hdf5).

  • real_tol (float, optional) – Tolerance below which imaginary components of eigenmodes are discarded. Default is 1e-10.

Returns:

  • freqs (ndarray of shape (M,)) – Phonon frequencies in THz.

  • modes (ndarray of shape (M, Nat, 3)) – Phonon eigenmodes. Converted to real if imaginary parts are negligible.

Raises:

ValueError – If the imaginary part of eigenmodes exceeds real_tol.

pypl.utils.parse_phonopy_yaml(fileName='band.yaml', real_tol=1e-10)[source]

Parse phonon frequencies and eigenmodes from a Phonopy band.yaml file.

Parameters:
  • fileName (str, optional) – Path to the band.yaml file. Default is "band.yaml".

  • real_tol (float, optional) – Tolerance for discarding imaginary parts of eigenmodes. Default is 1e-10.

Returns:

  • freqs (ndarray of shape (3N,)) – Phonon frequencies in THz.

  • modes (ndarray of shape (3N, N, 3)) – Phonon eigenmodes.

Raises:

ValueError – If the imaginary component of any eigenvector exceeds real_tol.