The Two-Component Polarizable-Embedding ... - ACS Publications

Subscriber access provided by - Access paid by the | UCSB Libraries

Article

Approaching phosphorescence lifetimes in solution: The twocomponent polarizable-embedding approximate coupled-cluster method Katharina Krause, Mirko Bauer, and Wim Klopper J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.6b00239 • Publication Date (Web): 09 May 2016 Downloaded from http://pubs.acs.org on May 10, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Theory and Computation is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Approaching phosphorescence lifetimes in solution: The two-component polarizable-embedding approximate coupled-cluster method Katharina Krause,† Mirko Bauer,‡ and Wim Klopper∗,† †Institute of Physical Chemistry, Theoretical Chemistry Group, Karlsruhe Institute of Technology, KIT Campus South, P. O. Box 6980, 76049 Karlsruhe, Germany ‡Mulliken Center for Theoretical Chemistry, Institute of Physical and Theoretical Chemistry, University of Bonn, Beringstr. 4, 53115 Bonn, Germany E-mail: [email protected] Fax: +49 721 60847225

Abstract The theoretical description of phosphorescence lifetimes in the condensed phase requires a method which takes into account both spin–orbit coupling and solvent– solute interactions. To obtain such a method, we have coupled our recently developed two-component coupled-cluster method with singles and approximated doubles to a polarizable environment. With this new method we investigate how different solvents effect the electronic phosphorescence energies and lifetimes of 4H-pyran-4-thione.

1 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Introduction

Many UV/Vis experiments are conducted in the condensed phase, even though one is usually interested in the absorption and emission of the solute only. However, the solvent may influence these processes. For example, by screening partial charges, states with a dipole moment are stabilized by polar solvents. Solvent induced lifetime differences of excited states may even lead to different reaction pathways in photochemical reactions. 1 Although there is a large variety of well-established methods for the computation of optical spectra in vacuum, 2 calculating those in solution still proves challenging. One reason is the increased system size. While the quantum-chemical calculation of properties of, for example, acetone requires the treatment of a system containing 32 electrons, placing the same molecule in solution, for example in a droplet of water molecules with a radius of only 10 Å, drastically increases the system size to about 1400 electrons. Thus, the straightforward quantum-mechanical (QM) description of the whole system is not feasible. Another aspect is the increased number of degrees of freedom. Therefore, the dynamics of the whole system play a more pronounced role than in the gas phase. For the prediction of optical transitions within a molecule in the gas phase, it is usually sufficient to focus on its equilibrium structure in the initial state. For ”liquid” systems, however, a representative sample of different structures of the solute coordinated by solvent molecules is required. To overcome the system size problem, several approaches have been developed in the past. The key ideas will be presented briefly, for a detailed discussion we refer to a recent review. 3 Based on the fragmentation of the solvent–solute supersystem into several subsystems, for instance each containing a single solvent or solute molecule, the frozen density embedding ansatz provides a quantum-chemical description of each subsystem yielding subsystem densities. 4 The interaction of the subsystems is then modeled by density functional theory. 5 Because of the division into subsystems, the contributions of the solute to a UV/Vis spectrum can be trivially separated from bands due to the solvent. 6 This advantage is still present if the solvent is treated only on a molecular mechanics (MM) 2 ACS Paragon Plus Environment

Page 2 of 26

Page 3 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


level of theory by means of a force field. The transition energies and moments are then computed for the solute which is embedded in a field of multipoles and polarizable sites obtained from the force field. These two approaches both rely on explicit structures of the solute coordinated by solvent molecules, whereas the explicit treatment of solvent molecules is abandoned if one applies continuum models. Here, the solute is replaced by a homogeneous dielectric continuum, which is defined by only a few parameters, for example its dielectric constant. Yet, effects due to hydrogen bonding or different dynamics are not accounted for. All of these approaches have been successfully combined with accurate wavefunction methods such as the coupled-cluster method with singles and approximate doubles (CC2) or the algebraic diagrammatic construction scheme up to second order (ADC(2)). 7–11 With the resulting methods, UV/Vis absorption spectra have been computed and fluorescence processes investigated. 12–19 In the present paper, we proceed to the computation of phosphorescence lifetimes in solution. Since this is a spin-forbidden process, it is usually slow compared to fluorescence. Therefore, structural relaxations of not only the electronically excited solute but also the surrounding solvent shell have to be considered. For this purpose, we combine our recently developed 2c-RI-CC2 (two-component resolution of the identity coupled-cluster singles and approximate doubles) method, 20 which has yielded promising results for triplet lifetimes in vacuum, with the simplified polarizable-embedding (sPE) ansatz 21 adapting the approximations suggested by Schwabe and coworkers 7 to complex wavefunctions. This QM/MM ansatz has the advantage that the structure of the solvent shell can be relaxed explicitly to model the solvated triplet excited state. Because the solvent is included in the computation of the solute in terms of multipole moments and polarizable sites only, the resulting method should be applicable to medium sized molecules in a sphere of several hundred solvent molecules. In the next section, the theory of the simplified polarizable embedding ansatz applied to the complex 2c-RI-CC2 wavefunction and its underlying generalized Hartree–Fock (GHF) wave-



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 26

