Skip to main content An official website of the United States government Here's how you know Here's how you know Official websites use .gov A .gov website belongs to an official government organization in the United States. Secure .gov websites use HTTPS A lock ( Lock Locked padlock icon ) or https:// means you've safely connected to the .gov website. Share sensitive information only on official, secure websites. Search Log in Dashboard Publications Account settings Log out Search… Search NCBI Primary site navigation Search Logged in as: Dashboard Publications Account settings Log in Search PMC Full-Text Archive Search in PMC Journal List User Guide PERMALINK Copy As a library, NLM provides access to scientific literature. Inclusion in an NLM database does not imply endorsement of, or agreement with, the contents by NLM or the National Institutes of Health. Learn more: PMC Disclaimer | PMC Copyright Notice J Chem Phys . Author manuscript; available in PMC: 2026 Apr 15. Published in final edited form as: J Chem Phys. 2026 Apr 14;164(14):144118. doi: 10.1063/5.0316578 Search in PMC Search in PubMed View in NLM Catalog Add to search LibppRPA : An Open-Source Library for Particle-Particle Random Phase Approximation Jincheng Yu Jincheng Yu 1) Department of Chemistry, Duke University, Durham, North Carolina 27708, United States 2) Department of Chemistry and Biochemistry, University of Maryland, College Park, Maryland 20742, United States Find articles by Jincheng Yu 1, 2, a) , Jiachen Li Jiachen Li 3) Department of Chemistry, Yale University, New Haven, Connecticut 06520, United States 1) Department of Chemistry, Duke University, Durham, North Carolina 27708, United States Find articles by Jiachen Li 3, 1, a) , Chaoqun Zhang Chaoqun Zhang 3) Department of Chemistry, Yale University, New Haven, Connecticut 06520, United States Find articles by Chaoqun Zhang 3 , Tianyu Zhu Tianyu Zhu 3) Department of Chemistry, Yale University, New Haven, Connecticut 06520, United States Find articles by Tianyu Zhu 3 , Weitao Yang Weitao Yang 1) Department of Chemistry, Duke University, Durham, North Carolina 27708, United States Find articles by Weitao Yang 1, b) Author information Copyright and License information 1) Department of Chemistry, Duke University, Durham, North Carolina 27708, United States 2) Department of Chemistry and Biochemistry, University of Maryland, College Park, Maryland 20742, United States 3) Department of Chemistry, Yale University, New Haven, Connecticut 06520, United States a) These two authors contributed equally b) Email: [email protected] PMC Copyright notice PMCID: PMC13078643 NIHMSID: NIHMS2147869 PMID: 41954251 The publisher's version of this article is available at J Chem Phys Abstract The accurate description of electron correlation and excitation energies remains a fundamental challenge in quantum chemistry. The particle-particle random phase approximation (ppRPA) has emerged as a promising method for capturing a broad range of excited-state properties. However, the implementation of ppRPA has been largely limited to in-house software, restricting its accessibility and usability. In this work, we present LibppRPA , an open-source and lightweight Python library designed for efficient and flexible ppRPA calculations of (1) electronic excitation energy and its associated analytical gradients and (2) the ground state correlation energy, and its associated analytical gradients. LibppRPA enables seamless integration with existing quantum chemistry packages, such as PySCF , by utilizing occupation numbers, molecular orbital coefficients, and three-center electron repulsion integrals. We implement both direct diagonalization and the iterative Davidson algorithm for solving the ppRPA equations, as well as active-space approximations, allowing users to balance accuracy and computational efficiency. We demonstrate the performance of LibppRPA through benchmark calculations on singlet-triplet gaps, double excitations, charge-transfer excitations, and valence/Rydberg excitations, showcasing its reliability across diverse molecular systems. The library provides a robust platform for studying electronic excitations and offers new opportunities for future developments in electronic structure theory. I. INTRODUCTION The accurate description of electron correlation and electronic excitations is a central task in theoretical chemistry. Over the past decades, significant efforts have been made to developing theoretical approaches that address these challenges, leading to a variety of methods with distinct strengths and limitations. Among the most widely used approaches are time-dependent density functional theory (TDDFT) ( 1 – 3 ) , the Bethe-Salpeter equation (BSE) formalism 4 – 6 , and wave-function-based methods 7 – 10 . TDDFT has become one of the most popular computational approaches for describing molecular systems due to its balance of computational efficiency and accuracy. The linear-response formulation of TDDFT has been widely implemented in modern quantum chemistry packages to calculate energies, structures, and other properties of excited states 11 – 17 . With the iterative Davidson algorithm 18 , the formal scaling of TDDFT is 𝒪 N 4 , where N is the size of the system. However, its performance can suffer for charge-transfer (CT) states and Rydberg states 19 , 20 , partly due to the incorrect long-range behavior for describing the potential energy surface of TDDFT with conventional density functional approximations (DFAs). Efforts to address these limitations have included the development of range-separated functionals 21 , 22 and the tuning of the Hartree-Fock (HF) exchange fraction in DFAs 23 , 24 . Furthermore, the performance of TDDFT is critically dependent on the selection of the exchange-correlation (XC) kernel 14 , which plays a central role in determining its accuracy. The BSE formalism 4 – 6 , derived from Green’s function theory or many-body perturbation theory, has also gained growing attention for computing optical spectrum of molecules, interfaces and solids. BSE is commonly performed on top of quasiparticle (QP) energies computed at the GW level, which is denoted as the BSE/ GW approach. By accurately capturing the long-range behavior and using the dynamic screening interaction for non-local electron correlations in real systems, the BSE/ GW approach provides reliable predictions of excitation energies for a wide range of systems 25 – 36 . However, there are several challenges for the BSE/ GW formalism. First, the accuracy of the BSE/ GW formalism heavily depends on the level of self-consistency in the GW calculation 30 , 37 – 39 , similar to the dependence on the DFA in TDDFT. Second, BSE/ GW suffers from the underscreening error, as a result of the missing vertex correction in the BSE kernel. Third, although BSE has the same 𝒪 N 4 scaling as the TDDFT, the proceeding GW calculation can be computationally demanding. Methods including obtaining QP energies from machine-learning frameworks 40 , 41 and combining BSE with generalized Kohn–Sham (KS) approaches, such as localized orbital scaling correction (LOSC) 42 and Koopmans–compliant functionals 43 , have been proposed to avoid the computational bottleneck. Additionally, the analytic gradient of BSE/ GW excitation energies is not available, where the relaxation of the excited-state structure can only be performed with approximations. Recently, the analytic gradients for BSE excitation energies from the Z-vector formalism have been developed 44 , 45 , which enable the rigorous relaxation of excited-state structures. Alongside the computationally affordable linear-response formalisms, highly accurate wave function methods are commonly used to predict excited energies of molecular systems. It has been shown that multireference methods such as complete active space second-order perturbation theory (CASPT2) 46 , multireference configuration interaction (MRCI) 47 can predict different excited states on the equal footing. However, the selection of the active space in multireference methods can be ambiguous. Another path is using single-reference wave function methods such as algebraic diagrammatic construction 48 and coupled cluster techniques 49 . These methods offer a systematically improvable path to increase the accuracy by including high-order expansions. However, due to their high computational cost, wave function methods are typically reserved for benchmark studies 9 , 50 , 51 . Recently the particle-particle random phase approximation (ppRPA) has gained increasing attentions for ground-state and excited-state properties of molecules and solid-state materials. ppRPA was originally developed to calculate the nuclear many-body correlation energy 52 , 53 , and has been extended to describe correlation in electronic systems by the Yang group 54 , 55 . It can be derived from different approaches, including the adiabatic connection using the pairing matrix fluctuations 54 , 55 TDDFT with the pairing field 56 , and the equation of motion (EOM) 53 , 57 , and particle-particle channel in BSE 58 . Commonly used particle-hole random phase approximation 59 , 60 (phRPA) contains information in the the particle-hole channel, which has been extensively implemented in quantum chemistry packages and open-source library such as LibRPA 61 for electron correlation energy. Compared to phRPA, ppRPA contains information in the particle-particle and the hole-hole channels, which associates with the fluctuation of the pairing matrix, which is the functional derivative of the pairing matrix with respect to the pairing field, taken at the limit of vanishing pairing field 54 . For the ground-state energy, ppRPA is equivalent to ladder couple cluster doubles when the Hartree-Fock approximation is used for the nonintereting reference system (CCD) 62 , 63 , which is exact up to the second order with the ladder partial summation. Recently, the extension of ppRPA based on a multireference character has been developed to describe the strong correlation in molecular systems 64 . ppRPA has also been extended to obtain charge-neutral excitation energies by the Yang group 65 , 66 : it calculates the two-electron addition and the two-electron removal energies, which can be considered as an approximation to double-electron-affinity or double-ionization-potential EOM-CCD. By taking the differences between the two-electron addition energies of the ( N - 2 ) -electron system 65 or the two-electron removal energies of the ( N + 2 ) -electron system 65 , 67 , 68 , excitation energies of the N -electron system can be obtained. ppRPA has gained great success in describing various ground-state and excited-state properties. For example, ppRPA can accurately predict a variety types of excitation energies including valence excitations 69 , singlet-triplet (S-T) energy gaps of diradicals 70 , 71 , charge-transfer (CT) excitations 72 , double excitations 65 , 73 , and Rydberg excitations 65 . Conical intersections 74 , oscillator strengths 65 , and excited-state geometries 75 can also be well described by ppRPA. More recently, with the development of the active-space formalism 76 , 77 , the computational cost of ppRPA can be significantly reduced without loss of accuracy. With this development, ppRPA has also been applied to predict accurate excitation energies of point defects 78 , 79 . In addition, ppRPA with the Tamm-Dancoff approximation (TDA) has been applied in the multireference DFT to predict dissociation energies and excitation energies 80 , 81 . In the context of the Green’s function theory, the ppRPA eigenstates are also utilized in the T-matrix self-energy 82 to calculate QP energies 83 – 85 , as a counterpart of phRPA used to construct the GW self-energy, The success of ppRPA can be attributed to several reasons. First, ppRPA can be viewed as an embedding approach in the Fock space. Two frontier electrons are treated in an exact way using subspace CI with a seamless integration of DFT for the remaining ( N - 2 ) electrons 71 , 77 . Thus, for diradicals and a large number of point defects, ppRPA describes the strong correlation of the two added electrons, free from the static correlation error and the spin contamination. Second, ppRPA contains the two-particle information by its construction. Therefore, ppRPA is naturally capable of describing double excitations 65 , 73 , which cannot be captured by linear-response particle-hole formalisms such as TDDFT and BSE with the adiabatic approximation. Third, ppRPA kernel presents the correct long-range asymptotic behavior, so it can accurately predict CT and Rydberg excitations 65 , 72 . The orbital energies involved for ground and excited states are all from either the unoccupied space in particle-particle channel or the occupied space in the hole-hole channel in the corresponding ( N - 2 ) -electron or ( N + 2 ) -electron reference system described with DFA. Thus the excitation energy for the N -electron system do not involve the HOMO-LUMO gaps, or the fundamental band gap, which has a major dependence on the DFA. In contrast, the HOMO-LUMO gap directly enters into the TDDFT equations. Therefore, the starting point dependence of ppRPA is much smaller compared to TDDFT. 69 , 72 , whose accuracy largely depends on the fundamental gap at the DFT level and the form of the XC kernel. In this work, we develop a reliable and stable implementation for ppRPA in LibppRPA , which will serve as an open-source and lightweight library conducting ppRPA calculations for electronic and magnetic properties of molecular and periodic systems ( Γ -point supercell approach). With the implementation completely in Python, we aim to provide LibppRPA users with great flexibility to perform ppRPA calculations based on data from any quantum chemistry packages of their choice. To demonstrate this we integrate LibppRPA with PySCF 86 , 87 . In principle, LibppRPA can take input generated from any Gaussian type orbital (GTO) based software, which provides immediate access for researchers to perform calculations with ppRPA. II. THEORY A. Pairing Matrix The ppRPA formalism can be derived from different approaches, including the equation of motion 53 , 57 , the adiabatic connection 54 , 55 , TDDFT with the pairing field 56 , and the BSE in the Green’s function theory 88 . As a counterpart to phRPA formulated with the fluctuation of the density-density response function, the ppRPA formalism is formulated with the fluctuation of the pairing matrix 52 – 55 κ x 1 , x 2 , t = Ψ 0 N ψ ˆ x 2 , t ψ ˆ x 1 , t Ψ 0 N (1) where x = ( r , σ ) is the combined space-spin variable, Ψ 0 N is the N -electron ground state, ψ ˆ † and ψ ˆ are the second quantization creation and annihilation operator. In the absence of an external pairing field, the pairing matrix κ x 1 , x 2 , t is zero, but its linear response to an external pairing field F ˆ x 1 , x 2 , t = f ( t ) ψ ˆ x 1 , t ψ ˆ x 2 , t is not zero, which is an external field coupled to the pairing matrix. Then the retarded dynamic pairing matrix fluctuation K ‾ t - t ′ is 54 K ¯ pqrs t - t ′ = - i θ t - t ′ Ψ 0 N a ˆ p ( t ) a ˆ q ( t ) , a ˆ s † t ′ a ˆ r † t ′ Ψ 0 N (2) which has poles at two-electron addition and removal energies that characterize double electron addition and ionization processes in the frequency space K ‾ ( ω ) pqrs = ∑ m Ψ 0 N a ˆ p a ˆ q Ψ m N + 2 Ψ m N + 2 a ˆ s † a ˆ r † Ψ 0 N ω - Ω m N + 2 + i η - ∑ m Ψ 0 N a ˆ s † a ˆ r † Ψ m N - 2 Ψ m N - 2 a ˆ p a ˆ q Ψ 0 N ω - Ω m N - 2 + i η (3) Here a ˆ p † and a ˆ p are the second quantization creation and annihilation operator for the orbital p , Ω N ± 2 is the two-electron addition/removal energy, η is a positive infinitesimal number. We use i , j , k , l for occupied orbitals, a , b , c , d for virtual orbitals, p , q , r , s for general orbitals, and m for the index of the two-electron addition/removal energy. In many-body perturbation theory, the time-ordered form is normally used for theoretical derivation. The time-ordered dynamic pairing matrix fluctuation is 54 K pqrs t - t ′ = - i Ψ 0 N T a ˆ p ( t ) a ˆ q ( t ) , a ˆ s † t ′ a ˆ r † t ′ Ψ 0 N (4) and in the frequency space is K ( ω ) pqrs = ∑ m Ψ 0 N a ˆ p a ˆ q Ψ m N + 2 Ψ m N + 2 a ˆ s † a ˆ r † Ψ 0 N ω - Ω m N + 2 + i η - ∑ m Ψ 0 N a ˆ s † a ˆ r † Ψ m N - 2 Ψ m N - 2 a ˆ p a ˆ q Ψ 0 N ω - Ω m N - 2 - i η (5) The dynamic paring matrix fluctuation, as shown in Eq. 5 clearly contains the double electron addition and removal excitation energies. It has also been shown that the dynamic paring matrix fluctuation provides a rigorous formulation for electron correlation in the particle-particle channel 54 , 55 , parallel to the electron correlation formulation in the particle-hole channel through the dynamic density density fluctuation, or the dynamic density response function 89 . To obtain the pairing matrix fluctuation K of the interacting system, the ppRPA approximates K in terms of the non-interacting K 0 by the Dyson equation 54 , 55 K = K 0 + K 0 V K (6) where the antisymmetrized interaction V pqrs = ⟨ p q | | r s ⟩ = ⟨ p q ∣ r s ⟩ - ⟨ p q ∣ s r ⟩ with ⟨ p q ∣ r s ⟩ = ∫ d x d x ′ ψ p * ( x ) ψ r ( x ) ψ q * x ′ ψ s x ′ r - r ′ . In Ref. 88 , the screened interaction is used in the Dyson equation, which leads to the BSE in the particle-particle channel. B. ppRPA Equation The Dyson equation in Eq. 6 can be cast into a generalized eigenvalue problem that gives two-electron addition and removal energies 54 , 69 A B B † C X Y = Ω I 0 0 - I X Y (7) with A a b , c d = δ a c δ b d ϵ a + ϵ b + ⟨ a b ‖ c d ⟩ (8) B a b , k l = ⟨ a b ‖ k l ⟩ (9) C i j , k l = - δ i k δ j l ϵ i + ϵ j + ⟨ i j ‖ k l ⟩ (10) where a < b , c < d , i < j , k < l , the eigenvalue Ω is the two-electron addition/removal energy, X and Y are the two-electron addition and removal eigenvectors. For closed-shell systems, the ppRPA equation in Eq. 7 can be cast into the spin-adapted form 90 . For singlet excitations the ppRPA matrix elements are A a b , c d s = δ a c δ b d ϵ a + ϵ b + 1 1 + δ a b 1 + δ c d ( ⟨ a b ∣ c d ⟩ + ⟨ a b ∣ d c ⟩ ) (11) B a b , k l s = 1 1 + δ a b 1 + δ k l ( ⟨ a b ∣ k l ⟩ + ⟨ a b ∣ l k ⟩ ) (12) C i j , k l s = - δ i k δ j l ϵ i + ϵ j + 1 1 + δ i j 1 + δ k l ( ⟨ i j ∣ k l ⟩ + ⟨ i j ∣ l k ⟩ ) (13) with a ≤ b , c ≤ d , i ≤ j and k ≤ l . And for triplet excitations the ppRPA matrix elements are A a b , c d t = δ a c δ b d ϵ a + ϵ b + ⟨ a b ‖ c d ⟩ (14) B a b , k l t = ⟨ a b ‖ k l ⟩ (15) C i j , k l t = - δ i k δ j l ϵ i + ϵ j + ⟨ i j ‖ k l ⟩ (16) with a < b , c < d , i < j and k < l . ppRPA eigenvalues obtained from Eq. 7 can be considered as an approximation to double-electron-affinity or double-ionization-potential equation-of-motion coupled-cluster doubles 65 , 66 . The ppRPA eigenvector is normalized as 54 , 55 X m , † X m - Y m , † Y m = ± 1 (17) where the upper sign is for two-electron addition excitations and the lower sign is for two-electron removal excitations. The neutral excitation energy of an N -electron system can be obtained from the energy difference between two-electron addition energies of the corresponding ( N - 2 ) -electron system, or the energy difference between two-electron removal energies of the corresponding ( N + 2 ) -electron system. Similar to TDDFT and BSE, the scaling of solving the ppRPA equation in Eq. 7 is 𝒪 N 4 with the Davidson algorithm 69 . In the recently developed activespace formalism, ppRPA excitation energies are shown to converge rapidly with canonical active-space orbitals for both molecular and periodic systems 76 , 78 , 79 , which significantly lowers the computational cost. C. ppRPA for Fractional Charge Systems With the extension of the Green’s function theory to fractional charge systems established in Ref. 91 , the ppRPA equation in Eq. 7 can be applied for fractional charge system with matrix elements defined as A a b , c d = δ a c δ b d ϵ a + ϵ b + 1 - n a 1 - n b 1 - n c 1 - n d ⟨ a b ‖ c d ⟩ (18) B a b , k l = 1 - n a 1 - n b n k n k ⟨ a b ‖ k l ⟩ (19) C i j , k l = - δ i k δ j l ϵ i + ϵ j + n i n j n k n k ⟨ i j ‖ k l ⟩ (20) where n is the occupation number. The total energy from Hartree-Fock energy with ppRPA correlation energy shows a linear energy behavior for systems with fractional charges and meets the flat-plane condition 54 . The fractional formulation of ppRPA also allows one to obtain the ppRPA chemical potential with the finite difference approach, which significantly outperforms phRPA for predicting ionization potentials and electron affinity 54 . D. ppRPA for Relativistic Hamiltonian Li and co-workers have extended particle-particle TDA method to treat relativistic two-component Hamiltonian 92 . ppRPA for relativistic Hamiltonian is also implemented in LibppRPA , compatible with relativistic generalized mean-field reference using spin orbital basis functions or j-adapted spinors. Various exact two-component theory (X2C) 93 – 96 Hamiltonians are supported, such as the X2C theory in its one-electron variant (X2C-1e) and the X2C-1e with atomic mean-field two-electron interaction (the X2CAMF scheme) 97 – 99 . Two-component ppRPA calculations can be run seamlessly using the interface to PySCF . We note that the relativistic Hamiltonian is treated variationally, i.e., added at the mean-field level. With the no-pair approximation, the ppRPA working equations remain the same as those for the non-relativistic calculations except for the lack of spin symmetry. E. Total Energy The ppRPA correlation energy can be obtained in terms of the solution to generalized eigenvalue problem in Eq. 7 as 54 , 55 E pp , c = ∑ Ω + 2 e - TrA = - ∑ Ω - 2 e - TrC (21) where T r means trace, Ω ± 2 e is two-electron addtion/removal energy. The expression in Eq. 21 is the same for fractional charge systems. For ground states, the ppRPA correlation energy can alternatively be interpreted as the sum of all ladder diagrams in the wavefunction theory, which is equivalent to the ladder-coupled-cluster doubles 62 , 63 . Then the ppRPA total energy is expressed as E pp = E HF + E pp , c (22) where E HF is the Hartree-Fock energy evaluated with the input orbitals. Alternatively, the ppRPA total energy can be obtained from the multireference DFT method 80 , 81 , which means E m pp , MR = E DFT , N ± 2 + Ω m ∓ 2 e (23) where E DFT , N ± 2 is the DFT energy of the ( N ± 2 ) -electron system and Ω ± 2 e is the two-electron addition/removal energy from ppRPA. Eq. 23 describes both ground and excited states on equal footing. As shown in Refs. 80 , 81 , the SCF solution of the multireference ppRPA energy can be achieved with the (generalized) optimized effective potential method 100 , 101 . F. Natural Transition Orbital The natural transition orbitals (NTOs) of ppRPA is developed to provide qualitative descriptions of electronic transitions 79 . For the NTOs of particle-hole formalisms, the dominant particle-hole pairs of an excited state are obtained from the SVD of the corresponding transition density matrix 102 . Similarly, NTOs in ppRPA convey information about particle-particle pairs and hole-hole pairs. For the m -th state, the two-electron addition eigenvector X m can be viewed as a triangular matrix of dimension N vir × N vir and the two-electron removal eigenvector Y m can be viewed as a triangular matrix with of dimension N occ × N occ . Therefore, the coefficients that transform molecular orbitals to natural transition orbitals of two particles (holes) can be obtained with the SVD of X m Y m X m = C p 1 , m λ p , m C p 2 , m † (24) Y m = C h 1 , m λ h , m C h 2 , m † (25) In Eq. 24 , the NTO coefficient matrices C p 1 and C p 2 of dimension N vir × N vir are associated with particle-particle pairs for adding two electrons, which are weighted with diagonal elements of matrix λ p . Similarly, in Eq. 25 , the NTO coefficient matrices C h 1 and C h 2 of dimension N occ × N occ are associated with hole-hole pairs for removing the two electrons, which are weighted with diagonal elements of matrix λ h . As a consequence of the normalization in Eq. 17 , the NTO weights satisfy the following relation ∑ a vir λ a p - ∑ i occ λ i h = ± 1 (26) where the upper sign is for two-electron addition excitations and the lower sign is for two-electron removal excitations. The resulting NTO weights ( λ a p and λ i h ) can thus be employed to qualitatively analyze the components and multireference character of the associated ground and excited states. G. Density Matrix As shown in Ref. 103 , the density of a DFA can be obtained as the functional derivative of the energy with respect to the external potential ρ ( r ) = δ E δ v ext ( r ) (27) which is used for deriving phRPA and ppRPA density 27 . Here we generalize the expression from the density to the one-particle density matrix, whose occupied-occupied and the virtual-virtual blocks in the orbital space are defined as D i j = ⟨ Ψ | a ˆ i a ˆ j † | Ψ ⟩ (28) D a b = ⟨ Ψ | a ˆ a † a ˆ b | Ψ ⟩ (29) In ppRPA, the energy is obtained as the summation of the energy of the ( N ± 2 ) -electron reference system calculated by DFT and the two-electron addition/removal energy calculated by ppRPA. Therefore, the total density matrix of the m -th state of the N -electron system is divided into two parts: the non-interacting part of the ( N ± 2 ) -electron reference system, and the two-electron addition/removal part D N , m = D N ± 2 , KS + D 2 e , m (30) which has a zero occupied-virtual block and agrees with the ppRPA density derived in Ref. 103 in the real-space diagonal limit. In Eq.30 , the two-electron addition/removal part is constructed with the ppRPA eigenvectors solved from Eq.7 , which contains information of adding and removing two electrons. For the m -th excitation, the virtual-virtual block for adding two electrons is D a b 2 e , m = X m X m , † a b + X m , † X m a b (31) and the occupied-occupied block for removing two electrons is D i j 2 e , m = - Y m Y m , † i j - Y m , † Y m i j (32) Then the two-electron addition/removal density matrix stratifies the following relation T r D 2 e , m = ± 2 (33) where the upper sign is for two-electron addition excitations and the lower sign is for two-electron removal excitations in Eq. 33 . H. Analytic Gradient To obtain the energy gradient in ppRPA, the total energy is expressed in the multireference DFT manner by Eq.23 , then the gradient can be calculated as 75 ∂ E m pp , MR ∂ λ = ∂ E DFT , N ± 2 ∂ λ + ∂ Ω m ∓ 2 e ∂ λ (34) where λ is the nuclear coordinate. By virtue of the Hellmann–Feynman theorem, the gradient of the two-electron addition/removal energy with respect to the nuclear coordinate is 75 ∂ Ω m ∓ 2 e ∂ λ = X m † Y m † ∂ ∂ λ A B B † C X m Y m (35) where more details can be found in Ref. 75 . III. IMPLEMENTATION DETAILS The LibppRPA library is an open-source and pure-Python package for performing ppRPA calculations within the Gaussian type orbital (GTO) framework. The workflow is illustrated in Fig.1 . LibppRPA takes three main input variables: a) occupation numbers, b) molecular orbital (MO) energies and c) the three-center density-fitting integrals in the MO space, which is constructed from the resolution of identity technique (Coulomb or overlap matrix norm) or the Cholesky decomposition of the two-electron integrals, in either a molecular or periodic (supercell with Γ -point sampling) calculation. The input variables are passed to LibppRPA in the numpy.ndarray format, which can be generated on the fly or read via HDF5 binary data file from popular quantum chemistry software and in-house programs. LibppRPA also offers a convenient interface that directly initializes the ppRPA calculation from a mean-field calculation by PySCF 86 , 87 . FIG. 1. Open in a new tab The architecture of LibppRPA . After loading the input variables, LibppRPA provides two routines to solve the ppRPA equation: direct diagonalization and the Davidson algorithm 69 . Both routines use the density-fitting (resolution-of-identity) technique to compute electron repulsion integrals and can be combined with the active-space formalism 76 . For symmetry preserved spin-restricted ppRPA, the ppRPA matrix is constructed in a spin-adapted manner 90 . For spin-unrestricted ppRPA, the ppRPA matrix is constructed and diagonalized in three subspaces ( α , α ; α , α ) , ( α , α ; β , β ) , and ( β , β ; β , β ) 90 . As a lightweight package, the core routines in LibppRPA only depend on standard Python packages such as Numpy 104 and Scipy 105 for efficient mathematical operations, which are parallelized with the OpenMP scheme. After solving the ppRPA equation, ppRPA two-electron addition/removal energies and eigenvectors in the numpy.ndarray format are returned as the output. LibppRPA provides several analysis tools for ground-state and excited-state properties such as NTOs, density matrix and oscillator strengths. For NTOs and density matrix, the results can be saved in the Gaussian Cube format 106 for further visualizations. LibppRPA can also take the output to calculate analytic gradients for geometry optimizations, which is integrated with PySCF to obtained required integrals. IV. EXAMPLE USAGE In this section the usage of LibppRPA for calculating excitation energies is demonstrated. We begin with the particle-particle channel. In Listing 1 , the calculation of charge-neutral excitation energies of H 2 O in the particle-particle channel is shown. The mean-field calculation of the ( N - 2 ) -electron system is performed by PySCF first. Then the input variables for ppRPA (the occupation number, orbital energies and the density-fitting matrix) are obtained from the mean-field object in PySCF . Finally the spin-adapted ppRPA calculations in the particle-particle channel are carried out. Listing 1. Example of using LibppRPA to calculate excitation energies of H 2 O in the particle-particle channel. Open in a new tab Similar to the ppRPA calculation in the particle-particle channel, the example for the ppRPA calculation in the hole-hole channel of H 2 O is shown in the comment blocks of Listing 1 . The mean-field calculation of the ( N + 2 ) -electron system is performed by PySCF to generate input variables for LibppRPA . Then the spin-adapted ppRPA calculations are performed in the hole-hole channel. V. RESULTS A. Computational Details To demonstrate the performance of the LibppRPA library, we reproduced a series of calculations in the literature: S-T gaps of diradical systems 70 , double excitation of molecular systems 73 , charge transfer excitation of the Stein’s set 72 , Rydberg excitation energies of atomic systems 76 , valence excitation energies of the Thiel’s set 69 , vertical excitation energies of point defects 78 , 79 , reaction barriers in DBH24 set 90 , dissociation curves of H 2 and Ar 2 54 , potential energy surfaces of small molecules 75 . The input SCF results and integrals for LibppRPA can be generated from any quantum chemistry package. In this work, we performed all ground-state SCF calculations in Gaussian basis sets with Gaussian density fitting using the PySCF quantum chemistry software package 86 , 87 . B. Excitation Energy The mean absolute errors (MAEs) of ppRPA for predicting excitation energies of different characters, including 24 singlet-triplet gaps of diradicals, 29 double excitations in Loos’s set 107 , 12 charge-transfer excitations in Stein’s set 21 , 8 Rydberg excitations of atoms, 38 valence excitations in Thiel’s set 108 , and 13 defect excitations are presented in Table I . It shows that ppRPA is capable with predicting accurate excitation energies for a broad range of systems and has a small starting point dependence on the chosen XC functional. In particular, for singlet-triplet gaps, double excitations and defect excitations, ppRPA significantly outperforms TDDFT 70 , 73 , 78 , 79 . For singlet-triplet gaps of diradicals, ppRPA provides errors of 3 to 5 kcal/mol compared to experiment references 70 , which outperforms spin-flipped TDDFT with errors of 3 to 10 kcal/mol 109 . For double excitations, ppRPA provides errors around 0.4 eV compared to the theoretical best estimations, similar to the wavefunction method CCSDT with around 0.3 eV errors 107 , 110 . For defect excitations, ppRPA gives errors around 0.2 eV compared to experiment references 78 , 79 , where TDDFT based on conventional functionals gives larger errors around 0.4 eV 15 . For charge-transfer and valence excitations, the accuracy of ppRPA is similar or slightly better than TDDFT based on conventional functionals 69 , 72 . In addition to the good accuracy, the computational cost can be largely reduced by combining with the active-space formalism, which is implemented in LibppRPA for both full diagonalization and Davidson algorithm. As shown in Ref. 76 , ppRPA can be solved with a small active space without loss of accuracy. With 200 occupied and 200 virtual orbitals in the active space, the ppRPA calculation of the NV − center in diamond with the 216-atom supercell only takes 1757 seconds on a 48-CPU node (Intel Xeon Gold 6442Y, 4000MHz), where the ground-state B3LYP calculation takes 35951 seconds TABLE I. Mean absolute errors (MAEs) of ppRPA based on HF, PBE and B3LYP for predicting different types of excitation energies. ppRPA@HF ppRPA@PBE ppRPA@B3LYP singlet-triplet gap (kcal/mol) 17.5 3.9 4.7 double excitation (eV) 0.38 0.39 charge-transfer excitation (eV) 0.51 0.86 0.72 valence excitation (eV) 0.79 0.39 0.37 Rydberg excitation (eV) 0.08 2.13 2.53 defect excitation (eV) 0.20 0.12 Open in a new tab C. Atomic Zero-Field Splitting Relativistic ppRPA is a natural choice for calculating spin properties for systems with two open-shell systems. We report the calculated atomic zero-field splittings (ZFS) for P 3 states of carbon group elements using X2CAMF-ppRPA method with uncontracted ANO-RCC basis sets 111 . Two-electron Coulomb interaction and Gaunt term have been included in X2CAMF Hamiltonian. All the calculations start with N-2 reference state. The occupied valence n s orbital and virtual orbitals below 10 hartree are included in the ppRPA calculations. The results are summarized in Table II . ppRPA@HF gives the best agreement with the experimental values. Relative errors are within 10%. However, ppRPA@PBE and ppRPA@B3LYP tend to overestimate the ZFS by approximately 30%. This might be due to the use of collinear functional. TABLE II. Calculated and experimental atomic zero-field splittings (ZFS) for P 3 states of carbon group elements using X2CAMF-ppRPA based on HF, PBE and B3LYP. The numbers are reported as the relative energy levels for J=1/J=2 states with respect to J=0 state. ZFS ( c m - 1 ) ppRPA@HF ppRPA@PBE ppRPA@B3LYP exp C 14/42 21/62 20/59 16/43 Si 73/211 104/300 99/287 77/223 Ge 510/1295 806/1958 752/1842 557/1410 Sn 1556/3212 2429/4688 2279/4446 1692/3428 Pb 7191/10024 10671/14237 10051/13518 7819/10650 Open in a new tab D. Total Energy In this section, the performance of LibppRPA for the total energy is shown. Dissociation curves of single-bond breaking in H 2 molecule and weakly interacting Ar 2 molecules obtained from ppRPA are shown on the left-hand side of Fig.2 . It shows that ppRPA@HF provides a similar dissociation curve to phRPA@HF for the H 2 dissociation, which overestimates the dissociation energies compared to the CCSD reference. Because H 2 only has two electrons, MR-ppRPA is exact and gives the same results as the CCSD reference. For other single-bond breaking problems, MR-ppRPA treats the two valence electrons in a subspace configuration interaction fashion and is free from the static correlation error, which has been shown to predict accurate dissociation energies and equilibrium bond lengths 80 , 81 . FIG. 2. Open in a new tab Left: dissociation curve of H 2 obtained from phRPA@HF, ppRPA@HF, ppRPA@HF (multireference DFT) and CCSD. Right: behaviors of PBE, phRPA@HF and ppRPA@HF total energies of Be as a function of the electron number. ppRPA@HF and MR-ppRPA@HF calculations were performed with LibppRPA . The cc-pVDZ basis set was used. On the right-hand side of Fig.2 , total energy obtained from ppRPA for fractional charge systems is shown. In the context of DFT, ppRPA is the first known functional that captures the energy derivative discontinuity in strongly correlated systems. As shown in Fig.2 , both phRPA and conventional functionals like PBE give convex curves for total energies between integer electron numbers, which originate from large delocalization errors and lead to huge errors for predicting chemical potentials 112 , 113 . ppRPA has no delocalization error with a nearly linear energy behavior for systems with fractional charges. It captures the derivative discontinuity of the energy at integer electron numbers and meets the flat-plane condition 114 , which gives good accuracy for predicting ionization potentials and electron affinity using the finite-difference approach 54 . E. Analysis and Visualization The ppRPA results can be further visualized and analyzed by NTOs and electron densities. NTOs of the triplet ground state in ⋅ CH 2 CH 2 4 C CH 3 H ⋅ obtained from ppRPA@B3LYP are shown in Fig.3 . ppRPA provides good descriptions for this long disjoint diradicals. It shows that two nonbonding electrons are added to carbon atoms on each chain end connected by a CH 2 n bridge. The corresponding NTO weight is 0.998, which means it can be properly described by a single-determinant approach and indicates the strong diradical character in ⋅ CH 2 CH 2 4 C CH 3 H ⋅ . FIG. 3. Open in a new tab Natural transition orbitals (NTOs) of the triplet ground state in ⋅ CH 2 CH 2 4 C CH 3 H ⋅ obtained from ppRPA@B3LYP. The NTO weight is 0.998. Top: NTO of adding the first nonbonding electron. Bottom: NTO of adding the second nonbonding electron. The aug-cc-pVTZ basis set was used. Geometry was taken from Ref. 70 . The isosurface value is 0.04 a.u. Then we visualize the charge transfer nature of the excitation in the anthracene-tetracyanoethylene (TCNE) system with the ppRPA density matrix. The electron density difference between the first singlet excited state and the ground state of anthracene-TCNE obtained from ppRPA@B3LYP in Fig.4 . It shows that the electron density flows from the aromatic donor anthracene to the TCNE acceptor. FIG. 4. Open in a new tab Electron density difference between the first singlet excited state and the ground state of anthracene-TCNE obtained from ppRPA@B3LYP. Yellow and blue indicate negative and positive electron densities, respectively. The cc-pVDZ basis set was used. Geometry was taken from Ref. 115 . The isosurface value is 0.003. F. Geometrical Gradient and Structural Optimization With analytical gradient techniques 75 , LibppRPA can perform structural optimization using ppRPA methods. The ∑ g − 3 , Δ g 1 , and ∑ g + 1 states of oxygen molecule have been optimized to demonstrate the applicability of the ppRPA analytical gradient methods. ppRPA@HF method was used with cc-pVDZ basis sets and N-2 reference state. The optimized bond lengths and the corresponding potential energy surfaces are plotted in Fig. 5 . The new implementation features several advances beyond the original publication 75 . The details will be reported in the future work. FIG. 5. Open in a new tab Optimized structures for ∑ g − 3 , Δ g 1 , and ∑ g + 1 states of oxygen molecule together with the corresponding potential energy surfaces. The optimized structures are denoted as red cross. The equilibrium bond lengths are 1.1663, 1.1681, and 1.1708 Å, for ∑ g − 3 , Δ g 1 , and ∑ g + 1 states, respectively. VI. CONCLUSIONS In this work, we develop LibppRPA , an open-source and lightweight library for conducting ppRPA calculations. Implemented entirely in Python, LibppRPA provides a flexible and userfriendly framework that can seamlessly integrates with existing quantum chemistry packages such as PySCF . By leveraging efficient matrix solvers, including direct diagonalization and the Davidson algorithm, the library enables accurate and scalable ppRPA calculations. Additionally, the incorporation of active-space approximations significantly reduces computational cost while maintaining high accuracy, making ppRPA more accessible for a wider range of molecular systems. Through extensive benchmark studies, we have demonstrated the reliability and versatility of LibppRPA in predicting excitation energies, including S-T gaps, double excitations, charge-transfer excitations, and Rydberg excitations. Our results highlight the advantages of ppRPA over conventional particle-hole approaches, particularly in describing multireference character and strong correlation effects. With its open-source nature and ease of integration, LibppRPA provides a robust platform for researchers to explore electronic excitations and advance the development of electronic structure methods. We anticipate that LibppRPA will serve as a valuable tool for the quantum chemistry community, facilitating both fundamental research and practical applications in excited-state calculations. Supplementary Material Supporting Information NIHMS2147869-supplement-Supporting_Information.pdf (173.7KB, pdf) See the Supporting Information for numerical results of excitation energies and groundstate energies obtained from ppRPA. ACKNOWLEDGMENTS J.Y. and W.Y. acknowledge the support from the National Institutes of Health (R35GM158181). T.Z., J.L., and C. Z. are supported by the National Science Foundation (Grant No. OAC-2513473). J.L. also acknowledges support from the Tony Massini Postdoctoral Fellowship in Data Science from Yale University. DATA AVAILABILITY STATEMENT Data and scripts pertaining to this work have been archived in the Duke Research Data Repository 116 . REFERENCES 1. Runge E and Gross EKU, Phys. Rev. Lett 52, 997 (1984). [ Google Scholar ] 2. Casida ME, in Recent Advances in Density Functional Methods, Recent Advances in Computational Chemistry, Vol. Volume 1 (WORLD SCIENTIFIC, 1995) pp. 155–192. [ Google Scholar ] 3. Ullrich CA, Time-Dependent Density-Functional Theory: Concepts and Applications (OUP Oxford, 2011). [ Google Scholar ] 4. Salpeter EE and Bethe HA, Phys. Rev 84, 1232 (1951). [ Google Scholar ] 5. Sham LJ and Rice TM, Phys. Rev 144, 708 (1966). [ Google Scholar ] 6. Hanke W and Sham LJ, Phys. Rev. Lett 43, 387 (1979). [ Google Scholar ] 7. Watson TJJ, Lotrich VF, Szalay PG, Perera A, and Bartlett RJ, J. Phys. Chem. A 117, 2569 (2013). [ DOI ] [ PubMed ] [ Google Scholar ] 8. Loos P-F, Matthews DA, Lipparini F, and Jacquemin D, J. Chem. Phys 154, 221103 (2021). [ DOI ] [ PubMed ] [ Google Scholar ] 9. Véril M, Scemama A, Caffarel M, Lipparini F, Boggio-Pasqua M, Jacquemin D, and Loos P-F, WIREs Comput. Mol. Sci 11, e1517 (2021). [ Google Scholar ] 10. Bartlett RJ, Phys. Chem. Chem. Phys 26, 8013 (2024). [ DOI ] [ PubMed ] [ Google Scholar ] 11. Jacquemin D, Perpète EA, Scuseria GE, Ciofini I, and Adamo C, J. Chem. Theory Comput 4, 123 (2008). [ DOI ] [ PubMed ] [ Google Scholar ] 12. Casida ME, J. Mol. Struct: THEOCHEM Time-Dependent Density-Functional Theory for Molecules and Molecular Solids, 914, 3 (2009). [ Google Scholar ] 13. Yuen-Zhou J, Tempel DG, Rodríguez-Rosario CA, and Aspuru-Guzik A, Phys. Rev. Lett 104, 043001 (2010). [ DOI ] [ PubMed ] [ Google Scholar ] 14. Laurent AD and Jacquemin D, Int. J. Quantum Chem 113, 2019 (2013). [ Google Scholar ] 15. Jin Y, Yu V. W.-z., Govoni M, Xu AC, and Galli G, J. Chem. Theory Comput 19, 8689 (2023). [ DOI ] [ PubMed ] [ Google Scholar ] 16. Knysh I, Raimbault D, Duchemin I, Blase X, and Jacquemin D, J. Chem. Phys 160, 144115 (2024). [ DOI ] [ PubMed ] [ Google Scholar ] 17. Poh YR, Morozov D, Kazmierczak NP, Hadt RG, Groenhof G, and Yuen-Zhou J, J. Am. Chem. Soc 146, 15549 (2024). [ DOI ] [ PubMed ] [ Google Scholar ] 18. Stratmann RE, Scuseria GE, and Frisch MJ, J. Chem. Phys 109, 8218 (1998). [ Google Scholar ] 19. Tozer DJ, J. Chem. Phys 119, 12697 (2003). [ Google Scholar ] 20. Dreuw A, Weisman JL, and Head-Gordon M, J. Chem. Phys 119, 2943 (2003). [ Google Scholar ] 21. Stein T, Kronik L, and Baer R, J. Chem. Phys 131, 244119 (2009). [ DOI ] [ PubMed ] [ Google Scholar ] 22. Refaely-Abramson S, Baer R, and Kronik L, Phys. Rev. B 84, 075144 (2011). [ Google Scholar ] 23. Brückner C and Engels B, Chem. Phys Electrons and Nuclei in Motion - Correlation and Dynamics in Molecules (on the Occasion of the 70th Birthday of Lorenz S. Cederbaum), 482, 319 (2017). [ Google Scholar ] 24. Gong J, Lam JWY, and Zhong Tang B, Phys. Chem. Chem. Phys 22, 18035 (2020). [ DOI ] [ PubMed ] [ Google Scholar ] 25. Blase X, Duchemin I, Jacquemin D, and Loos P-F, J. Phys. Chem. Lett 11, 7371 (2020). [ DOI ] [ PubMed ] [ Google Scholar ] 26. Monino E and Loos P-F, J. Chem. Theory Comput 17, 2852 (2021). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 27. Jiang X, Zheng Q, Lan Z, Saidi WA, Ren X, and Zhao J, Sci. Adv 7, eabf3759 (2021). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 28. Knysh I, Duchemin I, Blase X, and Jacquemin D, J. Chem. Phys 157, 194102 (2022). [ DOI ] [ PubMed ] [ Google Scholar ] 29. Cho Y, Bintrim SJ, and Berkelbach TC, J. Chem. Theory Comput 18, 3438 (2022). [ DOI ] [ PubMed ] [ Google Scholar ] 30. Li J, Golze D, and Yang W, J. Chem. Theory Comput 18, 6637 (2022). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 31. Ud Din N and Liu Z-F, J. Phys. Chem. C 126, 20694 (2022). [ Google Scholar ] 32. Wu J, Hou B, Li W, He Y, and Qiu DY, Phys. Rev. B 110, 075133 (2024). [ Google Scholar ] 33. Bhattacharya S, Li J, Yang W, and Kanai Y, J. Phys. Chem. A 128, 6072 (2024). [ DOI ] [ PubMed ] [ Google Scholar ] 34. Hillenbrand C, Li J, and Zhu T, J. Chem. Phys 162, 174117 (2025). [ DOI ] [ PubMed ] [ Google Scholar ] 35. Zhou R, Yao Y, Blum V, Ren X, and Kanai Y, J. Chem. Theory Comput 21, 291 (2025). [ DOI ] [ PubMed ] [ Google Scholar ] 36. Liu Z-F, ACS Nano 19, 5861 (2025). [ DOI ] [ PubMed ] [ Google Scholar ] 37. Jacquemin D, Duchemin I, and Blase X, Mol. Phys 114, 957 (2016). [ Google Scholar ] 38. Hung L, da Jornada FH, Souto-Casares J, Chelikowsky JR, Louie SG, and Öğüt S, Phys. Rev. B 94, 085125 (2016). [ Google Scholar ] 39. Förster A and Visscher L, J. Chem. Theory Comput 18, 6779 (2022). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 40. Venturella C, Hillenbrand C, Li J, and Zhu T, J. Chem. Theory Comput 20, 143 (2024). [ DOI ] [ PubMed ] [ Google Scholar ] 41. Venturella C, Li J, Hillenbrand C, Leyva Peralta X, Liu J, and Zhu T, Nat Comput Sci, 1 (2025). [ Google Scholar ] 42. Li J, Jin Y, Su NQ, and Yang W, J. Chem. Phys 156, 154101 (2022). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 43. Elliott JD, Colonna N, Marsili M, Marzari N, and Umari P, J. Chem. Theory Comput 15, 3710 (2019). [ DOI ] [ PubMed ] [ Google Scholar ] 44. Villalobos-Castro J, Knysh I, Jacquemin D, Duchemin I, and Blase X, J. Chem. Phys 159, 024116 (2023). [ DOI ] [ PubMed ] [ Google Scholar ] 45. Tölle J, Kitsaras M-P, and Loos P-F, J. Phys. Chem. Lett 16, 11134 (2025). [ DOI ] [ PubMed ] [ Google Scholar ] 46. Andersson K, Malmqvist P-Å, and Roos BO, J. Chem. Phys 96, 1218 (1992). [ Google Scholar ] 47. Siegbahn PEM, Almlöf J, Heiberg A, and Roos BO, J. Chem. Phys 74, 2384 (1981). [ Google Scholar ] 48. von Niessen W, Schirmer J, and Cederbaum LS, Comput. Phys. Rep 1, 57 (1984). [ Google Scholar ] 49. Koch H, Jensen HJA, Jo/rgensen P, and Helgaker T, J. Chem. Phys 93, 3345 (1990). [ Google Scholar ] 50. Loos P-F, Scemama A, Blondel A, Garniron Y, Caffarel M, and Jacquemin D, J. Chem. Theory Comput 14, 4360 (2018). [ DOI ] [ PubMed ] [ Google Scholar ] 51. Loos P-F, Lipparini F, Boggio-Pasqua M, Scemama A, and Jacquemin D, J. Chem. Theory Comput 16, 1711 (2020). [ DOI ] [ PubMed ] [ Google Scholar ] 52. Ripka SRPG, Blaizot J-P, and Ripka G, Quantum Theory of Finite Systems (MIT Press, 1986). [ Google Scholar ] 53. Ring P and Schuck P, The Nuclear Many-Body Problem, softcover repri edition ed. (Springer, Berlin Heidelberg, 2004). [ Google Scholar ] 54. van Aggelen H, Yang Y, and Yang W, Phys. Rev. A 88, 030501 (2013). [ Google Scholar ] 55. van Aggelen H, Yang Y, and Yang W, J. Chem. Phys 140, 18A511 (2014). [ Google Scholar ] 56. Peng D, van Aggelen H, Yang Y, and Yang W, J. Chem. Phys 140, 18A522 (2014). [ Google Scholar ] 57. ROWE DJ, Rev. Mod. Phys 40, 153 (1968). [ Google Scholar ] 58. Marie A, Romaniello P, Blase X, and Loos P-F, J. Chem. Phys 162, 134105 (2025). [ DOI ] [ PubMed ] [ Google Scholar ] 59. Bohm D and Pines D, Phys. Rev 82, 625 (1951). [ Google Scholar ] 60. Ren X, Rinke P, Joas C, and Scheffler M, J Mater Sci 47, 7447 (2012). [ Google Scholar ] 61. Shi R, Zhang M-Y, Lin P, He L, and Ren X, Computer Physics Communications 309, 109496 (2025). [ Google Scholar ] 62. Peng D, Steinmann SN, van Aggelen H, and Yang W, J. Chem. Phys 139, 104112 (2013). [ DOI ] [ PubMed ] [ Google Scholar ] 63. Scuseria GE, Henderson TM, and Bulik IW, J. Chem. Phys 139, 104113 (2013). [ DOI ] [ PubMed ] [ Google Scholar ] 64. Wang Y, Fang W-H, and Li Z, J. Chem. Phys 163, 174101 (2025). [ DOI ] [ PubMed ] [ Google Scholar ] 65. Yang Y, van Aggelen H, and Yang W, J. Chem. Phys 139, 224105 (2013). [ DOI ] [ PubMed ] [ Google Scholar ] 66. Berkelbach TC, J. Chem. Phys 149, 041103 (2018). [ DOI ] [ PubMed ] [ Google Scholar ] 67. Bannwarth C, Yu JK, Hohenstein EG, and Martínez TJ, J. Chem. Phys 153, 024110 (2020). [ DOI ] [ PubMed ] [ Google Scholar ] 68. Yu JK, Bannwarth C, Hohenstein EG, and Martínez TJ, J. Chem. Theory Comput 16, 5499 (2020). [ DOI ] [ PubMed ] [ Google Scholar ] 69. Yang Y, Peng D, Lu J, and Yang W, J. Chem. Phys 141, 124104 (2014). [ DOI ] [ PubMed ] [ Google Scholar ] 70. Yang Y, Peng D, Davidson ER, and Yang W, J. Phys. Chem. A 119, 4923 (2015). [ DOI ] [ PubMed ] [ Google Scholar ] 71. Yang Y, Davidson ER, and Yang W, Proc. Natl. Acad. Sci 113, E5098 (2016). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 72. Yang Y, Dominguez A, Zhang D, Lutsker V, Niehaus TA, Frauenheim T, and Yang W, J. Chem. Phys 146, 124104 (2017). [ DOI ] [ PubMed ] [ Google Scholar ] 73. Yu J, Li J, Zhu T, and Yang W, J. Chem. Phys 162, 094101 (2025). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 74. Yang Y, Shen L, Zhang D, and Yang W, J. Phys. Chem. Lett 7, 2407 (2016). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 75. Zhang D, Peng D, Zhang P, and Yang W, Phys. Chem. Chem. Phys 17, 1025 (2014). [ DOI ] [ PubMed ] [ Google Scholar ] 76. Li J, Yu J, Chen Z, and Yang W, J. Phys. Chem. A 127, 7811 (2023). [ DOI ] [ PubMed ] [ Google Scholar ] 77. Zhang D and Yang W, J. Chem. Phys 145, 144105 (2016). [ DOI ] [ PubMed ] [ Google Scholar ] 78. Li J, Jin Y, Yu J, Yang W, and Zhu T, J. Phys. Chem. Lett 15, 2757 (2024). [ DOI ] [ PubMed ] [ Google Scholar ] 79. Li J, Jin Y, Yu J, Yang W, and Zhu T, J. Chem. Theory Comput 20, 7979 (2024). [ Google Scholar ] 80. Chen Z, Zhang D, Jin Y, Yang Y, Su NQ, and Yang W, J. Phys. Chem. Lett 8, 4479 (2017). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 81. Li J, Chen Z, and Yang W, J. Phys. Chem. Lett 13, 894 (2022). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 82. Springer M, Aryasetiawan F, and Karlsson K, Phys. Rev. Lett 80, 2389 (1998). [ Google Scholar ] 83. Zhang D, Su NQ, and Yang W, J. Phys. Chem. Lett 8, 3223 (2017). [ DOI ] [ PubMed ] [ Google Scholar ] 84. Li J, Chen Z, and Yang W, J. Phys. Chem. Lett 12, 6203 (2021). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 85. Loos P-F and Romaniello P, J. Chem. Phys 156, 164101 (2022). [ DOI ] [ PubMed ] [ Google Scholar ] 86. Sun Q, Berkelbach TC, Blunt NS, Booth GH, Guo S, Li Z, Liu J, McClain JD, Sayfutyarova ER, Sharma S, Wouters S, and Chan GK-L, WIREs Comput. Mol. Sci 8, e1340 (2018). [ Google Scholar ] 87. Sun Q, Zhang X, Banerjee S, Bao P, Barbry M, Blunt NS, Bogdanov NA, Booth GH, Chen J, Cui Z-H, Eriksen JJ, Gao Y, Guo S, Hermann J, Hermes MR, Koh K, Koval P, Lehtola S, Li Z, Liu J, Mardirossian N, McClain JD, Motta M, Mussard B, Pham HQ, Pulkin A, Purwanto W, Robinson PJ, Ronca E, Sayfutyarova ER, Scheurer M, Schurkus HF, Smith JET, Sun C, Sun S-N, Upad-hyay S, Wagner LK, Wang X, White A, Whitfield JD, Williamson MJ, Wouters S, Yang J, Yu JM, Zhu T, Berkelbach TC, Sharma S, Sokolov AY, and Chan GK-L, J. Chem. Phys 153, 024109 (2020). [ DOI ] [ PubMed ] [ Google Scholar ] 88. Marie A, Romaniello P, and Loos P-F, Phys. Rev. B 110, 115155 (2024). [ Google Scholar ] 89. Langreth DC and Perdew JP, Phys. Rev. B 15, 2884 (1977). [ Google Scholar ] 90. Yang Y, van Aggelen H, Steinmann SN, Peng D, and Yang W, J. Chem. Phys 139, 174110 (2013). [ DOI ] [ PubMed ] [ Google Scholar ] 91. Yang W, Mori-Sánchez P, and Cohen AJ, J. Chem. Phys 139, 104114 (2013). [ DOI ] [ PubMed ] [ Google Scholar ] 92. Williams-Young D, Egidi F, and Li X, J. Chem. Theory Comput 12, 5379 (2016). [ DOI ] [ PubMed ] [ Google Scholar ] 93. Dyall KG, J. Chem. Phys 106, 9618 (1997). [ Google Scholar ] 94. Kutzelnigg W and Liu W, J. Chem. Phys 123, 241102 (2005). [ DOI ] [ PubMed ] [ Google Scholar ] 95. Iliaš M and Saue T, J. Chem. Phys 126, 064102 (2007). [ DOI ] [ PubMed ] [ Google Scholar ] 96. Liu W and Peng D, J. Chem. Phys 131, 031104 (2009). [ DOI ] [ PubMed ] [ Google Scholar ] 97. Heß BA, Marian CM, Wahlgren U, and Gropen O, Chem. Phys. Lett 251, 365 (1996). [ Google Scholar ] 98. Liu J and Cheng L, J. Chem. Phys 148, 144108 (2018). [ DOI ] [ PubMed ] [ Google Scholar ] 99. Zhang C and Cheng L, J. Phys. Chem. A 126, 4537 (2022). [ DOI ] [ PubMed ] [ Google Scholar ] 100. Yang W and Wu Q, Phys. Rev. Lett 89, 143002 (2002). [ DOI ] [ PubMed ] [ Google Scholar ] 101. Jin Y, Zhang D, Chen Z, Su NQ, and Yang W, J. Phys. Chem. Lett 8, 4746 (2017). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 102. Martin RL, J. Chem. Phys 118, 4775 (2003). [ Google Scholar ] 103. Jin Y, Su NQ, Chen Z, and Yang W, Faraday Discuss. 224, 9 (2020). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 104. Harris CR, Millman KJ, van der Walt SJ, Gommers R, Virtanen P, Cournapeau D, Wieser E, Taylor J, Berg S, Smith NJ, Kern R, Picus M, Hoyer S, van Kerkwijk MH, Brett M, Haldane A, del Río JF, Wiebe M, Peterson P, Gérard-Marchant P, Sheppard K, Reddy T, Weckesser W, Abbasi H, Gohlke C, and Oliphant TE, Nature 585, 357 (2020). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 105. Virtanen P, Gommers R, Oliphant TE, Haberland M, Reddy T, Cournapeau D, Burovski E, Peterson P, Weckesser W, Bright J, van der Walt SJ, Brett M, Wilson J, Millman KJ, Mayorov N, Nelson ARJ, Jones E, Kern R, Larson E, Carey CJ, Polat I, Feng Y, Moore EW, VanderPlas J, Laxalde D, Perktold J, Cimrman R, Henriksen I, Quintero EA, Harris CR, Archibald AM, Ribeiro AH, Pedregosa F, and van Mulbregt P, Nat Methods 17, 261 (2020). [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 106. Gaussian Cube File Format, https://paulbourke.net/dataformats/cube/ . 107. Kossoski F, Boggio-Pasqua M, Loos P-F, and Jacquemin D, J. Chem. Theory Comput 20, 5655 (2024). [ DOI ] [ PubMed ] [ Google Scholar ] 108. Silva-Junior MR, Schreiber M, Sauer SPA, and Thiel W, J. Chem. Phys 129, 104103 (2008). [ DOI ] [ PubMed ] [ Google Scholar ] 109. Bernard YA, Shao Y, and Krylov AI, J. Chem. Phys 136, 204103 (2012). [ DOI ] [ PubMed ] [ Google Scholar ] 110. Loos P-F, Boggio-Pasqua M, Scemama A, Caffarel M, and Jacquemin D, J. Chem. Theory Comput 15, 1939 (2019). [ DOI ] [ PubMed ] [ Google Scholar ] 111. Roos BO, Lindh R, Malmqvist P-Å, Veryazov V, and Widmark P-O, J. Phys. Chem. A 108, 2851 (2004). [ Google Scholar ] 112. Cohen AJ, Mori-Sánchez P, and Yang W, Science 321, 792 (2008). [ DOI ] [ PubMed ] [ Google Scholar ] 113. Mori-Sánchez P, Cohen AJ, and Yang W, Phys. Rev. A 85, 042507 (2012). [ Google Scholar ] 114. Mori-Sánchez P, Cohen AJ, and Yang W, Phys. Rev. Lett 102, 066403 (2009). [ DOI ] [ PubMed ] [ Google Scholar ] 115. Stein T, Kronik L, and Baer R, J. Am. Chem. Soc 131, 2818 (2009). [ DOI ] [ PubMed ] [ Google Scholar ] 116. Duke data research repository, 10.7924/r4c82k88g. [ DOI ] [ Google Scholar ] Associated Data This section collects any data citations, data availability statements, or supplementary materials included in this article. Supplementary Materials Supporting Information NIHMS2147869-supplement-Supporting_Information.pdf (173.7KB, pdf) Data Availability Statement Data and scripts pertaining to this work have been archived in the Duke Research Data Repository 116 . ACTIONS View on publisher site PDF (1.1 MB) Cite Collections Permalink PERMALINK Copy RESOURCES Similar articles Cited by other articles Links to NCBI Databases Cite Copy Download .nbib .nbib Format: AMA APA MLA NLM Add to Collections Create a new collection Add to an existing collection Name your collection * Choose a collection Unable to load your collection due to an error Please try again Add Cancel Follow NCBI NCBI on X (formerly known as Twitter) NCBI on Facebook NCBI on LinkedIn NCBI on GitHub NCBI RSS feed Connect with NLM NLM on X (formerly known as Twitter) NLM on Facebook NLM on YouTube National Library of Medicine 8600 Rockville Pike Bethesda, MD 20894 Web Policies FOIA HHS Vulnerability Disclosure Help Accessibility Careers NLM NIH HHS USA.gov Back to Top