Communication : Configuration interaction singles has a large systematic bias against charge-transfer states

We show that standard configuration interaction singles (CIS) has a systematic bias against chargetransfer (CT) states, wherein the computed vertical excitation energies for CT states are disproportionately too high (by >1 eV) as compared with non-CT states. We demonstrate this bias empirically for a set of chemical problems with both interand intra-molecular electron transfer, and then, for a small analytical model, we prove that this large difference in accuracy stems from the massive changes in electronic structure that must accompany long-range charge transfer. Thus far, the conclusion from this research is that, even in the context of wave function theory, CIS alone is insufficient for offering a balanced description of excited state surfaces (both CT and non-CT) and explicit electron-electron correlation must be included additionally (e.g., via CIS(D)) for charge-transfer applications. © 2011 American Institute of Physics. [doi:10.1063/1.3627152]


I. INTRODUCTION
The simplest means for computing molecular excited states is configuration interaction singles (CIS), where the wave function is given the form | CI S = ia t a i | a i .Following standard protocol, in this communication, we denote occupied orbitals ij and virtual orbitals ab.For CIS theory, the one electron difference density matrix can be be written as the sum of electron attachment ("the particle") ρ att ab = i t a i t b i and electron detachment ("the hole") ρ det ij = a t a i t a j . 1 As a tool for computing excited state energies, CIS has several advantages: first, it is computationally cheap; second it recovers both non-charge-transfer (CT) and CT states, including the correct −1/r asymptotic behavior of CT states that comes about because of the Coulombic attraction between electron attachments and detachments. 2The disadvantage of CIS is that, because the method makes no attempt to allow for electron-electron correlation in the excited state, 3 it often does not recover accurate vertical excitation energies.
In recent years, time-dependent density functional theory 4 (TD-DFT) has emerged as a convenient alternative to CIS that usually computes excitation energies better than CIS and that, like CIS, can also be applied to large systems of chemical interest.Excelling at short-range interactions, TD-DFT does remarkably well at computing excited states without charge-transfer character.3][14][15][16][17][18][19][20][21][22][23][24] These so-called longa) Electronic mail: subotnik@sas.upenn.edu.ranged corrected (LRC) DFT functionals have a slight flavor of CIS and attempt to find a meaningful balance between DFT functionals that work for ground-state properties and the correct asymptotic form for CT states.Progress is being made for these LRC-TD-DFT functionals, but as with all DFT development, there is no systematic way to improve accuracy.
In this article, we return to basic, systematically improvable wave-function methods for excited states.We ask a very simple question: knowing that CIS recovers both non-CT and CT states (including the correct −1/r asymptotic form for CT states), can we say that the relative CIS energies are unbiased between CT and non-CT states?In other words, because CIS is a simple wave function approach without any approximate DFT functional, might it be true that the CIS error will not discriminate for or against CT states (vs.non-CT states)?In the remainder of this article, we will answer these questions emphatically in the negative.
We will show below that usually, compared with non-CT states, CIS systematically gives worse excitation energies for CT states, and that the error is not small; they are too high usually by more than 1 eV.While Tozer and co-workers 25,26 have argued that TD-DFT errors can be correlated with a measure of charge-transfer (though this is not always true 27 ), we now make surprisingly similar claims for CIS itself.In Sec.II, we will give empirical evidence of this assertion for two meaningful chemical systems-with interand intra-molecular electron transfer-and in Sec.III, we will provide analytical evidence to support our claim.Throughout our discussion, we will use a very simple form of perturbation theory, CIS(D), 28,29 to compute a correction to CIS, E = E CI S(D) − E CI S .By calculating CIS and CIS(D) energies at hundreds of different nuclear geometries and for many electronic states each, we will see clear trends in E versus a CT coordinate that suggest a large failure of CIS when computing excited state surfaces.For small systems, we will verify these results with equation of motion coupled-cluster singles doubles (EOM-CCSD). 30

