Effective Fragment

Feature Article pubs.acs.org/JPCA

Hybrid Equation-of-Motion Coupled-Cluster/Effective Fragment Potential Method: A Route toward Understanding Photoprocesses in the Condensed Phase Debashree Ghosh* Physical and Materials Chemistry Division, CSIR-National Chemical Laboratory, Pune, India 411008 ABSTRACT: The prediction of accurate solvatochromic shifts to the electronic excited states of chromophores is a challenge, especially in the complex biological phase, due to the importance of long-range electrostatic interactions. Hybrid quantum mechanical/ molecular mechanical (QM/MM) methods are generally employed for the calculation of quantum mechanical properties in complex systems. To be predictive, there is need for an accurate quantum mechanical method that can depict the charge transfer states correctly and incorporate higher than single excited determinants in its linear response ansatz. On the contrary, for the correct depiction of the environment interactions (MM region), one needs to account for polarizability in a balanced manner. These two challenges are successfully addressed by the recently developed hybrid quantum mechanical/effective fragment potential (QM/EFP) methods, with equation-of-motion coupled-cluster (EOM-CC) as the QM method of choice. The result is an efficient method to estimate excitation energy, ionization energy, electron affinity, and redox potential in the condensed phase. It has further been extended to biological systems.

■

INTRODUCTION

Here it should be noted that the recent development of computational capabilities and methodologies have triggered intense activity in the theoretical community to understand the processes in DNA damage, light harvesting, etc.17−22 Dominant pathways for the photoprotection and photodamage of isolated DNA bases have been studied extensively. There has been similar success from the theoretical perspective in various other photoprocesses of small molecules.23−32 Quantum mechanical (QM) methods are, however, computationally expensive and the computational scaling of these methods are O(N5)−O(N6) or higher. Thus, although it is possible to treat small to medium sized systems (tens of atoms) with QM methods, such calculations are unfeasible for large biological systems or systems in condensed phase with hundreds or thousands of atoms. Because the electrostatic interactions are long-range in nature, consideration of only the chromophore is not enough to understand the excited states in condensed phase and, therefore, one needs to consider the effect of environment interactions with the chromophore. For these situations, hybrid quantum mechanical/molecular mechanical (QM/MM) methods33−36 have been developed over the last few decades. The hybrid QM/MM methods combine the strengths of both QM methods (accuracy) and MM methods (speed). In these methods, the important part of the

Light matter interaction is the cornerstone of the functioning of many processes in biology, such as light harvesting, photosynthesis, vision, and photoinduced catalysis.1−6 Understanding such processes involves a synergy between time dependent spectroscopy and computational modeling of the processes, which requires the understanding of the electronic structure of the excited states. Over the last few decades, quantum mechanical methods have been well developed for the closed shell and ground state systems. This development has been fueled by the improvements in algorithm and computer speed. However, predictive calculations of excited states and open shell species are still nontrivial, especially for large systems. This is mainly becasue, in wave-function-based approaches, the orbitals are optimized for the accurate estimation of ground states and, therefore, in the excited states one finds large numbers of near-degeneracies in the electronic configurations of these systems. These neardegeracies are not treated adequately by the independent orbital (electron) model or the single reference assumption (the assumption that each state can be depicted by a single configuration of the molecular orbitals). Equation-of-motion coupled-cluster (EOM-CC) is a robust and sophisticated method to treat small to medium sized open shell systems and excited states, based on the consideration of important parts of the Fock space.7−15 Reference 16 gives a detailed overview of these methods. © 2016 American Chemical Society

Received: August 16, 2016 Revised: October 5, 2016 Published: October 24, 2016 741

DOI: 10.1021/acs.jpca.6b08263 J. Phys. Chem. A 2017, 121, 741−752

Feature Article

The Journal of Physical Chemistry A

the effective fragments (EFs) consist of electrostatic, polarization, dispersion, and exchange−repulsion interactions, whereas the interactions between EFs and the QM region consist of electrostatics and polarization components. To understand the excited state phenomena, hybrid QM/EFP, specifically hybrid equation-of-motion coupled-cluster/effective fragment potential (EOM-CC/EFP), has been developed.81−84 It is important to note that hybrid QM/EFP schemes have been developed with other QM methods such as CIS/EFP,85 CIS(D)/EFP, and TD-DFT/EFP.86 In this paper, we will discuss the EOM-CCSD and EFP methods and the hybrid formalism for EOM-CCSD/EFP in details. The capability of this approach is illustrated by the various excited stated properties that have been calculated with this hybrid formalism in our group.

system (active part) is treated quantum mechanically, i.e., the subsystem where the reaction or electron transfer happens. The rest of system is treated molecular mechanically, i.e., the subsystem that provides an effective perturbation/field (spectator part) on the QM part. Thus, the Hamiltonian can be written as, Ĥ = Ĥ QM + ĤMM + Ĥ QM/MM (1) ̂ ̂ where Ĥ QM and Ĥ MM are the Hamiltonians for the quantum mechanical (QM) and molecular mechanical (MM) parts of the system that are treated with two different kinds of physics (quantum and classical). Ĥ̂ QM/MM is the interaction between the QM and MM part. Alternative directions taken by various groups toward lower computational scaling have been in local orbital-based linear scaling methods,37−41 ONIOM, and fragmentation schemes.42−50 The simplest form of MM methods consists of point charges on the atoms in the MM subsystem,51−58 where the interaction between QM and MM is treated by Coulomb’s law,

EQM/MM =

∑∫ i

QM

ρQM ( r ⃗)qi | r ⃗ − ri |⃗

■

THEORY Effective Fragment Potential. EFP is a sophisticated, discrete molecular mechanics method; i.e., each solvent/ environment molecule is treated as a separate entity or fragment (Figure 1). It is an empirical parameter-free (based

dr (2)

where qi is the charge on the ith MM point and ρQM(r)⃗ is the electron density on the QM part at position r.⃗ However, the traditional QM/MM with simple electrostatics is not accurate enough to predict solvatochromic shifts of excited state processes. Reference 59 shows the need to go beyond electrostatics for even qualitative prediction of the red/blue shift of the UV−vis spectra in condensed phase. Using Rayleigh−Schrodinger perturbation theory (RSPT), it can be shown that the first-order perturbation gives rise to electrostatic or Coulomb interactions, i.e, the interaction energy due to the static charges of the two interacting bodies.60,61 Polarization and dispersion interactions are due to the secondorder perturbation corrections. Polarization of the MM part (spectator and inactive) by the QM part (active) changes the effective field that is felt by the QM part and consequently the excited state energies and properties. Reference 36 establishes the need for polarization along with electrostatics as crucial to the quantitative estimation of the solvatochromic shifts. The absence of polarization is the root cause behind the inadequacy of traditional point charge QM/MM method. Furthermore, there is need for self-consistency in estimating the effect of polarization because the QM part is affected by the MM part charges, and the MM charges are in turn modified by the interaction with the QM charge density and so on. Several approaches such as polarizable embedding62,63 and various flavors of continuum models64−69 have dealt with the accurate treatment of polarization in the hybrid scheme and it continues to be an active field of research. Thus, it would be ideal if a force field, based on ab initio calculations, i.e., without empirical parameters, is capable of depicting the electrostatics and polarization of the environment accurately. Effective fragment potential (EFP)70−75 has been developed by Gordon and co-workers, to mainly treat the effects of water solvation very accurately, originally called EFP1. In this method, the potentials due to each fragment (typically a molecule) have been calculated from an a priori ab initio calculation. This approach has been later extended for any molecule.76 EFP provides us with a sophisticated and empirical parameter-free polarizable force field. The details of the methods can be found in refs 77−80. The interactions between

Figure 1. Effective fragment potential as a discrete molecular mechanics method. The figure shows the solute molecule (thymine) considered as a quantum mechanical region in the environment of water molecules, each of which are considered as EFs.

on a priori ab initio calculation) method. The first implementation of EFP in GAMESS87,88 was targeted at understanding solvation in water and later generalized to describe other solvent molecules. The interactions between two EFs can be divided into long-range (those that depend on 1/ RN) and short-range interactions (that depends exponentially on R) terms. The long-range terms can be divided into electrostatic, polarization, and dispersion interactions. The short-range term deals with the exchange−repulsion and screening or damping terms. The total interaction energy between two effective fragments is given by E EFP − EFP = ECoul + Epol + Edisp + Eexch − rep

(3)

where ECoul, Epol, Edisp, and Eexch−rep are the Coulomb, polarization, dispersion, and exchange−repulsion energies, respectively. The charge transfer interaction energies can also be included in this description. However, in our implementa742


Feature Article

The Journal of Physical Chemistry A tion in Q-CHEM,89 charge transfer interaction energies are omitted. The Coulomb or electrostatic interaction is the leading term in the interactions between two fragments or molecules. The most accurate way to treat it would be to calculate the Coulomb interactions between charge densities. However, to reduce computational cost, one can approximate the charge densities as an expansion of multipoles (charges, dipoles, quadrupoles, etc.) over multiple points. In EFP, Stone’s distributed multipole analysis (DMA)90,91 up to octupoles, which involves distribution of multipoles based on the basis set space, is used. Thus, the Coulomb potential is expanded as V Coul ≈ V charge + V dipole + V quadrupole + V octupole

Coul EQM/EFP =

k

+

C6 R6

μak ,

k Θab ,

∑ μak a

Fak R k3

+

∑ Θkab a,b

k Fab

3R k5

k ⎫ ⎪ Fabc ⎬ρ ( r ) dr 7 ⎪ QM ⃗ 5R k ⎭

(7)

k Ωabc

where q , and are the charges, dipoles, quadrupoles, and octupoles at the multipole points k and a, b, c denote the x, y, z coordinates. The Fka, Fkab, Fkabc terms refer to the solute electric field, field gradient, and field Hessian, respectively, and Rk denotes the position of the kth multipole point. Polarization is the other component of the QM/MM interaction and, in fact, is a crucial component in many properties, such as enzyme catalysis, redox potentials, ionization energies, etc. The polarization of the EFs and the QM region induces an effective induced dipole moment on each EF. The effect of this induced dipole can be calculated in a manner very similar to that for the effect of the permanent dipole, as shown in eq 7,

(4)

pol EQM/EFP =

1 2

x ,y,z

∑ ∑ [μã p Fap,(ai) − μã p Fap,(mult + nuc)] p

a

(8)

where p denotes the polarizability points on the LMO centers and μ̃pa and μ̃ã p denotes the induced dipole moment and conjugate induced dipole moment. The Fp denotes the field at the point p. The computation of induced dipole moments, μ̃ pa, is performed self-consistently. The details of the implementations of various flavors of QM/EFP can be found in ref 86. Equation-of-Motion Coupled-Cluster/Effective Fragment Potential. An extremely successful approach toward excited states of molecular systems is the linear response or equation-of-motion formalism.92−96 In this approach, the excited state is built as a linear combination of derivatives of the ground state wave function. Time dependent density functional theory (TD-DFT) and configuration interaction singles (CIS) are the most well-known examples of this formalism. The equation-of-motion (EOM) is a similar formalism to the configuration interaction (CI). It solves for excited states by diagonalizing the similarity transformed Hamiltonian, H¯ = e−THeT. The similarity transformation involved in a coupled cluster wave function incorporates the dynamic correlation of the predominantly single reference state, typically the closed shell ground state. The EOM formalism solves for the excited state or target state by using the Fock space around the reference state (CCSD). Therefore, its success stands on finding a good reference state from which low order (singles, doubles) excitations can describe the target state. EOM-CCSD can rigorously treat excited states with considerable multireference nature as long as there is a closed shell well behaved reference state. There are various flavors of EOM-CCSD depending on the reference state and the excited state that is targeted, e.g., EOM-IP-CCSD, EOM-EE-CCSD, EOM-SF-CCSD, etc. However, all of the EOM-CCSD methods have steep computational scaling (O(N6)). To alleviate this problem, there have been development of perturbative approximations to the EOM-CCSD method by truncating the Baker−Campell−Hausdorff expansion at various levels. The second-order perturbative approximation to the EOM-CCSD is denoted as EOM-MP2 or EOM-CCSD(2) and has been the most useful approximation because it reduces the

(5)

where μ̃ ⃗ is the induced dipole moment that is created due to the electric field F⃗, given that α is the polarizability tensor. The dispersion interaction can be expanded in the London dispersion series of the form R−N, where N = 6, 8, 10, .... In EFP, the dispersion term is truncated at the first term and takes the form Edisp =

∑ Ωkabc a,b,c

k

The distributed multipole points for the EFs are placed on the atoms and the bond centers. Although this approximation is accurate when the fragments are far away from each other, it starts to show considerable error at short distances, due to the charge penetration effect, which is not accounted for in the distributed multipole analysis and, therefore, an exponentially decaying (short-range) damping term (1 − e−αR) is added to the electrostatic interaction terms. Apart from the effect of the static multipoles, an important component of the interactions is polarization, i.e., the selfconsistent change in the charge density of one fragment due to the presence of the other fragment. It is expressed in the form of an induced dipole formed on fragment A due to the permanent multipoles present on the fragment B. This induced dipole is calculated as a tensor sum on the localized molecular orbitals (LMOs) in a fragment (distributed polarizability) to improve convergence and accuracy.

μ ̃ ⃗ = αF ⃗

⎧ k ⎪q ⎨ k + QM ⎪ ⎩R

∑∫

(6)

where C6 is expressed as sum over induced dipole−induced dipole interactions calculated at the LMOs. The short-range exchange−repulsion term is derived as an expansion over overlap integrals calculated in the LMO basis and truncated at the quadratic term. Hybrid Quantum Mechanics/Effective Fragment Potential. Hybrid QM/MM involves expressing the Hamiltonian as a sum of Hamiltonians describing the quantum mechanical and molecular mechanical (classical mechanical) parts as well as the interactions between them (QM and MM), as given in eq 1. Though the purely quantum mechanical (QM) and pure classical mechanical part of the system are conceptually easy to solve, it is the interaction between the QM and MM part that gives rise to the different approaches of QM/MM. The QM/ MM interaction in EFP accounts for the Coulomb and polarization interactions. The Coulomb interaction terms are treated in a way analogous to the point charge MM methods, i.e. 743


Feature Article

The Journal of Physical Chemistry A scaling of some of the flavors of EOM-CCSD, viz. EOM-IPCCSD.97−100 However, the computational scaling of the perturbative approximation is still quite steep and can be applicable to smallmedium sized molecules. Thus, for the EOM-CCSD methods to be applicable to large systems in the condensed phase, the hybrid EOM-CC/EFP formalism has been developed. This excited state QM/EFP is different from its ground state counterpart, because the electron density of the QM part in the excited state (ρexc(r)⃗ ) is significantly different from the ground state electron density (ρgr(r)⃗ ). The differences in electron densities cause different electric fields on the EFs due to the excited and ground states, and thus, the induced dipoles that are formed on the EFs will be different for these states. To account for this differential effect of the ground and excited states on the EFs, there are two approaches: (i) solving the selfconsistent polarizabilities for each state separately in a state specific way; (ii) perturbative correction due to the effect of the induced dipole due to different states depending on their electron densities. Both these methods have been tested and implemented.85 It has been found that solving the polarization response of the environment self-consistently for each state can be tedious and prone to convergence problems, and therefore, for the rest of the review, we will discuss the results as obtained from the perturbative approximation.

Both these benchmark papers on noncovalent interactions show that the accuracy of EFP interaction energies is comparable to ab initio methods such as MP2 and SCS-MP2. Furthermore, the decomposition of interaction energies shows that the components are similar to SAPT-based methods. On the contrary, computationally EFP is significantly more affordable than these ab initio methods.

■

SOLVATOCHROMIC SHIFT The effect of the solvent (condensed phase) on the solute can in most cases be regarded as the effect of differential stabilization of the solute in the ground and excited states. This gives rise to the resultant solvatochromic shifts. Figure 2

■

