Electron Repulsion Integrals¶
- dopyqo.eri(c_ip: ndarray, b: ndarray, mill: ndarray, fft_grid: ndarray | None = None, use_gpu: bool = True) ndarray[source]¶
Calculate Electron Repulsion Integrals (ERIs) via pair densities in the physicists’ index order We calculate
\[h_{ijkl}=4\pi \sum_{p \neq 0} \rho_{il}(-p)\rho_{jk}(p)/|p|² =4\pi \sum_{p \neq 0} \rho*_{li}(p)\rho_{jk}(p)/|p|²\]where \(i,j,k,l\) index spatial orbitals. Since the momenta \(p\) and reciprocal-space wavefunctions \(\psi_i(p)\) (given as
c_ip) are given as a list, we need to transfer \(\psi_i(p)\) to a 3D grid. Then \(\psi_i(p)\) is represented on a 3D reciprocal-space grid. We calculate \(1/|\mathbf{p}|^2\) on all grid points, resulting in an infinite value at the center of the grid. This infinite value is set to zero; therefore, technically we later perform a sum over all momenta \(p\) except the zero momentum. Note that there are different techniques for handling this singularity, e.g. using a resolution-of-identity.Note that \(h_{ijkl}\) does not explicitly depend on the k-point but only implicitly depends on the k-point via the plane-wave coefficients
c_ip. Any explicit dependency is cancelled out since the k-point contributes a phase to each KS-orbital \(e^{i\mathbf{k}\cdot\mathbf{r}}\), but since the pair densities are a product of two orbitals where one is complex conjugated and the other is not, the phases cancel: \((e^{i\mathbf{k}\cdot\mathbf{r}})^* e^{i\mathbf{k}\cdot\mathbf{r}} = e^{-i\mathbf{k}\cdot\mathbf{r}} e^{i\mathbf{k}\cdot\mathbf{r}} = 1\).- Parameters:
c_ip (np.ndarray) – Array of coefficients describing the Kohn-Sham orbitals in the plane wave basis, shape (#bands, #waves)
b (np.ndarray) – Array of reciprocal lattice vectors, shape (3, 3)
mill (np.ndarray) – Array of Miller indices, shape (#waves, 3)
fft_grid_size (np.ndarray | None) – Edge lengths of the use FFT grid. Has to be larger or equal to grid where wavefunctions are defined on. Defaults to None, for which the wavefunction grid is used.
- Returns:
ERIs in reciprocal space
- Return type:
np.ndarray
- dopyqo.get_frozen_core_energy_pp(p: ndarray, c_ip_core: ndarray, b: ndarray, mill: ndarray, cell_volume: float, atom_positions: ndarray, atomic_numbers: ndarray, occupations_core: ndarray, pseudopots: list[Pseudopot], fft_grid: ndarray | None = None, use_gpu: bool = True) float[source]¶
Calculate frozen core energy ERI part
\[2\sum_i^{\text{frozen}} h_{ii} + \sum_{ij}^{\text{frozen}} (2h_{iijj} - h_{ijji})\]- Parameters:
p (np.ndarray) – Array of momentum vectors, shape (#waves, 3)
c_ip_core (np.ndarray) – Array of coefficients describing the core Kohn-Sham orbitals in the plane wave basis, shape (#bands, #waves)
b (np.ndarray) – Array of reciprocal lattice vectors, shape (3, 3)
mill (np.ndarray) – Array of Miller indices, shape (#waves, 3)
cell_volume (float) – Cell volume of the computational real space cell.
atom_positions (np.ndarray) – 1D array of positions of each atom
atomic_numbers (np.ndarray) – 1D array of atomic number of each atom
occupations_core (np.ndarray) – Array of occupations of the core orbitals.
- Returns:
Frozen core energy
- Return type:
- dopyqo.get_frozen_core_energy_given_pp(p: ndarray, c_ip_core: ndarray, b: ndarray, mill: ndarray, cell_volume: float, pp_core: ndarray, occupations_core: ndarray, fft_grid: ndarray | None = None, use_gpu: bool = True) float[source]¶
Calculate frozen core energy ERI part
\[2\sum_i^{\text{frozen}} h_{ii} + \sum_{ij}^{\text{frozen}} (2h_{iijj} - h_{ijji})\]- Parameters:
p (np.ndarray) – Array of momentum vectors, shape (#waves, 3)
c_ip_core (np.ndarray) – Array of coefficients describing the core Kohn-Sham orbitals in the plane wave basis, shape (#bands, #waves)
b (np.ndarray) – Array of reciprocal lattice vectors, shape (3, 3)
mill (np.ndarray) – Array of Miller indices, shape (#waves, 3)
cell_volume (float) – Cell volume of the computational real space cell.
pp_core (np.ndarray) – Pseudopotential matrix for the core orbitals
occupations_core (np.ndarray) – Array of occupations of the core orbitals.
- Returns:
Frozen core energy
- Return type:
- dopyqo.get_frozen_core_pot_and_energy_given_pp(c_ip_core: ndarray, c_ip_active: ndarray, b: ndarray, mill: ndarray, p: ndarray, cell_volume: float, pp_core: ndarray, occupations_core: ndarray, fft_grid: ndarray | None = None, use_gpu: bool = True) tuple[ndarray, float][source]¶
- dopyqo.get_frozen_core_pot(c_ip_core: ndarray, c_ip_active: ndarray, b: ndarray, mill: ndarray, fft_grid: ndarray | None = None, calc_eri_energy: bool = False, cell_volume: float | None = None, use_gpu: bool = True) ndarray | tuple[ndarray, float][source]¶
Calculate frozen core effective single-particle potential
\[V_{tu} = \sum_{i}^{\text{frozen}} (2h_{tuii} - h_{tiiu})\](chemists’ index order, notation from Saad Yalouz et al. 2021 Quantum Sci. Technol. 6 024004, Appendix B). Here t, u are active orbital indices and i are core orbital indices.
In physicists’ index order this becomes
\[V_{tu} = \sum_{i}^{\text{frozen}} (2h_{tiiu} - h_{tiui})\]- Parameters:
c_ip_core (np.ndarray) – Array of coefficients describing the core Kohn-Sham orbitals in the plane wave basis, shape (#bands, #waves)
c_ip_active (np.ndarray) – Array of coefficients describing the active Kohn-Sham orbitals in the plane wave basis, shape (#bands, #waves)
b (np.ndarray) – Array of reciprocal lattice vectors, shape (3, 3)
mill (np.ndarray) – Array of Miller indices, shape (#waves, 3)
calc_eri_energy (bool) – Defaults to False.
cell_volume (float | None) – Has to be given if calc_eri_energy is True. Defaults to None.
- Returns:
- ERIs in reciprocal space or
tuple of (ERIs in reciprocal space, frozen core energy) if calc_eri_energy is True
- Return type: