A corrected-linear response formalism for the calculation of electronic excitation energies of solvated molecules with the CCSD-PCM method

Computational and Theoretical Chemistry xxx (2014) xxx–xxx

Contents lists available at ScienceDirect

Computational and Theoretical Chemistry journal homepage: www.elsevier.com/locate/comptc

A corrected-linear response formalism for the calculation of electronic excitation energies of solvated molecules with the CCSD-PCM method Marco Caricato ⇑ Gaussian, Inc., 340 Quinnipiac St. Bldg. 40, Wallingford, CT 06492, USA

a r t i c l e

i n f o

Article history: Received 23 December 2013 Received in revised form 4 February 2014 Accepted 4 February 2014 Available online xxxx Keywords: CCSD-PCM Corrected linear response Excitation energy Solution

a b s t r a c t In this work, we present an efficient approximation to compute excitation energies in solution when coupled cluster (CC) methods are combined with a polarizable solvation model. Two formalisms exist to compute excited state energies with polarizable solvation models: state-specific (SS) and linear-response (LR). The former more accurately describes the solute–solvent polarization in the excited state, but is computationally intensive. The LR formalism is efficient, but lacks proper relaxation effects. An approximate method, called corrected-LR (cLR), was originally formulated in the context of time-dependent density functional theory and the polarizable continuum model of solvation (PCM), and was shown to be able to recover most of the relaxation contributions of the SS formalism at a cost similar to LR. We have expanded the cLR idea to CC theory, and introduced an extra approximation that further reduces the computational effort with negligible loss of accuracy. The test cases reported in this contribution clearly show that the cLR-CC-PCM method is able to estimate transition energies in very close agreement with the SS formalism at a cost that is similar (in fact, slightly smaller) than the LR formalism. The average SS–LR difference is of the order of 0.10–0.20 eV for nonequilibrium calculations, and 0.30–0.55 eV for equilibrium calculations. The SS–cLR average difference is, on the other hand, 0.01–0.02 eV for nonequilibrium calculations, and 0.05–0.15 eV for equilibrium calculations. Therefore, the cLR approach is a promising alternative for computing excited state energies in solution with computationally intensive CC methods. Ó 2014 Elsevier B.V. All rights reserved.

1. Introduction Electronic excited states are central for research in many areas, from materials to energy science, from photoactivated catalysis to biochemistry. Most of these experiments are conducted in condensed phase, often in solution. Therefore, a realistic simulation of the photochemistry of chromophores surrounded by a solvent is of paramount importance to interpret experimental findings and predict the behavior of novel chromophoric species. Such simulations are, however, particularly challenging because they combine the difficulty of describing excited states with that of describing the effect of a large number of surrounding molecules (the solvent) on the excitation process. Accurate methods as those belonging to the coupled cluster (CC) family [1–6] are often necessary to obtain reliable photochemical information on a particular system, but accounting for the environmental effects would require a computational effort that is, in many cases, too large. Efficient computational strategies are thus necessary to preserve the

⇑ Tel.: +1 2032842501. E-mail address: [email protected]

accuracy of these methods while reducing the computational cost to a manageable level. To achieve this goal, a multilayer approach can be applied where the chromophore is treated at quantum mechanical (QM) level while the solvent is described by a classical polarizable field. The interaction between the two levels of theory is achieved by introducing the effect of the classical field in the QM Hamiltonian as a one-electron operator, which polarizes the wave function of the chromophore. These solvation models can use an atomistic description of the solvent (explicit models) or a continuum dielectric representation (implicit models). The former allow the description of specific solute–solvent interactions (e.g., hydrogen bonds) at the price of requiring multiple QM calculations to reach a statistical averaging over the possible solvent configurations. The latter implicitly include this averaging by using the macroscopic permittivity of the solvent (related to the static dielectric constant e), but are not able to describe specific interactions. Excited state solvation also introduces two complications: the time scale of the solvent response [7–10], and the choice between state-specific (SS) and linear-response (LR) formalisms [11,12]. The first issue distinguishes the calculation of vertical excitations, where the solvent molecules do not have enough time to equili-

http://dx.doi.org/10.1016/j.comptc.2014.02.001 2210-271X/Ó 2014 Elsevier B.V. All rights reserved.

Please cite this article in press as: M. Caricato, Comput. Theoret. Chem. (2014), http://dx.doi.org/10.1016/j.comptc.2014.02.001

2

M. Caricato / Computational and Theoretical Chemistry xxx (2014) xxx–xxx

brate with the new electronic distribution of the solute, from the calculation of excited state stationary points (e.g., minima of the excited potential energy surface of the solute), where both solute and solvent relax simultaneously. Vertical excitation energies are thus computed in a nonequilibrium solvation regime, while relaxed quantities are computed in an equilibrium solvation regime as for the ground state. In continuum solvation models, nonequilibrium is obtained by dividing the solvent polarization in a fast (or dynamic) contribution (related to the optical dielectric constant e1 ), and a slow (or inertial) contribution, which is kept frozen in equilibrium with the initial state of the solute. The difference between SS and LR formalisms resides in how the solvent polarization in the excited state is computed [11,12]. For SS, the excited state wave function is considered explicitly, and the solvent responds to the excited state electronic density. For LR, the solvent interaction term is evaluated as the linear response of the ground state solvent operator to an external perturbation (i.e., the solvent responds to a transition electronic density). These two formulations provide different results, and also require a different computational effort. SS contains a more complete account of the mutual solute–solvent polarization in the excited state, but is computationally more demanding since it involves an iterative cycle over these mutual polarization effects, and different electronic states must be considered separately. LR, on the other hand, lacks such polarization but the solvent effect can be computed for multiple states at once. Therefore, the LR formalism is comparable in cost to a gas phase calculation. For both formalisms, equilibrium and nonequilibrium regimes can be described. The goal of this work is to reduce the gap between these two formalisms when solvation models are used with CC methods, by resorting to an approach called corrected-LR (cLR) [13]. The cLR approach was developed in the context of time-dependent density functional theory (TDDFT), and is based on ‘‘correcting’’ the linear response transition energy computed in the frozen ground state reaction field with a state-specific energy term. cLR can also be seen as the first iteration of the SS formalism. This approach was shown to reduce considerably the difference between SS and LR energies while maintaining a computational cost similar to LR [13]. In this contribution, we consider the cLR approach in the framework of CC theory, and compare the cLR-CC and LR-CC results against full SS-CC calculations. Two formulations of the cLR corrections are considered, one that includes the relaxation of the CC T amplitudes and one that neglects this term. The latter is computationally very efficient since it allows the calculation of the cLR energy term directly for multiple excited states, thus resulting in a computational effort that is even slightly smaller than that of the LR formalism. We stress that, although the following discussion focus on implicit solvation models (in particular the polarizable continuum model, PCM [10]), the cLR approach can be also applied to explicit polarizable models. The article is organized as follows. The theory and computational considerations are presented in Section 2. Numerical results on a set of test compounds are reported in Section 3. A discussion of these findings and conclusive remarks are presented in Section 4. 2. Theory This section presents a theoretical description of the cLR approach [13] for the calculation of excited state energies of molecules in solution applied to the combination of CC methods for the electronic structure of the solute and PCM for the solvent effect. A detailed description of CC methods [1–6] and PCM [10,14] can be found elsewhere. We also refer the reader to the available literature on how to combine CC methods with continuum solvation models [15–30].

The calculation of electronic excitation energies of molecules in condensed phase, where the latter is described with a polarizable solvation model, can be performed using two formalisms, namely state-specific (SS) and linear-response (LR). The difference stems from the Hartree partition of the solute–solvent state that allows to treat the latter classically. A comparison between these two formalisms for exact states was performed by Cammi et al. [11] and Corni et al. [12], and by this author in the context of CC methods [31]. The reader is referred to these articles for more details. In summary, in the SS formalism the solute–solvent mutual polarization is treated at the same footing in the ground and in the excited states, whereas in the LR only a response of the ground state reaction field due to the excitation process is computed. Therefore, the SS provides a better description of the solvent polarization due to the excited state solute electronic density. However, this approach is computationally more demanding since it requires a self-consistent solution of the excited state equations in the presence of the solvent to account for this mutual polarization, and this can be only performed separately for each state. Furthermore, the ground and excited states are no longer orthogonal, which renders the evaluation of transition moments more cumbersome. On the other hand, the LR formalism is a much more efficient computational paradigm, since the cost of the explicit solvent contribution is negligible, multiple states can be considered at once, and transition properties are well defined. In fact, the computational effort for a LR calculation in solution is comparable with that in gas phase. Nonequilibrium solvation formula can be derived in both formalisms for the calculation of vertical excitation energies [11,12,29,31]. From a numerical standpoint, these two formalisms provide different values of vertical excitation energies that can amount to 0.1–0.2 eV already in small systems, although, relative solvatochromic shifts between different solvents are usually in good agreement [13,31]. This difference can become of the order of 0.5–0.6 eV in the equilibrium solvation regime. Additionally, the sign of the difference between the SS and LR excitation energies is not always the same (say, SS transition energies are not always larger than LR ones). It was shown that the relative magnitude of the excitation energies is related to the relative magnitude of the dipole difference between states (representing SS) and the transition dipole (representing LR) [12]. 2.1. corrected-LR (cLR) In order to bridge the gap between LR and SS, a corrected-LR (cLR) approach was proposed, and implemented in the context of time-dependent density functional theory (TDDFT) [13]. The cLR transition energies between the ground and the Kth state are computed as:

1 2

xK ¼ x0 þ VD  Q D 1 2

D D xneq K ¼ x0 þ V  Q dyn

ð1Þ ð2Þ

where x0 is the excitation energy computed in the presence of the ground state reaction field, V is the solute electrostatic potential on the cavity surface, and Q are the corresponding charges. The overbar indicates that the operators are contracted with an appropriate density, and the superscript D indicates the density difference between ground and excited states. xK in Eq. (1) is the transition energy to the Kth state in the equilibrium solvation regime (which assumes that the solvent has had enough time to reorganize around the new solute charge distribution), while xneq in Eq. (2) is the verK tical excitation energy. The dynamic charges Q dyn provide the response of the solvent electrons to changes in the solute electronic distribution, and are computed replacing the static dielectric

Please cite this article in press as: M. Caricato, Comput. Theoret. Chem. (2014), http://dx.doi.org/10.1016/j.comptc.2014.02.001

M. Caricato / Computational and Theoretical Chemistry xxx (2014) xxx–xxx

constant e with the optical dielectric constant e1 in the PCM equations [13]. x0 is the same for equilibrium and nonequilibrium regimes since it is computed with the frozen ground state charges. For DFT, the ground state charges Q 0 are computed during the solution of the self-consistent field (SCF) equations (also called self-consistent reaction field equations, SCRF). For a derivation of Eqs. (1) and (2), see Ref. [13]. In that article, the D density is computed by considering the orbital response to the excitation of interest, which is obtained as the solution of the Kohn–Sham Z-vector equations (relaxed density). Therefore, the energy correction due to the Q D is state-specific. The cLR approach can be also seen as a SS calculation terminated at the first iteration. However, contrary to a full SS iterative procedure, transition properties are well defined since the transition energy is computed in a frozen reaction field, and the excitation amplitudes are still orthogonal to the ground state. Therefore, cLR improves upon the LR formalisms while maintaining its computational advantages. A downside of the cLR approach is that the excitation energy is no longer variationally optimized with respect to the excitation amplitudes, which makes the derivation of analytic gradients more cumbersome.

The main scope of this work is to extend the the cLR approach to CC methods in solution. Although the following formulas are applicable to any level of truncation in the CC expansion, the actual implementation is limited to singles and doubles excitations (CCSD) [1,28,29,32]. An important difference between CC- and SCF-PCM methods is in the definition of the ground state response field. In CC, the ground state charges are separated into two contributions: the reference charges, Q 0 , which are computed with the reference (SCF) one-particle density matrix (1PDM), and the correlation charges Q N , which are computed with the CC reduced 1PDM [22,25]:

c0pq ¼ hU0 jð1 þ KÞeT fpy qgeT jU0 i

and Q N þ Q LN for the left-hand diagonalization. The transition charges Q RN and Q LN are computed with the transition densities: T T c0K pq ¼ hU0 jð1 þ KÞ½e fp; qge ; RK jU0 i K0 T T cpq ¼ hU0 jLK e fp; qge jU0 i

ð3Þ

where U0 is the reference wave function, T and K are the CC excitation and de-excitation operators [5,6], respectively, and p; q refer to generic molecular orbitals (MOs). The charges Q 0 do not change during the solution of the CC equations since the orbitals are not reoptimized as in typical gas phase calculations, and they will be implicitly included in the Fock operator. The charges Q N , on the other hand, enter the CC equations explicitly as a one-electron operator, and couple the T and K equations (contrary to gas phase calculations). Approximations have been proposed to reduce the computational cost of the CC-PCM method [26–29], although, these will not be considered here. For excited states, SS and LR formalisms for CCSD-PCM are compared in Ref. [31]. In summary, the main difference resides in what charges are used in the calculation. In the SS formalism, the excited state charges Q KN are computed with the Kth state reduced 1PDM:

cKpq ¼ hU0 jLK ½eT fpy qgeT ; RK jU0 i þ hU0 jð1 þ KK ÞeT fpy qgeT jU0 i ð4Þ where LK and RK are the Kth left- and right-hand eigenvectors of the CC similarity transformed Hamiltonian, respectively, and KK is the excited state equivalent of the ground state CC de-excitation operator. With this charges, the excited state free energy Lagrangian is minimized with respect to changes in the four sets of amplitudes: T; LK ; RK , and KK . The solvent energy term depends quadratically on cK , thus it couples all the amplitudes equations similarly to the case of the ground state. In the LR formalism, the charges that enter the excited state equations are the sum of the ground state and transition charges: Q N þ Q RN for the right-hand diagonalization,

ð5Þ

respectively. Since the right-hand charges do not depend on the LK amplitudes (and similarly for the left-hand charges), the amplitudes equations are decoupled as in a gas phase diagonalization. As for TDDFT, the scope of the cLR approach is to preserve the computational efficiency of the LR formalism while correcting the final transition energy with a state-specific contribution. The cLR transition energies for CCSD-PCM can be computed with the same expressions in Eqs. (1) and (2), where the difference with DFT is in how the x0 and Q D terms are evaluated. x0 is obtained by diagonalization of the similarity transformed Hamiltonian in the presence of the frozen reaction field of the ground state, which is now composed by the two sets of charges Q 0 and Q N . This is again in common for equilibrium and non-equilibrium regimes. The D charges Q DN are computed with a density that is the difference between the excited state and ground state reduced densities, Eqs. (3) and (4):

cD ¼ cK  c0

2.2. corrected-LR (cLR) for CC methods

3

ð6Þ

Dynamic or equilibrium D charges are obtained by using e1 or e, respectively. An important difference between this formulation and that used in Ref. [13] is that here no orbital relaxation is considered. This is for consistency with the SS formulation of CCSD-PCM [28,29]. The cLR-CCSD-PCM method maintains the amplitudes equations decoupled as in the LR formalism, but introduces a SS energy correction at the end of the excited state calculation. Additionally, transition moments can be computed since the frozen reaction field in the excited state calculation does not break the orthogonality condition between states. Transition properties, like oscillator and rotatory strengths, can be then computed by using these transition moments and the corrected transition energies. The cost of the cLR calculation can be further reduced if the term that depends on KK in Eq. (4) is neglected. This corresponds to neglecting the T amplitudes response to the excitation perturbation. The approximate density has the form [2]:

c~Kpq ¼ hU0 jLK eT fpy qgeT RK jU0 i

ð7Þ

that would replace cK in Eq. (6). A similar approximation was shown to be very effective in the calculation of transition properties with LR-CCSD-PCM [31]. In this context, this approximation allows to compute cLR corrections for multiple states in a single calculation without adding any significant cost. In fact, all necessary sets of amplitudes are available at the end of the similarity transformed ~K in Eq. (7) Hamiltonian diagonalization, and the evaluation of c scales as OðN 5 Þ compared to the OðN 6 Þ scaling of CCSD. The validity of this approximation is also tested in Section 3. A limitation of the cLR approach, as mentioned in the previous section, is that analytic gradients are more complex to evaluate compared to those in iterative procedures that minimize the excited state energy with respect to the LK and RK amplitudes (as in both SS and LR formalisms). However, this is a minor issue because energy gradients are meaningless in a nonequilibrium calculation, while in an equilibrium calculation one is usually interested in only one state (e.g., in an excited state geometry optimization) and the overhead of the full SS formalism is less accentuated. The appeal of the cLR approach in the equilibrium regime is to provide preliminary information, for instance to determine the ordering of the states before selecting the proper one whose geometry can be optimized with the full SS method. The cLR-CCSD-PCM approach was implemented in a development version of the Gaussian suite of programs [33].

Please cite this article in press as: M. Caricato, Comput. Theoret. Chem. (2014), http://dx.doi.org/10.1016/j.comptc.2014.02.001

4

M. Caricato / Computational and Theoretical Chemistry xxx (2014) xxx–xxx