A. Water clusters
For our first set of model problems, we calculate excitation energies for water dimers, trimers, and tetramers.For each set of molecules, we consider 500 different geometries, chosen from a 5 ps long ground-state trajectory from the TINKER molecular dynamics program 31 at room temperature.For each cluster, we compute all excited states up to ∼18 eV, which corresponds, respectively, to 15 excited states (dimer), 25 states (trimer), and 35 states (tetramer).We compute excitation energies using standard CIS and also perturbationcorrected CIS(D), in a 6-31g** basis.
For these clusters with multiple charge centers, we choose to "measure" charge-transfer character according to the Coulombic interaction between electron attachment and detachment, On the one hand, when ζ is large (and negative), the attachment and detachment are nearby in space and there is little charge transfer character.On the other hand, when ζ is small, the attachment and detachment are far apart in space and there is a great deal of charge transfer character.Note that ζ is not a perfect measure of charge-transfer character, however, as the value depends also on size, shape, and symmetry of the electronic densities; nevertheless, ζ roughly measures CT character.In Fig. 1, we make a scatter plot of E = E CI S(D) − E CI S and the magnitude of ζ .The data is unambiguous: we see that the correlation energy systematically gets more negative for excited states with more charge-transfer character (smaller ζ ).Moreover, the scale here is large; note that the y-axis in Fig. 1 is in units of eV.

B. PYCM
Our conclusions for water clusters are confirmed and amplified when we now consider the chargetransfer states of the molecule 2-(4-(propan-2-ylidene) cyclohexylidene)malononitrile shown in the inset of Fig. 2, abbreviated as PYCM.PYCM has been studied previously experimentally 32 and theoretically 33 in the context of radiative decay as modulated by solvent.The lifetime of the PYCM CT state depends strongly and nonlinearly on the ambient solvent, polar vs. nonpolar, suggesting mixing between adiabatic states.As such, a proper treatment of this molecule requires a balanced treatment of both CT and non-CT excited states.
In Fig. 2, we make a scatter plot of the correction energy ( E) for the first 10 excited states of PYCM.This includes all vertical excitations up to (roughly) 8.5 eV, and included in this set is at least one long-range charge-transfer state where the di-methyl alkene functional group donates an electron to the dicyano functional group.For this case of only two charge centers, we choose to "measure" charge transfer by the magnitude of the relative dipole moment of the CIS excited state, | μ |.Here, the x-component of μ is CI S |μ x | CI S = ab t a i t b i x ab − ij t a i t a j x ij , and similarly for y and z.As for the case of water clusters, we find that the CIS(D) correction energy is far larger and more negative for the CT states (by 1-2 eV).In fact, here we see something stronger: E appears to be linearly related to the dipole moment (for large dipole moments).

III. ANALYTICAL TREATMENT: HE DIMER
Having demonstrated empirically that, on average, CIS excitation energies for CT states are too large relative to non-CT states, we now explore why this is so for a small model problem.Consider a helium dimer, where the nuclei are separated by a distance r, and we choose a minimal basis with two atomic s orbitals per center (for a total of 4 orbitals).For convenience, we insert an external charge of 0.8 a.u.far away from the helium to break the spatial symmetry of the problem.
Without loss of generality, we assume that the helium atoms are so far apart that electrons on the two centers do not interact.In this case, we can write the Hartree-Fock ground state as | H F = |L LR R .Here, we denote an occupied molecular orbital on the left by "L" and on the right by "R"; we let L * and R * be virtual orbitals on the left and right helium centers.The Hartree-Fock orbitals must diagonalize the Fock operator (F ), and ignoring all matrix operators between different centers, this implies, e.g., where "O" represents all one-electron operators and we are using physics notation for the antisymmetrized two-electron integrals, pq||uv .
Starting from the ground state | H F , the four singlet CIS states will be approximately (2) (3) (4) It is now straightforward to see that the CT states (| (3) CIS and | (4) CIS ) have much stronger couplings to the doubles excitations than do the non-CT states, (| (1) CIS and | (2) CIS ).For instance, on the one hand, for non-CT states, the couplings to double excitations look like Here, we have ignored all couplings between left and right centers and applied Eq. ( 2).In Eq. ( 8), the total expression will be small under the assumption that, in some average sense, the charge densities and Coulombic potentials of orbitals L and L * are similar.This assumption will never be exactly true, but at the very least, if these charge densities are at all similar, we do expect that the two terms in Eq. ( 8) should be of the same sign, so there will be partial cancellation.
On the other hand, however, for CT states, we find singles-doubles couplings of the form, e.g., To derive Eq. ( 10), we ignore all couplings between left and right centers and invoke Eq. ( 2).We now find only one Coulombic term in the matrix element, rather than two partially cancelling terms.Assuming that orbitals L and L * do not have different symmetries so that all these integrals do not vanish (as they would, say, for H 2 dimers), the end result of Eq. ( 10) is that singles-double couplings are very large for CT states, and CT states cannot be treated adequately using CIS theory alone.The reasoning here is very clear: when an α electron in orbital R moves to the excited orbital L * , it is essential to allow the β electron in orbital L to relax partially into orbital L * .Without this relaxation, CIS vastly overestimates the vertical energy of this excitation.Similar orbital relaxation effects should exist for most molecular CT states provided there are no unusual symmetry effects.All of the assertions above are realized in Fig. 3.In (a), we plot the CIS energies for the four CIS excited states: the lower two are non-CT and the upper two are CT.In (b), we plot the CIS(D) energy correction for each excited state.As would be expected from the arguments above, the correction is much stronger (< −1eV) for the CT states.In (c), (d), we show that these corrections do not depend on the correlation method: they are also found with EOM-CCSD, 30 as shown for excited states 1 and 4. Finally, note that, in Fig. 3(b), E depends linearly on r for small r (and also μ , not shown).This is exactly the behavior we saw above for PYCM, now as captured by a very small and tractable model problem.Thus, we have a simple explanation for why CIS overestimates CT vertical excitation energies.