function is discussed. Details on the implementation of both the PE-GHF and 2c-sPERI-CC2 methods are given in section 3. With the new 2c-sPERI-CC2 method, the solvent effect on phosphorescence frequencies of an exemplary organic molecule is investigated in section 4.

2 2.1

Theory Polarizable-embedding generalized Hartree–Fock

In the following, the polarizable-embedding Hartree–Fock method will be presented for the case of a GHF wavefunction ansatz. The derivation closely follows the presentations in References 7 and 22. ˆ 2c-QM |GHF⟩ of the solute in vacuum is obtained as the expecThe energy E GHF = ⟨GHF|H ˆ 2c-QM which is set up according to tation value of the Hamiltonian H  ˆ 2c-QM = H

N ∑ i

 ˆ αα

ˆ αβ

 h (i) h (i)    + gˆ + Vnuc . ˆ βα (i) h ˆ ββ (i) h

(1)

ˆ N refers to the number of electrons of the solute, h(i) is a Dirac-like one-electron operator introducing spin–orbit coupling, gˆ is the non-relativistic electron–electron repulsion and Vnuc the repulsion of the nuclei. If this solute is placed in an environment of M molecular mechanics sites consisting of multipole moments and a dipole polarizability tensor each, the total energy of the system will contain two additional terms:

E PE-GHF = E GHF + E es + E pol .

(2)

The first one is the sum of electrostatic interactions of the multipole moments at the MM sites with the whole system, which can be separated into electronic contributions and those


Page 5 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


due to nuclei and multipole moments of the MM sites (ne, non-electronic),

ˆ es |GHF⟩ + E es,ne . E es = ⟨GHF|G

(3)

The latter term of Eq. (2), E pol , contains interaction terms due to the mutual polarization of the QM and MM subsystems. These terms are obtained based on the approximation that the induced dipole µind u at the MM site at r u is linearly proportional to the field F (r u ) ≡ F u at that site, 22

