where Ψ 0 and E 0 are the ground state wavefunction and energy, respectively, and n is the index of the excited states. Indeed, since ∆E HL = E 1 -E 0 is commonly much smaller than the energy gap of higher excited states, the first term of the sum in Eq. 1 provides a first-order approximation to the infinite series and, hence, there could exist an inversely proportional relationship between ∆E HL and α, i.e., α ∝ (∆E HL ) -1 . Equation 1 can only be analytically evaluated for simple model systems (such as the hydrogen atom or a quantum Drude oscillator). For many-electron systems, the sum can only be evaluated numerically and requires including bound-bound and bound-continuum transi-tion dipoles . Modeling atoms or larger coarse-grained fragments with QDOs and solving the dipole-dipole screening equations is known to be an effective method to predict polarizability, and it is also the basis for the Many-Body Dispersion (MBD) method . Since the response properties of all atoms and molecules can be represented by QDOs by carefully setting the three parameters {charge q, frequency ω, mass µ} of the model, the analysis of the polarizability of the QDO Hamiltonian should generally be transferable to any systems. Therefore, we devote some attention to this model.