IV. FUTURE DIRECTIONS
We have argued that CIS provides an inadequate and unbalanced view of charge-transfer states: Even though CIS predicts the correct −1/r asymptotic energy as charges separate, the method grossly overestimates the vertical excitation energies of charge-transfer states (by >1 eV) compared to noncharge-transfer states.In some sense, these results may not be surprising because quantum chemists have long known that, for CIS accuracy, the electronic correlations in the excited state must resemble those in the ground state-which is clearly false for CT states. 34What is surprising, however, is simply the scale and systematic nature of the CIS error in all the figures above, which partially resemble the failures of TD-DFT. 25,26  the future, it will be essential to explore the size and extent of this CIS error for benchmarking purposes.Having done so, we anticipate at least two options for computing CT states.(1) If we seek wave functions for describing charge transfer, CIS must be supplemented with fast correlation methods to account for double or higher order excitations In order to prove that these results are not an artifact of the CIS(D) correction, in (c) and (d), we plot the energies of states S 1 and S 4 for EOM-CCSD and our results are nearly unchanged.Note how much more inaccurate CIS energies are for CT states S 3 , S 4 compared with non-CT states S 1 , S 2 .For CT states, μ ∝ r, so that for small r, (b) reflects the same physics as in Fig. 2.
and new techniques are likely needed.(2) As an alternative to wave function techniques, if we seek to compute excited CT states with TD-DFT, it will be essential to check the behavior of long-range corrected TD-DFT functionals given our new information about CIS's failures.After all, LRC-DFT functionals recover the −1/r asymptotic energy of CT states by introducing exact HF exchange, similar to CIS.Perhaps we will find that because standard TD-DFT underestimates CT vertical excitation energies and CIS overestimates CT vertical excitation energies, LRC-TD-DFT will recover the correct answers. 3,35 r perhaps not.The problem of computing CT states in a fast, reliable, and balanced manner is a key roadbloack for future modeling of excited state charge transfer dynamics.

FIG. 2 .
FIG. 2. (a) A scatter plot of E versus the magnitude of the relative dipole moment, | μ | for the molecule shown in (b).The formal name for this molecule is 2-(4-(propan-2-ylidene)cyclohexylidene)malononitrile, abbreviated here as PYCM.For each of 500 nuclear geometries, we plot the first 10 electronic states.

FIG. 3 .
FIG. 3. (a)CIS and (b) E curves for the He 2 dimer separated by a distance r in a VDZ basis (with an external charge far away to break symmetry).In order to prove that these results are not an artifact of the CIS(D) correction, in (c) and (d), we plot the energies of states S 1 and S 4 for EOM-CCSD and our results are nearly unchanged.Note how much more inaccurate CIS energies are for CT states S 3 , S 4 compared with non-CT states S 1 , S 2 .For CT states, μ ∝ r, so that for small r, (b) reflects the same physics as in Fig.2.