1 ∑ ( tot GHF )⊤ GHF ). ) Ruv F tot F u (D v (D 2 uv M

E pol = −

(4)

The total field is the sum of the electronic (el) and non-electronic (ne) contributions,

GHF GHF F tot ) = F el ) + F ne u (D u (D u

(5)

The classical response matrix R, which is real and symmetric, is obtained from the blockdiagonal matrix of the anisotropic dipole polarizability tensor α and the matrix of the dipole–dipole interaction tensor T dip,dip ( )−1 R = α−1 − T dip,dip .

(6)

For a detailed derivation of the polarization term, we refer to Reference 22. From the energy expression (2), the Fock operator of the polarizable-embedding generalized Hartree–Fock method is derived, ˆ PE =F ˆ 2c-QM + G ˆ es F −

M ∑ (

GHF F tot ) u (D

)⊤

(7) ˆv , Ruv F

uv



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 26

ˆ v is defined in second quantization with p, q referring where the operator of the electric field F ⃗ p (r 1 ): to general 2c-spinors ϕ ˆv =− F

∑ pq

=−

∑

⃗ p (r 1 )| ⟨ϕ

rv − r1 ⃗ |ϕq (r 1 )⟩ a ˆ†p a ˆq 3 |r v − r 1 | (8)

θ v,pq a ˆ†p a ˆq

.

pq

With this, the electric field due to an arbitrary density D can be written as

F el u (D) = −

∑

θ u,pq Dpq .

(9)

pq

Since the Hartree–Fock density D GHF is a complex but Hermitian matrix, the electric field F el (D GHF ) obtained from it is real. On a side note, we want to point out that even though some of the additional terms due to the polarizable environment in Eq. (7) are formally similar to the nuclear attraction potential, they are not treated on an equal footing. While the nuclear attraction potential is included in the decoupling procedure of the four-component one-electron Dirac equation which yields ˆ the operator h(i) given in Eq. (1), the transformation of the additional one-electron terms appearing in equation (7) is neglected. This approximation is justified since the resulting terms are smaller by at least one order of magnitude than the averaged nuclear–electron attraction.

2.2

Polarizable-embedding 2c-RI-CC2

The polarizable-embedding generalized Hartree–Fock state serves as the reference wavefunction for the two-component polarizable-embedding coupled-cluster singles and approximated doubles method. We define its Lagrangian L based on the simplified polarizable-embedding


Page 7 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


coupled-cluster method, 21 LsPECC =E PE-GHF

(10)

ˆ pol (D ∆CC )|CC⟩ . + ⟨Λ|ˆ gN + FˆNPE + G N The subscript N indicates that the operator is normal ordered with respect to the Hartree– Fock state. The elements of the difference density matrix D ∆CC are defined as GHF ∆CC , ˆq |CC⟩ − Dpq = ⟨GHF|ˆ a†p a Dpq

(11)

hence, the density matrix is complex and non-Hermitian. This subject is further addressed ∆CC in the appendix. From this density matrix the fields F el ) due to the coupled-cluster v (D

wavefunction are obtained via Eq. (9). In Eq. (10), we introduced the electronic polarization ˆ pol (D ∆CC ), which is computed from these fields according to operator G ˆ pol (D) = − 1 G 2

M ( ∑

û F

)⊤

Ruv F el v (D) .

(12)

uv

To obtain the 2c-sPECC2 method, the coupled-cluster expansion is truncated to single excitations Tˆ1 and double excitations Tˆ2 and all terms in second or higher order in the doubles equation are neglected. 23 To further reduce the computational demand of the resulting method, the PERI-CC2 approximations 7 are applied, namely, doubles contributions to the projected coupled-cluster density matrix are neglected, GHF ∆CCS , ˆq eT1 |GHF⟩ − Dpq = ⟨GHF|e−T1 a ˆ†p a Dpq ˆ

ˆ


(13)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 26

as well as all explicit coupling terms of double excitations with the polarization operator ˆ pol ≡ G ˆ pol (D ∆CCS ). Finally, the 2c-sPERI-CC2 Lagrangian reads G L =E PE-GHF 1 ˆ pol Tˆ1 |GHF⟩ + ⟨GHF|ˆ gN (Tˆ2 + Tˆ12 ) + G 2 ∑ ˆ˜ pol |GHF⟩ + t¯µ1 ⟨µ1 |[FˆNPE , Tˆ1 ] + gˆÑ + G µ1

+

∑

(14)

t¯µ1 ⟨µ1 |[gˆÑ , Tˆ2 ]|GHF⟩

µ1

+

∑

t¯µ2 ⟨µ2 |gˆÑ + [FˆNPE , Tˆ2 ]|GHF⟩ ,

µ2

where µ1 indicates a single excitation and µ2 a double excitation. The singles similarity ˆ ˆ transformation was abbreviated by the shorthand notation xˆ˜ = e−T1 xêT1 . Similarly, we

introduce the following notation for the matrix representation of an operator in Λ basis, ⃗ p |e−T1 xêT1 |ϕ⃗q ⟩ . x˜pq ≡ ⟨ϕ ˆ

ˆ

(15)

Since the coupling terms of double excitations with the polarization operator are neglected, explicit polarization contributions are obtained only for derivatives with respect to singles amplitudes or multipliers. The explicit expressions of all derivatives which are relevant for the computation of vertical excitation energies and oscillator strengths are given in Table 1, where we have introduced the permutation operator

Pîjab f (a, i, b, j) = f (a, i, b, j) + f (b, j, a, i) .


(16)

Page 9 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 1: Explicit polarization contributions to important quantities derived from the 2csPERI-CC2 Lagrangian. The indices i, j, k (a, b, c) refer to active occupied (virtual) 2c∆CC spinors. Again, Gpol ) is abbreviated by Gpol pq (D pq . Derivative Abbr. Polarization Contribution ∂L ˜ pol Ωai G ai ∂ t¯ia ∂L ηia 2Gpol ia ∂tia t¯=0 M ∑ ∂ 2L ˜ ⊤ Ruv θ v,jb ˜ pol δij − G ˜ pol δab − 1 Aia,bj G θ u,ai ji ab ∂tbj ∂ t¯ia 2 ∂ L ∂tai ∂tbj 2

3

Fai,bj

uv

M ∑

¯ pol ¯ pol θ⊤ u,ia Ruv θ v,jb − tib Gja − tja Gib uv ( ) M ∑ 1 ∑ ˆ ab ∑ ¯ ˜ ⊤ ⊤ ˜ Ruv θ v,ia − P tjc θ u,cb − t¯kb θ u,jk 2 uv ij c k −

Implementation

The PE-GHF method was implemented in the ridft module of the Turbomole program package 24 based on the RI-GHF implementation. 25 The electrostatic interactions are added ˆ prior to the SCF iterations. The to the transformed two-dimensional one-electron operator h electrostatic fields due to the Hartree–Fock density are updated in each SCF iteration. From these, the contribution to the Fock operator is calculated and added to both the αα and ββ blocks of the Fock matrix. The 2c-sPERI-CC2 method was implemented as an extension of the 2c-RI-CC2 method 20 in the ricc2 module. 26 In contrast to the original PERI-CC2 method, 7 the amplitude equations are decoupled from the multipliers due to the simplified ansatz. Therefore, no macro iterations are necessary to update the coupled-cluster density D ∆CCS . Following the original PERI-CC2 implementation, the additional polarization terms are computed by means of effective density matrices D eff via 1 eff ˜ pol G pq (D ) = − 2

M ∑

˜ ⊤ Ruv F el (D eff ) . θ u,pq v

uv


(17)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 26

With truncated coupled-cluster theory, problems occur in the vicinity of conical intersections of two excitations of the same irreducible representation. Instead of two nearly degenerate excitation energies, a pair of eigenvalues is obtained, which are each other’s complex conjugate. The corresponding eigenvectors approach linear dependence. 27 This may be the case if phosphorescence energies of an organic molecule in an inhomogeneous environment are computed. By the interaction with the asymmetric environment, the symmetry of the wavefunction is reduced to C1 . Consequently, the energetically nearly degenerate triplet components of the solute, which are of different irreducible representations in vacuum, all transform according to A1 in the polarizable-embedding computation. Based on the work by KÃűhn et al., 27 we have implemented an algorithm to correct both energies Em , En and eigenvectors C nS , C m S of two nearly degenerate states n and m, if both the difference of the real parts, |ℜ(Em − En )|, and the sum of the imaginary parts of their energies, |ℑ(Em + En )|, are below 10−5 Eh . This correction is applied in each iteration of the computation of the excitation vectors. To correct for the near linear dependence of the vectors, a corrected overlap Σ is assumed. It is defined as suggested in Reference 27 as ( Σ = Smax tanh

S

) (18)

Smax

if the energy difference λ = Em − En is "real", or more precisely if |ℑ(λ)| < |ℜ(λ)|, and according to

( Σ = Smax tanh

1 S Smax

) (19)

else. Both for real and imaginary energy differences, the corrected overlap depends on the overlap of the uncorrected eigenvectors S = |(C nS )† C m S |. To obtain vectors whose overlap is equal to Σ and which are most closely resembling the original eigenvectors, the vectors C nS , C m S are symmetrically orthogonalized in a first step. Then, they are linearly combined to fulfill the overlap condition Σ. The maximal overlap Smax was set to 0.2 as was suggested


Page 11 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


in Reference 27. The energies of the two states are set to

E± =

    1 (Em + En ) ±

√

   1 (Em + En ) ±

|λS| √ 1−S 2

2

2

|λ| 1−S 2

if λ is "real", (20) else.

Hence, in the case of a complex conjugated pair of energies, that is in the case of an imaginary splitting in the vicinity of conical intersections, the corrected energies are real quantities. For small values of S ≪ Smax , the corrections of both the energies and vectors are negligible.

4

Solvent dependent phosphorescence energies and frequencies

With our newly developed 2c-sPERI-CC2 method, the solvent effect on the electronic phosphorescence frequencies of 4H-pyran-4-thione, which are well described in vacuum by the 2c-RI-CC2 method, 20 was investigated for three different solvents: methanol (MeOH), 2propanol (2-PrOH) and tetrachloromethane (CCl4 ).

4.1

Molecular dynamics simulations

Since PE computations require explicit structures of the embedding environment, N pT molecular dynamics (MD) simulations (1.0 bar, 295.15 K) were performed modeling 4Hpyran-4-thione in its triplet state surrounded by solvent molecules. The three solutions were modeled by cubic boxes with average lengths of 23–24 Å with periodic boundary conditions containing a single 4H-pyran-4-thione molecule and 96–229 solvent molecules to reproduce the solvent’s macroscopic density. These calculations were performed with the QMSIM program 28 which is based on a quantum mechanically derived force field (QMDFF). 29 The QMDFF parameters had been obtained from PBE0-D3(BJ)/def2-QZVP 30–33 computations of the three free solvent molecules in their ground states and of 4H-pyran-4-thione in its



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

triplet state. A total time of 5 ns was simulated with integration steps of 0.5 fs. Every 0.25 ps, the coordinates were written to disk.

4.2

Generation of the embedding potential

The embedding potential due to the solvent is approximated by the potential caused by static dipole and quadrupole moments as well as induced dipole moments, that is, polarizability tensors. Starting from a supermolecular structure, each embedding solvent molecule is replaced by a single MM site at its center of mass. Instead of computing the multipole moments and polarizability tensor of each MM site, that is, of each embedding solvent molecule, they are computed only once for each solvent type, namely for a free solvent molecule in its equilibrium structure. The corresponding orbital-relaxed RI-CC2 computations 26,34,35 of the free molecules were performed in a center of nuclear mass coordinate system utilizing the aug-cc-pVQZ 36–38 basis set and neglecting correlation effects of core orbitals. The explicit parameters of different solvent molecules are given in the Supporting Information. For each MM site, these parameters are transformed according to the orientation of the free molecule relative to the solvent molecule in the supermolecular structure.

4.3

Details of the 2c-sPERI-CC2 computations

For each solvent, 2c-PE-GHF and 2c-sPERI-CC2 computations were performed based on the last twenty MD snapshots which had been written out. The same computational setup as for the vacuum computation 20 was chosen: the 2c-X2C decoupling scheme 39 was applied to obtain the one-electron spin–orbit operator. The def2-TZVPP basis sets and auxiliary basis sets 33 were used after decontraction of the p shells. Excitations originating from the 22 core orbitals were neglected. The only differences to the vacuum computations were the embedding potential and the explicit structure of 4H-pyran-4-thione which was directly taken from each MD snapshot in the present case. The embedding potentials were obtained as described 12 ACS Paragon Plus Environment

Page 12 of 26

Page 13 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


in Section 4.2 by replacing the solvent molecules of the MD snapshots with MM sites. Of the resulting MM sites, only those with distances less than 25 a0 to the center of mass of 4H-pyran-4-thione were kept. Due to the different particle densities of the three solutions, the resulting number of MM sites strongly depends on the solvent. In case of methanol, 161–172 MM sites were taken into account in the 2c-PE computations, 91–96 in case of 2-propanol, and only 60–67 in case of tetrachloromethane. To avoid unphysical polarization due to the induced-point-dipole model, 40,41 the exponential Thole damping model was applied. 42,43 Since these computations were based on an embedding potential with multipole moments k up to second order and an anisotropic polarizability tensor p, they are referred to as k2p2 following the nomenclature of Reference 22. Additionally, computations applying modified embedding potentials were performed to investigate the impact of the different contributions to the frequencies. In one set of computations, both electrostatic interactions and the anisotropy of the polarizability tensors were neglected for methanol and 2-propanol (k0p1) resulting in a potential similar to tetrachloromethane. In a third set of computations, no embedding potential was applied (k0p0), hence the differences in the results can directly be attributed to the differences in the structures of 4H-pyran-4-thione. For each combination of solvent and potential type, the results obtained for the twenty structures were averaged and the standard deviation was computed.

Table 2: Vertical 2c-sPERI-CC2 and 2c-CC2 emission energies (∆E in eV) of the state 3 A2 of 4H-pyran-4-thione in different environments and based on different embedding potentials. For each solvent, the average ± standard deviation of 20 separate computations is given. k0p0 k0p1 k2p2 2-PrOH 2.35 ± 0.00 2.51 ± 0.04 2.68 ± 0.24 MeOH 2.32 ± 0.00 2.45 ± 0.03 2.22 ± 0.17 CCl4 2.00 ± 0.00 2.18 ± 0.06 2.18 ± 0.06 vacuum 2.05



0.7 0.6 0.5

total geometrical polarization electrostatic

0.4 solvent shift / eV

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

0.3 0.2 0.1 0 -0.1 -0.2 -0.3 CCl4

MeOH

2-PrOH

Figure 1: Solvent shifts of averaged phosphorescence energies of 4H-pyran-4-thione and their respective geometrical δE geo , polarization δE pol and electrostatic as well as anisotropic δE elec contributions.

Table 3: Vertical 2c-sPERI-CC2 and 2c-CC2 phosphorescence frequencies (k in ms−1 ) of the state 3 A2 of 4H-pyran-4-thione in different environments. For each solvent, the average ± standard deviation of 20 separate computations is given.

2-PrOH MeOH CCl4 vacuum

k(⟨3 A2 ⟩) k(A1 ) k(B2 ) k(B1 ) 12 ± 6 39 ± 17 0.22 ± 0.35 0.051 ± 0.025 4.8 ± 2.3 15 ± 7 0.09 ± 0.25 0.033 ± 0.024 9.5 ± 2.5 32 ± 7 0.024 ± 0.010 0.021 ± 0.009 4.3 15 0.025 0.0033


Page 14 of 26

Page 15 of 26

5

vacuum 2-PrOH MeOH CCl4

4 oscillator strength / 10-3

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


3

2

1

0 1.5

2.0

2.5

3.0

3.5

E / eV

Figure 2: Phosphorescence emission spectra of 4H-pyran-4-thione in vacuum and in three different solutions: methanol (MeOH), 2-propanol (2-PrOH) and tetrachloromethane (CCl4 ). The spectra have been obtained from 60 bands per solution which each are broadened by Gaussian functions (σ 2 = 0.08 eV2 ) normalized to their oscillator strengths.

4.4

Phosphorescence energies

In Table 2, the averaged phosphorescence energies ∆E(⟨3 A2 ⟩) and their standard deviations obtained for the three solvents as well as for the vacuum case are listed. The averaged phosphorescence energies obtained by applying the full embedding potential are blue shifted relative to the vacuum result. The shifts due to the 2-propanol environment are by far the largest (0.63 eV). They are four times as large as the shift due to the methanol environment (0.17 eV). The smallest shift is obtained for tetrachloromethane (0.13 eV). To investigate the influence of the different types of solute–solvent interactions, the solvent shift δE relativ to the vacuum emission energy ∆E(vac) was defined as sum of three contributions,

δE = δE geo + δE pol + δE elec ,

(21)

which are due to the change of the solute’s geometry compared to the free molecule δE geo , the polarization contribution δE pol , and due to the electrostatic interactions and the anisotropy



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 26

of the polarizability tensor δE elec ,

δE geo = ∆E(k0p0) − ∆E(vac) ,

(22a)

δE pol = ∆E(k0p1) − ∆E(k0p0) ,

(22b)

δE elec = ∆E(k2p2) − ∆E(k0p1) .

(22c)

The last contribution vanishes in case of tetrachloromethane, since neither multipole moments nor anisotropic contributions to the polarizability tensor are present in the force field due to the Td symmetry of the reference MM solvent molecule. For all three solvents, these contributions are visualized in Figure 1. Pronounced geometrical contributions to the solvent shifts were obtained for both methanol and 2-propanol, which are of similar value (0.27 eV and 0.30 eV). The polarization contributions are all positive and increase with increasing isotropic averaged polarizability of the solvent, that is, the largest shift was obtained for tetrachloromethane (70.0 a30 ), then for 2-propanol (45.7 a30 ), and the smallest shift was obtained for methanol (21.4 a30 ). However, since the averaged densities of MM sites differ for the three solvent types, the polarization contribution to the solvent shift is not linearly proportional to the polarizability of the free solute. The electrostatic and anisotropic contributions to the solvent shifts significantly differ for methanol and 2-propanol. While it is of the same sign and order of magnitude as the polarization contribution in case of 2-propanol, it is negative in case of methanol and nearly twice as large in magnitude as the polarization contribution. A significant variation of the solvent shifts obtained from different MD snapshots and thus dependence of the energies on the explicit structure of the environment was only obtained when electrostatic interactions were included.


20

20

15

15

15

10

5

0

frequency / (ms)-1

20

frequency / (ms)-1

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


frequency / (ms)-1

Page 17 of 26

10

5

va CC Me 2cu l OH PrO um 4 H

0

(a) k0p0

10

5


(b) k0p1

0


(c) k2p2

Figure 3: Averaged phosphorescence frequencies of 4H-pyran-4-thione obtained for different types of embedding potentials. The standard deviations are given in terms of error bars.

4.5

Phosphorescence frequencies

The phosphorescence frequency of each triplet component,

kn =

2 (En − E0 )2 fn , c3

(23)

was computed from its oscillator strength fn and phosphorescence energy En − E0 . Here, c is the speed of light. From this, the total phosphorescence frequency k(⟨3 A2 ⟩) at 295.15 K was computed by weighting each component with its Boltzmann factor. For all three solvents, the averages and corresponding standard deviations obtained by taking into account the full potential are given in Table 3. Since the embedding potential reduces the symmetry of 4H-pyran-4-thione from C2v to C1 , all triplet components transform according to A1 . To distinguish between the states, they are assumed to originate from the same states as in the vacuum case giving rise to the nomenclature A1 , B2 and B1 with decreasing emission frequency. Again, the largest impact of the solvent is observed in the case of 2-propanol. Its total



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

frequency is three times as large as in the vacuum case, while the frequency is only doubled in the case of tetrachloromethane. Surprisingly, no significant solvent shift of the total phosphorescence frequency is observed for methanol. In case of both methanol and 2-propanol, the structure of the environment strongly influences the phosphorescence frequency. This becomes evident in the large relative standard deviation of approximately 50 %. Obviously, the phosphorescence process is always dominated by a single transition (A1 →A1 ) which exhibits the same environment dependence as the total frequency discussed before. In case of methanol and 2-propanol, the influence of the explicit structure of the environment on k(B1 ) and k(B2 ) is even stronger pronounced than on k(A1 ), resulting in standard deviations of more than 100 %. In contrast to the vacuum values, the frequencies k(B1 ) and k(B2 ) obtained for tetrachloromethane are similar in magnitude while their ratio is only of the order of 0.1 in the vacuum case. For both methanol and propanol, this ratio is 0.3–0.4. This increased ratio is – at least partly – attributed to the inability of truncated coupled-cluster methods to reproduce almost degenerate states of the same irreducible representation. These two emitting states are almost degenerate with an energy difference of the order of 10 µeV which leads to linear dependencies as pointed out in Reference 27. Even though the vectors are corrected to reproduce a reasonable overlap as described in Section 3, information still is lost for numerical reasons. Still, the resulting error of the total frequency is assumed to be low, since the total frequency is dominated by k(A1 ) which is at least two orders of magnitude larger than k(B2 ) in all three cases. The resulting emission spectra were plotted (Figure 2).

The separate effects of the polarization and the electrostatic interactions become evident in Figure 3. Here, the averaged phosphorescence frequencies obtained for the three types of potentials k0p0, k0p1 and k2p2 are visualized. The geometry of the solute has only a minor influence on the emission frequency (Figure 3(a)) while it has a significant impact on the solvent shift of the energy as pointed out before. If polarization interactions are taken into


Page 18 of 26

Page 19 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


account (Figure 3(b)), the averaged frequencies are doubled in all three cases. However, since the MM site density decreases with increasing polarization in the present example, similar frequencies are obtained for all three solvents. If additionally electrostatic interactions are included in the calculations (Figure 3(c)), the frequency is further increased by 5 % in case of 2-propanol, while it is reduced by 50 % in case of methanol. This effect is similar to the one observed for the solvent shifts of the emission energies. Again, the large variations of the separate frequencies obtained for different MD snapshots are due to electrostatic interactions and the anisotropic character of the polarizability tensor. These investigations illustrate the dependence of the transition moments on the embedding potential. Furthermore they indicate that for comparisons with experimental values more elaborate potentials should be used.

5

Conclusion

We have generalized the PE-HF and sPERI-CC2 methods to two-component wavefunctions. Since the resulting 2c-sPERI-CC2 method takes into account both spin–orbit coupling and interactions with the environment, for example solvent molecules, this new method opens up the opportunity to compute phosphorescence lifetimes of small to medium sized molecules in solution. The polarizable-embedding ansatz incorporates information on the explicit structure of the environment, therefore, effects due to the dynamics of the solvent can be investigated. With this new method, we have computed electronic phosphorescence energies and emission frequencies of 4H-pyran-4-thione in three different solutions: methanol, 2-propanol and tetrachloromethane. Our results indicate that the averaged electronic phosphorescence frequencies are both significantly influenced by the polarization of the environment and the interactions with the static charge distribution due to the solvent. In case of both methanol and 2-propanol, a strong dependence on the explicit structure of the surrounding solvent shell was observed.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 26

Acknowledgement K.K. and W.K. gratefully acknowledge support by the Deutsche Forschungsgemeinschaft (DFG) through SFB/Transregio 88 “3MET” (Project No. C1).

Supporting Information Available • Force-field parameters of methanol, 2-propanol, and tetrachloromethane; CCS-based energies, transition frequencies and shifts This material is available free of charge via the Internet at http://pubs.acs.org/.

A

Definition of the density in complex polarizable-embedding coupled-cluster theory

In Sec. 2.2, the explicit equations of the 2c-sPERI-CC2 method have been presented based on a complex, non-Hermitian coupled-cluster density. At first glance, it seems to be more reasonable to define a real density matrix which ensures that physical meaningful properties are obtained. In the present section, we therefore discuss this alternative definition of the coupled-cluster density matrix. For the sake of clarity, we refrain from applying the CC2 approximation in the following equations. For real wavefunctions, the density is defined as the Λ density, ∆Λ GHF Dpq = ⟨Λ|ˆ a†p a ˆq |CC⟩ − Dpq ,

(24)

which is used in the original PECC method, 22 or as the projected density, GHF ∆CC , ˆq |CC⟩ − Dpq = ⟨GHF|ˆ a†p a Dpq


(25)

Page 21 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


which yields the simplified PECC method. 21 However, in the case of complex wavefunctions, both matrices are non-Hermitian. Thus, the resulting fields may have imaginary components. This unphysical result is a direct consequence of the non-Hermitian coupled-cluster ansatz. To overcome this problem, one could symmetrize the density, ( ) ¯ = 1 D + D† . D 2

(26)

The real Lagrangian, here based on the simplified polarizable-embedding coupled-cluster method, then reads LrPE =E PE-GHF 1( ¯ ∆CC )|CC⟩ ˆ pol (D + ⟨Λ|FˆNPE + gˆN + G N 2 ) ¯ ∆CC )|Λ⟩ . ˆ pol (D + ⟨CC|Fˆ PE + gˆN + G N

(27)

N

The main drawback of the symmetrized density is the resulting coupling of amplitudes t and Lagrangian multipliers t¯ with their complex conjugated forms. In contrast to the Lagrangian of the QM system in vacuum, the Lagrangian LrPE (27) is not separable into two complex terms, which depend either only on t and t¯ or only on t∗ and t¯∗ . Due to this coupling, ”mixed” derivatives of the Lagrangian do not necessarily vanish. The most prominent example is the CC Jacobian which contains non-zero off-diagonal matrix elements B and B ∗ . Consequently, its eigenvalue equation reads 



 n



 n

 A B  C   C   ωn   = En B ∗ A∗ En


(28)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 26

with the excitation energy ωn and the matrix elements

Aµν

∂ 2 LrPE = = ∂tν ∂ t¯µ

Bµν

∂ 2 LrPE = = ∂t∗ν ∂ t¯µ

( (

∂ 2 LrPE ∂t∗ν ∂ t¯∗µ ∂ 2 LrPE ∂tν ∂ t¯∗µ

)∗ (29a) )∗ .

(29b)

µ, ν refer to arbitrary excitations. Since one possible solution of the eigenvalue equation is E n = (C n )∗ , it seems to be sufficient to solve the reduced eigenvalue equation, AC n + B(C n )∗ = C n ωn ,

(30)

which is similar to the procedure in the vacuum case where B vanishes. However, the scaling of Eq. (28) with an arbitrary phase factor eiφ gives a surprising result, 



 iφ

n



 iφ

n

 A B  e C   e C    =  ωn . B ∗ A∗ eiφ E n eiφ E n

(31)

The resulting eigenvector of the second row eiφ E n = eiφ (C n )∗ is not the complex conjugate of the first row eigenvector, eiφ C n . Hence, the result E n = (C n )∗ holds only in the special case of a real phase factor eiφ = ±1 – even in the absence of spin–orbit interactions. This is an artefact of the artificial symmetrization of the non-Hermitian coupled-cluster ansatz. To conclude, subsequently symmetrizing the density is not advisable in the context of polarizable-embedding coupled-cluster theory.

References (1) Sebej, P.; Lim, B. H.; Park, B. S.; Givens, R. S.; Klan, P. Org. Lett. 2011, 13, 644. (2) González, L.; Escudero, D.; Serrano-Andrés, L. ChemPhysChem 2012, 13, 28.


Page 23 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(3) Gomes, A. S. P.; Jacob, C. R. Annu. Rep. Prog. Chem., Sect. C: Phys. Chem. 2012, 108, 222. (4) Höfener, S.; Visscher, L. J. Chem. Phys. 2012, 137, 204120. (5) Wesolowski, T. A.; Warshel, A. J. Phys. Chem. 1993, 97, 8050. (6) Neugebauer, J. Phys. Rep. 2010, 489, 1. (7) Schwabe, T.; Sneskov, K.; Olsen, J. M. H.; Kongsted, J.; Christiansen, O.; Hättig, C. J. Chem. Theory Comput. 2012, 8, 3274. (8) Höfener, S. J. Comput. Chem. 2014, 35, 1716. (9) Lunkenheimer, B.; Köhn, A. J Chem. Theory Comput. 2013, 9, 977. (10) Cammi, R. Int. J. Quantum Chem. 2010, 110, 3040. (11) Cammi, R. Int. J. Quantum Chem. 2012, 112, 2547. (12) Paasche, A.; Schirmeister, T.; Engels, B. J Chem. Theory Comput. 2013, 9, 1765. (13) Hrsak, D.; Holmegaard, L.; Poulsen, A. S.; List, N. H.; Kongsted, J.; Denofrio, M. P.; Erra-Balsells, R.; Cabrerizo, F. M.; Christiansen, O.; Ogilby, P. R. Phys. Chem. Chem. Phys. 2015, 17, 12090. (14) Schwabe, T.; Beerepoot, M. T. P.; Olsen, J. M. H.; Kongsted, J. Phys. Chem. Chem. Phys. 2015, 17, 2582. (15) Georgieva, I.; Aquino, A. J. A.; Plasser, F.; Trendafilova, N.; Köhn, A.; Lischka, H. J. Phys. Chem. A 2015, 119, 6232. (16) Brisker-Klaiman, D.; Dreuw, A. Chem. Phys. Chem. 2015, 16, 1695. (17) Höfener, S.; Gomes, A. S. P.; Visscher, L. J. Chem. Phys. 2013, 139, 104106.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(18) Caricato, M. J. Chem. Theory Comput. 2012, 8, 4494. (19) Caricato, M. J. Chem. Phys. 2013, 139, 044116. (20) Krause, K.; Klopper, W. J, Chem. Phys. 2015, 142, 104109. (21) Krause, K.; Klopper, W. J. Chem. Phys. 2016, 144, 041101. (22) Sneskov, K.; Schwabe, T.; Kongsted, J.; Christiansen, O. J. Chem. Phys. 2011, 134, 104108. (23) Christiansen, O.; Koch, H.; Jørgensen, P. Chem. Phys. Lett. 1995, 243, 409. (24) Furche, F.; Ahlrichs, R.; Hättig, C.; Klopper, W.; Sierka, M.; Weigend, F. WIREs: Comput. Mol. Sci. 2014, 4, 91. (25) Armbruster, M. K.; Weigend, F.; van Wüllen, C.; Klopper, W. Phys. Chem. Chem. Phys. 2008, 10, 1748. (26) Hättig, C.; Weigend, F. J. Chem. Phys. 2000, 113, 5154. (27) Köhn, A.; Tajti, A. J. Chem. Phys. 2007, 127, 044105. (28) Shushkov, P.; Grimme, S. QMSIM. Universität Bonn, 2015. (29) Grimme, S. J. Chem. Theory Comput. 2014, 10, 4497. (30) Adamo, C.; Barone, V. J. Chem. Phys. 1999, 110, 6158. (31) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. J. Chem. Phys. 2010, 132, 154104. (32) Grimme, S.; Ehrlich, S.; Goerigk, L. J. Comput. Chem. 2011, 32, 1456. (33) Weigend, F.; Ahlrichs, R. Phys. Chem. Chem. Phys. 2005, 7, 3297. (34) Winter, N. O. C.; Hättig, C. Chem. Phys. 2012, 401, 217.


Page 24 of 26

Page 25 of 26

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(35) Friese, D. H.; Winter, N. O. C.; Balzerowski, P.; Schwan, R.; Hättig, C. J. Chem. Phys. 2012, 136, 174106. (36) Dunning, T. H. J. Chem. Phys. 1989, 90, 1007. (37) Kendall, R. A.; Dunning, T. H.; Harrison, R. J. J. Chem. Phys. 1992, 96, 6796. (38) Weigend, F.; Köhn, A.; Hättig, C. J. Chem. Phys. 2002, 116, 3175. (39) Peng, D.; Middendorf, N.; Weigend, F.; Reiher, M. J. Chem. Phys. 2013, 138, 184105. (40) Elking, D.; Darden, T.; Woods, R. J. J. Comput. Chem. 2007, 28, 1261. (41) van Duijnen, P. T.; Swart, M. J. Phys. Chem. A 1998, 102, 2399. (42) Thole, B. T. Chem. Phys. 1981, 59, 341. (43) van Duijnen, P. T.; ; Swart, M. J. Phys. Chem. A 1998, 102, 2399.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Table of Contents Graphic

Figure 4: For Table of Contents Only


Page 26 of 26

The Two-Component Polarizable-Embedding ... - ACS Publications

Recommend Documents