Figure 2. Schematic of the cause behind the solvatochromic shift of excited state energies. The solvent interactions typically stabilize the ground and the excited states of a solute molecule (chromophore). It is the differential stabilization that gives rise to the solvatochromic shift.

NONCOVALENT INTERACTIONS Because the EFP method decomposes interaction energies between fragments (molecules) into electrostatic (Coulomb), polarization, dispersion, and exchange−repulsion energy terms, along with an accurate force field, it offers an energy decomposition scheme, similar in spirit to symmetry-adapted perturbation theory.101,102 This can also be used as a method to accurately estimate noncovalent interactions,103−106 i.e., interactions between fragments (typically molecules) of closed shell systems that do not involve significant electron sharing. These noncovalent interactions have been quite problematic for quantum chemical methods, because the small magnitude of the interactions requires high accuracy. Thus, noncovalent interactions can be properly captured only with high-level correlated methods, such as the coupled cluster with single, double, and perturbative triple excitations (CCSD(T)), combined with large augmented basis sets. We have shown that the EFP interaction energies of nucleic acid base (NAB) dimers (H-bonded and stacked) are in good agreement (within 1−3 kcal/mol) with the interaction energies calculated with high-level ab initio methods.78 However, because EFP is computationally affordable when compared to the high-level ab initio methods, system sizes up to 20 NABs could be studied to investigate the asymptotic behavior of different components of interaction energies. It is noticed that dispersion interactions are the leading component in the π−π interactions for all the oligomers. Slipchenko and co-workers have benchmarked the components of EFP interaction energies for the S22 data set.107 They find that the mean absolute deviation (MAD) for EFP is significantly lower than DFT-based methods (0.57 kcal/mol vs 4.5−6.2 kcal/mol) for dispersion-dominated interacting species. On the contrary, EFP performs significantly better than classical force fields, such as Amber, OPLSAA, and MMFF94, in the case of H-bonded species (1.82 kcal/mol vs 3.6−4.4 kcal/mol). Overall the MAD of EFP is 0.89 kcal/mol.

shows the schematic of the interactions that can give rise to red or blue solvatochromic shift. Furthermore, the stabilization depends on the QM/MM interaction energies from both Coulomb and polarization interactions. In the next subsections, we analyze the accuracy of QM/EFP treatment on vertical ionization, excitation energies, and redox potentials and extend the formalism from simple small molecule solvents to protein environments. Excitation Energy. Slipchenko and co-workers have used p-nitroaniline (pNA) as the test solute molecule to benchmark the effect of using the hybrid EOM-CCSD/EFP and CIS(D)/ EFP approach.84 Small clusters of pNA−water show that the QM/EFP approach is accurate to within 0.05 eV of the parent QM method. It is further noticed that Coulomb interactions contribute to a large part of the solvatochromic shift and second-order perturbation has a much smaller effect. Once the QM/EFP approach is benchmarked for microsolvated structures with respect to full QM method, the natural step is to employ this technique to predict the effect of bulk solvation and extract molecular understanding. CIS(D)/EFP calculations, performed on geometries obtained from classical MD simulations or EFP MD simulations, show that the CIS(D)/ EFP average vertical excitation energies (VEE) and the full width at half maxima (fwhm) are an excellent match to the experiment. The solvatochromic shift to VEE (ΔVEE) calculated by EFP is 1.0 eV compared to the experimental value of 1.0 eV, and the fwhm is estimated at 0.46 eV as compared to 0.6 eV from experiment (shown in Table 1). Here it should be noted that the prediction of the average properties need significantly less sampling (number of snapshots from MD simulation) than the prediction of the accurate line shape 744


Feature Article

The Journal of Physical Chemistry A Table 1. Comparison of Experimental and Hybrid QM/EFP Solvatochromic Shifts (eV) in Various Polar and Nonpolar Solvents experimental hybrid QM/EFP

water

1,4-dioxane

cyclohexane

1.0 1.0

0.72 0.44

0.41 0.20

or fwhm and the error in fwhm could have originated from insufficient sampling. However, QM/EFP shows large errors on the solvatochromic shifts in the case of solutes in nonpolar solvents, such as 1,4-dioxane and cyclohexane (

Effective Fragment

Recommend Documents