In this section, the excitation energies for four small molecules are computed in polar solvents within the equilibrium and nonequilibrium regimes. These test systems are methylencyclopropene (MCP) in methanol solution, acrolein, acetone, and water in water solution. These compounds were used in Ref. [31] for a comparison between the SS and LR formalisms. Here, the comparison is extended to the x0 and cLR transition energies. For the latter, the ~K in Eqs. choice of the complete density cK or T-frozen density c (4) and (7), respectively, is also tested. The results obtained with c~K are indicated in the following as ‘‘cLR0’’. The geometries are taken from Ref. [31], and the aug-cc-pVDZ basis set [34] is used throughout. The PCM cavity is built as a series of interlocking spheres centered on each atomic nucleus and using the solvation model density (SMD) radii [35]. The integral equation formalism PCM (IEF-PCM) [36–38] with the continuous surface charge (CSC) [39] formalism is applied. The lowest two excited states for each irrep are computed for each molecule. The results are reported as signed cumulative error bars (in eV) normalized for the number of states considered, using the SS results as reference. The absolute values of the excitation energies are reported in the supporting information. The results for MCP in methanol are reported in Fig. 1. Even from a cursory visual inspection, it is evident how the cLR results are much closer to those with the SS formalism than x0 and LR. The average error for x0 and LR is >0.1 eV in the nonequilibrium regime and >0.3 eV in the equilibrium regime. The fact that the former is smaller is due to the small values of e1 used in the calculation of the excited state dynamic charges. In other words, in the nonequilibrium regime the excited state charges are dominated

by their inertial contribution, which is frozen in the ground state response field. This inertial contribution is common to all of the approaches, thus the difference between them is small. Nevertheless, the results obtained with the cLR approach provide an excellent approximation of the SS results, with negligible errors (<0.01 eV) for the nonequilibrium case, and average errors of the order of 0.05 eV in the equilibrium case. Considering the response of the T amplitudes in the density used for the cLR approach seems to have very little effect on the results. In terms of states ordering, the SS formalism stabilizes the forth and fifth states (1A2 and 1B1 symmetry) from the x0 order to become the third and fourth states both in the nonequilibrium and equilibrium solvation regimes (see supporting material). The LR formalism does not recover this effect, while the cLR approach does. However, in cLR the 1A2 and 1B1 states become basically degenerate in the equilibrium regime (the difference is <0.01 eV) and switch order. This is of little consequence since the difference in energy between these two states in the SS formalism is about 0.04 eV, thus below the average error of cLR. The errors for acrolein in water are plotted in Fig. 2. Also in this case, the average errors for x0 and LR are much larger than for cLR. The former approximations show average errors around 0.1 eV for nonequilibrium, and about 0.35 eV for equilibrium solvation. On the other hand, the errors with both cLR approximations are about 0.01 eV for the nonequilibrium regime, and <0.15 eV for the equilibrium regime. The difference between cLR and cLR0 is also negligible. LR underestimates the lowest A0 state by about 0.25 eV and 0.37 eV for nonequilibrium and equilibrium solvation, respectively. All other approximations overestimate the reference values. The results for acetone in water are reported in Fig. 3. The difference in average error between x0 and LR compared to cLR is

Fig. 1. Cumulative normalized error (eV) compared to SS results for MCP in methanol. The normalization is with respect to the number of states considered.

Fig. 2. Cumulative normalized error (eV) compared to SS results for acrolein in water. The normalization is with respect to the number of states considered.

3. Numerical applications

Please cite this article in press as: M. Caricato, Comput. Theoret. Chem. (2014), http://dx.doi.org/10.1016/j.comptc.2014.02.001

M. Caricato / Computational and Theoretical Chemistry xxx (2014) xxx–xxx

5

Fig. 4. Cumulative normalized error (eV) compared to SS results for water in water. The normalization is with respect to the number of states considered. Fig. 3. Cumulative normalized error (eV) compared to SS results for acetone in water. The normalization is with respect to the number of states considered.

even more striking for this system. For nonequilibrium solvation, the former two approaches provide an average error around 0.18 eV while the cLR approach provides error below 0.02 eV with both choices of density. For equilibrium solvation, the average error is around 0.55 eV for x0 and around 0.50 eV for LR, while the error for both cLR approximations is around 0.13 eV. Additionally, the order of the fourth and fifth states (2A2 and 2B2) is inverted for x0 and LR in the equilibrium regime. Also in this case, the relaxation of the T amplitudes in cK density has a negligible effect on the cLR excitation energy. The results of the last test system, water in water solution, are reported in Fig. 4. The difference between LR and cLR is the largest for this case. The average error for LR is above 0.2 eV while it is negligible for both cLR approximations in the nonequilibrium calculations. The average error for x0 is slightly below 0.25 eV. The results for the equilibrium calculations further emphasize the difference between the LR and cLR results. The average error with LR is around 0.55 eV while it is less than 0.1 eV with cLR in both approximations. The error for x0 is about 0.65 eV, which shows that the state-specific corrections introduced with cLR recover most of the excited state solute–solvent polarization that cannot be recovered with the LR formalism. For this molecule, the order of the fifth and sixth states (1B2 and 2A1) in the equilibrium regime is inverted between SS and LR (and x0 ) while cLR follows the SS order. 4. Discussion and conclusions In this work, an efficient and accurate approximation is introduced to bridge the gap between the SS and LR solvation formal-

isms for the calculation of excitation energies in solution with the CCSD method and a polarizable solvation model. The approximate approach, introduced a few years ago for TDDFT [13] and called corrected-LR (cLR), has been extended here in the context of CC methods for equilibrium and nonequilibrium calculations.  D do not include the orbital Differently from Ref. [13], the charges Q relaxation for consistency with the usual CC formulation, where orbitals are not modified during the solution of the CC equations. Furthermore, the effect of the T amplitudes relaxation on the excitation energy is tested. The results in Section 3 for several test systems clearly show that the cLR approximation recovers much of the difference between the LR and SS formalisms at a negligible computational cost. In fact, this correction can be computed for multiple states at once, especially when the T relaxation is neglected. The latter seems indeed a very effective approximation since it avoids the calculation of the KK amplitudes (see Eq. (4)) at no apparent loss of accuracy. ~K density in Eq. With this approximation (i.e., with the use of the c (7)), the cLR correction is basically free (assuming that one has solved for both the LK and RK eigenvectors). The average difference between SS and LR for this test set in the nonequilibrium calculation varies in a range of 0.1–0.2 eV while it is an order of magnitude smaller with cLR. Note that using the x0 approximation (i.e., performing the excited state calculation in the frozen field of the ground state charges Q 0 þ Q N ) provides an average error of the same order of the LR formalism. Moreover, the x0 and cLR approximations consistently overestimate the SS results (albeit with different magnitude) while the LR results are in some cases below the SS values (i.e., MCP and acrolein, see Figs. 1 and 2). For the equilibrium regime, the difference is even more dramatic. The average discrepancy between SS and LR ranges between 0.3

Please cite this article in press as: M. Caricato, Comput. Theoret. Chem. (2014), http://dx.doi.org/10.1016/j.comptc.2014.02.001

6

M. Caricato / Computational and Theoretical Chemistry xxx (2014) xxx–xxx

and 0.55 eV, while the range is reduced to 0.05–0.15 eV for cLR. Also in this regime, the difference between the x0 approximation and the LR formalism is rather small, with the former being slightly worse (the largest difference is found for water, where the average error for x0 is 0.65 eV while it is 0.55 eV for LR). It is remarkable how a one-time correction to the x0 value is able to capture most of the polarization effects introduced in the SS formalism. The overall small difference between the LR and x0 transition energy is also interesting since it shows that the effect of the transition charges Q RN and Q LN on the transition energy is rather small even in the equilibrium regime. The results for the test cases reported in this work clearly show the potential for the cLR approach to obtain a very good approximation of the SS transition energies at a fraction of the computational cost. This is certainly a great advantage when accurate and computationally intensive methods are employed, like those belonging to the CC family. The level of agreement with the SS results is excellent for the calculation of vertical excitation energies. For equilibrium calculations, cLR can be successfully used for exploratory calculations on multiple states before a full SS calculation is performed on the state of interest to investigate, for instance, how the solvent affects the order of the states.

[14]

[15] [16]

[17]

[18]

[19]

[20]

[21]

[22]

Acknowledgement [23]

I wish to thank Prof. Benedetta Mennucci for fruitful discussions. Appendix A. Supplementary material Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.comptc.2014.02. 001. References [1] J.F. Stanton, R.J. Bartlett, The equation of motion coupled-cluster method – a systematic biorthogonal approach to molecular-excitation energies, transition-probabilities, and excited-state properties, J. Chem. Phys. 98 (9) (1993) 7029–7039. [2] J.F. Stanton, Many-body methods for excited-state potential-energy surfaces. 1. General-theory of energy gradients for the equation-of-motion coupledcluster method, J. Chem. Phys. 99 (11) (1993) 8840–8847. [3] J.F. Stanton, J. Gauss, Analytic energy gradients for the equation-of-motion coupled-cluster method – implementation and application to the HCN/HNC system, J. Chem. Phys. 100 (6) (1994) 4695–4698. [4] M. Kallay, J. Gauss, Calculation of excited-state properties using general coupled-cluster and configuration-interaction models, J. Chem. Phys. 121 (19) (2004) 9257–9269. [5] R.J. Bartlett, M. Musial, Coupled-cluster theory in quantum chemistry, Rev. Mod. Phys. 79 (1) (2007) 291–352. [6] I. Shavitt, R.J. Bartlett, Many-Body Methods in Chemistry and Physics, Cambridge University Press, Cambridge, 2009. [7] M. Aguilar, F.O. Del Valle, J. Tomasi, Nonequilibrium solvation: an ab initio quantum-mechanical method in the continuum cavity model approximation, J. Chem. Phys. 98 (1993) 7375. [8] H.J. Kim, J.T. Hynes, Equilibrium and nonequilibrium solvation and solute electronic structure. I. Formulation, J. Chem. Phys. 93 (7) (1990) 5194. [9] B. Mennucci, R. Cammi, J. Tomasi, Excited states and solvatochromic shifts within a nonequilibrium solvation approach: a new formulation of the integral equation formalism method at the self-consistent field, configuration interaction, and multiconfiguration self-consistent field level, J. Chem. Phys. 109 (7) (1998) 2798–2807. [10] J. Tomasi, B. Mennucci, R. Cammi, Quantum mechanical continuum solvation models, Chem. Rev. 105 (8) (2005) 2999–3093. [11] R. Cammi, S. Corni, B. Mennucci, J. Tomasi, Electronic excitation energies of molecules in solution: state specific and linear response methods for nonequilibrium continuum solvation models, J. Chem. Phys. 122 (10) (2005) 104513. [12] S. Corni, R. Cammi, B. Mennucci, J. Tomasi, Electronic excitation energies of molecules in solution within continuum solvation models: investigating the discrepancy between state-specific and linear-response methods, J. Chem. Phys. 123 (13) (2005) 134512. [13] M. Caricato, B. Mennucci, J. Tomasi, F. Ingrosso, R. Cammi, S. Corni, G. Scalmani, Formation and relaxation of excited states in solution: a new time dependent

[24]

[25]

[26] [27]

[28]

[29] [30]

[31]

[32] [33]

[34]

[35]

[36]

[37]

polarizable continuum model based on time dependent density functional theory, J. Chem. Phys. 124 (12) (2006) 124520. S. Miertus, E. Scrocco, J. Tomasi, Electrostatic interaction of a solute with a continuum – a direct utilization of abinitio molecular potentials for the prevision of solvent effects, Chem. Phys. 55 (1) (1981) 117–129. O. Christiansen, K.V. Mikkelsen, A coupled-cluster solvent reaction field method, J. Chem. Phys. 110 (3) (1999) 1365–1375. O. Christiansen, K.V. Mikkelsen, Coupled cluster response theory for solvated molecules in equilibrium and nonequilibrium solvation, J. Chem. Phys. 110 (17) (1999) 8348–8360. A. Osted, J. Kongsted, K.V. Mikkelsen, O. Christiansen, A CC2 dielectric continuum model and a CC2 molecular mechanics model, Mol. Phys. 101 (2003) 2055–2071. J. Kongsted, T.B. Pedersen, A. Osted, A.E. Hansen, K.V. Mikkelsen, O. Christiansen, Solvent effects on rotatory strength tensors. 1. Theory and application of the combined coupled cluster/dielectric continuum model, J. Phys. Chem. A 108 (2004) 3632–3641. M. Caricato, B. Mennucci, G. Scalmani, G.W. Trucks, M.J. Frisch, Electronic excitation energies in solution at equation of motion CCSD level within a state specific polarizable continuum model approach, J. Chem. Phys. 132 (8) (2010) 084102. R. Cammi, R. Fukuda, M. Ehara, H. Nakatsuji, Symmetry-adapted cluster and symmetry-adapted cluster-configuration interaction method in the polarizable continuum model: theory of the solvent effect on the electronic excitation of molecules in solution, J. Chem. Phys. 133 (2) (2010) 024104. R. Fukuda, M. Ehara, H. Nakatsuji, R. Cammi, Nonequilibrium solvation for vertical photoemission and photoabsorption processes using the symmetryadapted cluster–configuration interaction method in the polarizable continuum model, J. Chem. Phys. 134 (10) (2011) 104109. R. Cammi, Quantum cluster theory for the polarizable continuum model. I. The CCSD level with analytical first and second derivatives, J. Chem. Phys. 131 (16) (2009) 164104. R. Cammi, Coupled-cluster theories for the polarizable continuum model. II. Analytical gradients for excited states of molecular solutes by the equation of motion coupled-cluster method, Int. J. Quantum Chem. 110 (2010) 3040–3052. R. Cammi, Coupled-cluster theory for the polarizable continuum model. III. A response theory for molecules in solution, Int. J. Quantum Chem. 112 (2012) 2547–2560. M. Caricato, G. Scalmani, G.W. Trucks, M.J. Frisch, Coupled cluster calculations in solution with the polarizable continuum model of solvation, J. Phys. Chem. Lett. 1 (2010) 2369–2373. M. Caricato, CCSD-PCM: improving upon the reference reaction field approximation at no cost, J. Chem. Phys. 135 (2011) 074113. M. Caricato, G. Scalmani, M.J. Frisch, Brueckner doubles coupled cluster method with the polarizable continuum model of solvation, J. Chem. Phys. 134 (2011) 244113. M. Caricato, Exploring potential energy surfaces of electronic excited states in solution with the EOM-CCSD-PCM method, J. Chem. Theory Comput. 8 (2012) 5081–5091. M. Caricato, Absorption and emission spectra of solvated molecules with the EOM-CCSD-PCM method, J. Chem. Theory Comput. 8 (2012) 4494–4502. M. Caricato, F. Lipparini, G. Scalmani, C. Cappelli, V. Barone, Vertical electronic excitations in solution with the EOM-CCSD method combined with a polarizable explicit/implicit solvent model, J. Chem. Theory Comput. 9 (2013) 3035–3042. M. Caricato, A comparison between state-specific and linear-response formalisms for the calculation of vertical electronic transition energy in solution with the CCSD-PCM method, J. Chem. Phys. 139 (2013) 044116. G.D. Purvis, R.J. Bartlett, A full coupled-cluster singles and doubles model – the inclusion of disconnected triples, J. Chem. Phys. 76 (4) (1982) 1910–1918. M.J. Frisch, G.W. Trucks, H.B. Schlegel, G.E. Scuseria, M.A. Robb, J.R. Cheeseman, G. Scalmani, V. Barone, B. Mennucci, G.A. Petersson, H. Nakatsuji, M. Caricato, X. Li, H.P. Hratchian, A.F. Izmaylov, J. Bloino, G. Zheng, J.L. Sonnenberg, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, J.A. Montgomery Jr., J.E. Peralta, F. Ogliaro, M. Bearpark, J.J. Heyd, E. Brothers, K.N. Kudin, V.N. Staroverov, T. Keith, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J.C. Burant, S.S. Iyengar, J. Tomasi, M. Cossi, N. Rega, J.M. Millam, M. Klene, J.E. Knox, J.B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R.E. Stratmann, O. Yazyev, A.J. Austin, R. Cammi, C. Pomelli, J.W. Ochterski, R.L. Martin, K. Morokuma, V.G. Zakrzewski, G.A. Voth, P. Salvador, J.J. Dannenberg, S. Dapprich, P.V. Parandekar, N.J. Mayhall, A.D. Daniels, O. Farkas, J.B. Foresman, J.V. Ortiz, J. Cioslowski, D.J. Fox, Gaussian Development Version, Revision h.09+ Edition, Gaussian, Inc., Wallingford, CT, 2010. T.H. Dunning, Gaussian basis sets for use in correlated molecular calculations. I. The atoms boron through neon and hydrogen, J. Chem. Phys. 90 (2) (1989) 1007–1023. A.V. Marenich, C.J. Cramer, D.G. Truhlar, Universal solvation model based on solute electron density and on a continuum model of the solvent defined by the bulk dielectric constant and atomic surface tensions, J. Phys. Chem. B 113 (18) (2009) 6378–6396. E. Cances, B. Mennucci, J. Tomasi, A new integral equation formalism for the polarizable continuum model: theoretical background and applications to isotropic and anisotropic dielectrics, J. Chem. Phys. 107 (8) (1997) 3032–3041. B. Mennucci, E. Cances, J. Tomasi, Evaluation of solvent effects in isotropic and anisotropic dielectrics and in ionic solutions with a unified integral equation

Please cite this article in press as: M. Caricato, Comput. Theoret. Chem. (2014), http://dx.doi.org/10.1016/j.comptc.2014.02.001

M. Caricato / Computational and Theoretical Chemistry xxx (2014) xxx–xxx method: theoretical bases, computational implementation, and numerical applications, J. Phys. Chem. B 101 (49) (1997) 10506–10517. [38] F. Lipparini, G. Scalmani, B. Mennucci, E. Cancès, M. Caricato, M.J. Frisch, A variational formulation of the polarizable continuum model, J. Chem. Phys. 133 (1) (2010) 014106.

7

[39] G. Scalmani, M.J. Frisch, Continuous surface charge polarizable continuum models of solvation. I. General formalism, J. Chem. Phys. 132 (11) (2010) 114110.

Please cite this article in press as: M. Caricato, Comput. Theoret. Chem. (2014), http://dx.doi.org/10.1016/j.comptc.2014.02.001