Description of spin–orbit coupling in excited states with two-component methods based on approximate coupled-cluster theory Katharina Krause and Wim Klopper Citation: The Journal of Chemical Physics 142, 104109 (2015); doi: 10.1063/1.4908536 View online: http://dx.doi.org/10.1063/1.4908536 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/142/10?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Two-component hybrid time-dependent density functional theory within the Tamm-Dancoff approximation J. Chem. Phys. 142, 034116 (2015); 10.1063/1.4905829 Communication: Two-component ring-coupled-cluster computation of the correlation energy in the randomphase approximation J. Chem. Phys. 139, 191102 (2013); 10.1063/1.4832738 Analytic second derivatives in closed-shell coupled-cluster theory with spin-orbit coupling J. Chem. Phys. 131, 164113 (2009); 10.1063/1.3245954 Closed-shell coupled-cluster theory with spin-orbit coupling J. Chem. Phys. 129, 064113 (2008); 10.1063/1.2968136 Relativistically corrected hyperfine structure constants calculated with the regular approximation applied to correlation corrected ab initio theory J. Chem. Phys. 121, 5618 (2004); 10.1063/1.1785772

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Wed, 19 Aug 2015 16:05:38

THE JOURNAL OF CHEMICAL PHYSICS 142, 104109 (2015)

Description of spin–orbit coupling in excited states with two-component methods based on approximate coupled-cluster theory Katharina Krause and Wim Kloppera) Karlsruhe Institute of Technology (KIT), Institute of Physical Chemistry, Theoretical Chemistry Group, KIT Campus South, P.O. Box 6980, 76049 Karlsruhe, Germany

(Received 5 December 2014; accepted 5 February 2015; published online 13 March 2015) A generalization of the approximated coupled-cluster singles and doubles method and the algebraic diagrammatic construction scheme up to second order to two-component spinors obtained from a relativistic Hartree–Fock calculation is reported. Computational results for zero-field splittings of atoms and monoatomic cations, triplet lifetimes of two organic molecules, and the spin-forbidden part of the UV/Vis absorption spectrum of tris(ethylenediamine)cobalt(III) are presented. C 2015 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4908536]

I. INTRODUCTION

Over the last almost 20 years, the approximated coupledcluster singles and doubles (CC2) method1 has become a popular method for the computation of properties of excited states. Nowadays, it is applied routinely to compute excitation energies of molecules consisting of “light” elements or oscillator strengths of spin-conserving transitions. For very recent applications, see, for example, Refs. 2 and 3. Compared to CCSD (coupled-cluster singles and doubles), double excitations are treated only approximately with CC2. As a consequence, the formal scaling with the system size N is reduced by one order of magnitude to O(N 5). Thus, it is still feasible for medium-sized molecules. Moreover, since most low-lying excited states are dominated by single excitations, CC2 yields similar results compared to CCSD while being superior to CCS (coupled-cluster singles) in terms of accuracy.4,5 The same applies to the closely related ADC(2) method (algebraic diagrammatic construction scheme up to second order) which is derived from a polarization propagator approach. For a recent review, see Ref. 6. In order to extend the field of possible applications to zero-field splittings, i.e., vertical excitation energies of molecules containing “heavy” elements, or spin-forbidden transitions, spin–orbit (SO) interactions have to be included. A straightforward way is to include the spin–orbit operator in the initial Hamiltonian. As a consequence, a wavefunction is obtained that includes spin–orbit effects intrinsically. A fundamental way would be to apply a four-dimensional Dirac operator. This has been realized, for example, by Gao et al. in the context of time-dependent density functional theory7 or very recently by Pernpointner in combination with the polarization propagator for the calculation of excited states.8 In principle, this four-component approach yields an accurate description of relativistic effects. It suffers, however, from two major drawbacks, namely, the presence of the negative-energy a)Author to whom correspondence should be addressed. Electronic mail:

[email protected]. Fax: +49 721 60847225.

0021-9606/2015/142(10)/104109/7/$30.00

solutions and the high computational demand.9 A commonly used cure is to apply a decoupling scheme to separate the large and small components and to restrict the resulting twodimensional spin–orbit operator to an effective one-electron operator. The corresponding wavefunction is then setup by two-component (2c) spinors. Especially for “heavy” elements, two-component effective core potentials (ECPs) provide an efficient alternative to a quasi-relativistic all-electron calculation. Since only the valence electrons are treated explicitly, the computational demand is reduced significantly compared to an all-electron calculation. However, for molecules containing both heavy and light elements, spin–orbit matrix elements are retained only for the heavy atoms of the molecule if ECPs are employed, while in case of an all-electron calculation, SO matrix elements for all atoms are considered. The thorough discussion of alternatives to including the spin–orbit operator in the initial Hamiltonian is beyond the scope of the present paper. However, we would like to mention at least some of them briefly. If the spin–orbit operator is treated as a perturbation to the non-relativistic or spin-free relativistic state, properties are available from response theory10,11 or from a subsequent spin–orbit configuration interaction (SOCI) computation.12 Wang et al. have demonstrated that, in the context of CCSD, it is sufficient to include spin–orbit coupling only in the post-Hartree–Fock computation.13,14 In the following, we present a generalization of the CC2 and the ADC(2) methods to two-component generalized Hartree–Fock (GHF) reference wavefunctions. Spin–orbit interactions are treated by an effective one-electron operator, which is two-dimensional, while the explicit electron-electron interaction is modeled by the non-relativistic Coulomb operator. Our implementation in the T program package15 supports either an all-electron treatment in terms of a decoupled Dirac–Coulomb (DC) operator or the use of 2cECPs. In Sec. II, the underlying theory is presented briefly. The important aspects of our implementation are highlighted in Sec. III. Exemplary calculations of zero-field splittings,

142, 104109-1

© 2015 AIP Publishing LLC

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Wed, 19 Aug 2015 16:05:38

104109-2

K. Krause and W. Klopper

J. Chem. Phys. 142, 104109 (2015)

triplet lifetimes, and spin-forbidden contributions to a UV/Vis absorption spectrum are present in Sec. IV.

response function ( ) X M0← η µX Cµ (ω n ) + ξ µX M¯ µ (ω n ) , n =

(11)

µ X M0→ n =

II. THEORY

ˆ 2c for In the following, a quasi-relativistic Hamiltonian H an N-electron system is assumed which contains an effective ˆ two-dimensional one-electron operator h(i) while treating the two-electron interactions by a non-relativistic Coulomb operator: ˆ 2c = H

N  i

N  1 ˆ + h(i) + VNN. r i> j i j

Nbf  σ

cσα/β p ∈ C,

|CC2⟩ = eT1+T2|GHF⟩. ˆ

(2)

(3)

The cluster amplitudes t µ are determined by zeroing the vector functions Ω µ , ˆ˜ Tˆ ]|GHF⟩, 0 = Ω µ 1 = ⟨µ1| Hˆ˜ + [ H, 2 ˆ ˆ Tˆ ]|GHF⟩, 0 = Ω = ⟨µ | H˜ + [F, 2

2

(4) (5)

where {⟨µ1|} ({⟨µ2|}) denotes the single (double) excitation manifold. The Tˆ1-transformation of the Hamiltonian is indiˆ ˆ Tˆ1 cated by the shorthand notation Oˆ˜ = e−T1Oe . While singles equation (4) is the same as for the CCSD method, doubles equation (5) has been approximated up to first order in the fluctuation potential. Following the methodology of coupledcluster ground state response theory, the frequencies ω n , which are related to the vertical excitation energies according to En − E0 = ~ω n

(in S.I. units),

(6)

are obtained as eigenvalues of the CC2 Jacobian A, ˆ˜ Tˆ ], τˆ ]|GHF⟩, Aµ 1ν1 = ⟨µ1|[ Hˆ˜ + [ H, 2 ν1 ˆ ˜ τˆ ]|GHF⟩, A = ⟨µ |[ H, µ 1ν 2

1

ν2

ˆ˜ τˆ ]|GHF⟩, Aµ 2ν1 = ⟨µ2|[ H, ν1 ˆ Aµ ν = ⟨µ2|[[F, Tˆ2], τˆν ]|GHF⟩. 2 2

2

ξ µX = ⟨µ|e−T2 Xˆ˜ eT2|GHF⟩,  ˆ ˆ˜ ˆ η µX = ⟨GHF| Xˆ | µ⟩ + t¯ν ⟨ν|e−T2[ X, τˆµ ]eT2|GHF⟩. ˆ

ˆ

(13) (14)

ν

Their computation requires the determination of the left ¯ n )A = ω n C(ω ¯ n) C(ω

(15)

and right eigenvectors

where χσ (r) denotes the σ-th basis function and cσα/β p denotes the corresponding expansion coefficient for spinor p, the 2c-CC2 method has been derived following the derivation presented in Ref. 1. The most important equations will be summarized below. The 2c-CC2 wavefunction is defined by applying a cluster  operator consisting of single Tˆ1 = µ 1 t µ 1τˆµ 1 and double Tˆ2  = µ 2 t µ 2τˆµ 2 excitation operators,

µ2

(12)

AC(ω n ) = ω n C(ω n ) of A which form a biorthonormal set according to  C¯ µ (ω m )Cµ (ω n ) = δ mn .

(16)

(17)

µ

cα χσ (r) * σβ p + , ,cσ p -

ˆ

ξ µX C¯ µ (ω n ),

µ

(1)

Here, the indices i and j refer to electrons, r i j is the nonrelativistic distance between the two electrons i and j, and VNN is the nuclear repulsion energy. Starting from a GHF wavefunction, which is setup by two-component complex spinors φ p (r) =



(7)

In addition, the solutions of the following eigenvalue problems are required: ¯tA = −η, Tˆ1+Tˆ2

ˆ τˆµ ]e η µ = ⟨GHF|[ H,

(18) |GHF⟩,

(19)

yielding the ground state Lagrange multipliers ¯t, and for each excited state n,

Fµν

¯ n )A = F( ¯t )C(ω n ), M(ω  ∂ Aκν , = ⟨GHF| Hˆ τˆµ τˆν |GHF⟩ + t¯κ κ ∂t µ

(20) (21)

¯ n ) which can be interpreted as excited state Lagiving M(ω grange multipliers.16 The algebraic diagrammatic construction (ADC) scheme up to second order is derived from a polarization propagator approach.17 It is closely related to the CC2 method, which becomes apparent if the final working equations are compared. The excitation energies are obtained as eigenvalues of an Hermitian matrix AADC(2) which corresponds to the symmetrized CC2 Jacobian, if the singles amplitudes t µ 1 are zeroed,18 1 ˆ ˆ ˆ AADC(2) = τˆν1]eT2|GHF⟩ ⟨µ1|e−T2[ H, µ 1ν 1 2  ˆ ˆ ˆ + ⟨ν1|e−T2[ H, τˆµ 1]eT2|GHF⟩† . (22) Since the ADC(2) Jacobian is Hermitian, the computation of transition moments is simplified significantly, if contributions that are second order in the ground state amplitudes are neglected, ∗ X X M0→ n = M0← n  ∗ ˆ ˆ = ⟨µ|e−T2 Xˆ eT2|GHF⟩ Cµ (ω n ) . (23) µ

(8) (9) (10)

The transition moments, which are associated with the position ˆ are obtained from the residues of the linear operator X,

III. IMPLEMENTATION

The 2c-RI-CC2 and 2c-RI-ADC(2) methods have been implemented in the T program package in the module 2. The implementation follows the existing program

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Wed, 19 Aug 2015 16:05:38

104109-3

K. Krause and W. Klopper

J. Chem. Phys. 142, 104109 (2015)

structure closely. It is based on similar intermediates which have been discussed in detail in Refs. 5 and 19 based on a restricted Hartree–Fock (RHF) reference. The most important aspects of our implementation and the major differences to the RHF and UHF variants of CC2 or ADC(2) will be highlighted in the following. A. Integral transformation

All four-index two-electron integrals are computed by an approximated resolution of the identity, (pq|r s) =

N aux 

B pq,Q BQ,r s ,

(24)

Q

reducing both the computational costs for the calculation and transformation of these integrals and the required disk space.5 While in the case of UHF, the B intermediates can be classified as α or β according to the two sets of MO-coefficients cα and c β , they contain contributions from both sets in case of a GHF reference, BQ,r s =

Nbf N aux   ρσ

P

α cρr −1/2 VQP (P| ρσ)* β + ,cρr -



α *cσβ s + . ,cσ s -

(25)

C. Linear dependence of degenerate CC2 eigenstates

Since A is non-Hermitian in the case of CC2, its eigenvectors are not orthogonal. Especially in the presence of degenerate excited states, the simultaneous iterative optimization of the eigenvectors might lead to linear dependencies. To overcome this problem, the trial vectors corresponding to (almost) degenerate states are orthogonalized against each other by a Schmidt orthogonalization scheme. Currently, states are treated as being degenerate if their energies differ by less than 10−6 hartree. D. Visualization of excitation processes

In a simplified picture, that is, if the contributions of both double excitations and electron correlation are neglected, a vertical excitation can solely be described by the singles excitation vector CS (ω n ) ≡ C. The difference density at the point r of the ground state and a state n is approximated by ( Nvir N occ  ρ n (r) = ℜ φa (r)†φb (r) Cia ∗Cib



ab

i

N occ 

Nvir 

Due to the approximation in the doubles equation and the utilization of canonical spinors, the doubles–doubles block of the CC Jacobian is a diagonal matrix. As a consequence, the inversion of (A DD − ω n 1)−1 is trivial and the doubles part of the excitation vectors (26)

can be computed from the corresponding singles part efficiently. In a similar way, the doubles part of a general system of equation L (A − λ1) = −ζ

(27)

can be easily computed according to

CS (ω n ) = On Wn Vn† .

L D (λ) = − (ζ D + LS AS D ) (A DD − λ1) .

(30)

=

N occ 

O ∗j I cσα/β j ,

(31)

α/β Vb I cσb ,

(32)

j vir α/β c˘σ I

=

Nvir  b

where the suffix n has been dropped for O, V, and c˘ . While NTSs are complex quantities which depend on both spin functions, the approximated difference densities are real onecomponent quantities by definition. Furthermore, they are distinct for each excitation that might be dominated by more than one NTS. If more than one transition contribute to one band, e.g., in the case of (nearly) degenerate states, an averaged density is easily computed from the weighted sum,

(28)

This leads to a reduction of the dimensions of the linear systems of equations to Nsingles × Nsingles and of the required disk space, since all doubles contributions can be calculated on the fly.

(29)

The singular values Wn correspond to the occupation numbers of the NTSs. The spinor coefficients result from the eigenvectors according to

ρ n1−n2(r) = −1

.

The first (second) term of this formula can be interpreted as the increase (decrease) of the occupation of the virtual (occupied) natural transition spinor (NTS). These NTSs are obtained from a singular-value decomposition of the singles excitation vector CS (ω n ),

occ α/β c˘σ I

B. Doubles amplitudes

) Cia ∗C ja

a

ij

As a second major difference, the spinor coefficients are complex, whereas they are real in case of a RHF or UHF reference. Thus, most intermediates which are computed in the 2c-RICC2 or 2c-RI-ADC(2) methods are complex as well. Since the derivation is based on a reference wavefunction consisting of canonical Hartree–Fock spinors, no explicit spin–orbit integrals have to be reevaluated in the correlated calculation. Instead, the SO interaction is included implicitly through the spinor energies.

C D (ω n ) = −(A DD − ω n 1)−1A DS CS (ω n )

φi (r) φ j (r) †

Nn1−n2 =

1

n2 

Nn1−n2

m=n 1

n2 

S 0→ m .

S 0→ m ρ m (r),

(33)

(34)

m=n 1

Here, the weights are defined as the relative transition strengths.

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Wed, 19 Aug 2015 16:05:38

104109-4

K. Krause and W. Klopper

J. Chem. Phys. 142, 104109 (2015)

IV. COMPUTATIONAL DETAILS AND RESULTS A. Zero field splittings

The first excited singlet and triplet states of Cd, Hg, Tl+, Ag , and Au+ were computed with our newly implemented 2cRI-ADC(2) and 2c-RI-CC2 methods in two ways, by including an effective core potential and by applying an all-electron relativistic operator. The first type of calculations, which are abbreviated by “2c-ECP” in Table I, utilized the two-component Dirac– Hartree–Fock effective core potentials25,26 (dhf-ecp-2c) taken from the T repository in combination with polarized quadruple-ζ valence basis sets (dhf-QZVPP-2c). Corresponding auxiliary basis sets were used. Additionally, spin-free relativistic calculations (1c-ECP) have been performed using the one-component variants of the aforementioned ECPs and basis sets. In the second type of calculations, the X2C decoupling scheme (2c-X2C) was applied to the Dirac-like one-electron operator in the underlying Hartree–Fock computation, while the electron-electron interactions were treated non-relativistically. Decontracted relativistic atomic-natural-orbital basis sets27,28 were used. For the underlying Hartree–Fock calculations, the same universal auxiliary basis set as in Ref. 29 was used for all atoms. The post-Hartree–Fock calculations were performed using auxiliary basis sets taken from the T repository (aug-cc-pwCV5Z-PP). In both types of +

calculations, the cores [Kr] and [Xe]4 f 14 were kept frozen for Cd/Ag and Au/Hg/Tl, respectively. The results are listed in Table I. For all calculations, there is no significant difference between the results obtained with the CC2 method and the ADC(2) method. On average, the excitation energies obtained with the CC2 method are 0.04 eV higher than with ADC(2). The deviations of the computed excitation energies of the averaged triplet states to the experimentally obtained ones are below 0.1 eV in case of the uncharged systems. While the averaged triplet excitation energies of Tl+, Ag+, and Au+ deviate by up to 0.11 eV in case of the 2c-X2C approach, the deviation is larger (up to 0.22) when 2c-ECPs are employed. The zero-field splittings are reasonably well described by all methods. The largest errors are encountered for the 3D1 states of Ag+ (0.14 eV) and Au+ (0.19 eV) in case of the 2c-X2C ansatz. In general, computations based on X2C yield larger zero-field splittings than the ECP-based calculations. The first vertical singlet excitation energies deviate in average by 0.15 eV in both cases with the largest errors being 0.23 eV (2c-ECP) and 0.37 eV (2c-X2C) both for the 1D2 state of Au+. Surprisingly, the 2c-X2C ansatz yields a better agreement with the experimental values for the singlet states for Cd, Hg, Tl+, and Ag+, while the 2c-ECP ansatz yields smaller deviations in case of Au+. The average excitation energies computed with the scalarrelativistic approach are, as expected, higher than the average

TABLE I. Vertical excitation energies and zero-field splittings relative to the average in eV. ADC(2) Atom

State 3P

CC2

1c-ECP

2c-ECP

2c-X2C

1c-ECP

2c-ECP

2c-X2C

Expt.

Cd

Average of Splitting of 3P0 Splitting of 3P1 Splitting of 3P2 1P 1

3.86 ... ... ... 5.57

3.85 −0.13 −0.07 0.07 5.58

3.82 −0.14 −0.07 0.07 5.50

3.91 ... ... ... 5.57

3.90 −0.13 −0.07 0.07 5.58

3.88 −0.14 −0.07 0.07 5.52

3.87 −0.14 −0.07 0.07 5.42

(reference 20)

Hg

Average of 3P Splitting of 3P0 Splitting of 3P1 Splitting of 3P2 1P 1

5.26 ... ... ... 6.80

5.23 −0.48 −0.27 0.26 6.84

5.22 −0.50 −0.29 0.27 6.82

5.31 ... ... ... 6.82

5.27 −0.48 −0.27 0.26 6.85

5.27 −0.50 −0.29 0.27 6.84

5.18 −0.51 −0.29 0.28 6.70

(reference 21)

Tl+

Average of 3P Splitting of 3P0 Splitting of 3P1 Splitting of 3P2 1P 1

7.31 ... ... ... 9.42

7.21 −0.89 −0.55 0.51 9.53

7.14 −0.95 −0.59 0.54 9.50

7.37 ... ... ... 9.44

7.25 −0.89 −0.55 0.51 9.55

7.17 −0.96 −0.59 0.55 9.51

7.10 −0.97 −0.60 0.55 9.38

(reference 22)

Ag+

Average of 3D Splitting of 3D3 Splitting of 3D2 Splitting of 3D1 1D 2

5.18 ... ... ... 5.63

5.16 −0.16 0.02 0.35 5.79

4.93 −0.22 −0.01 0.52 5.68

5.20 ... ... ... 5.66

5.18 −0.16 0.02 0.35 5.81

5.03 −0.22 −0.01 0.52 5.78

5.03 −0.18 0.02 0.39 5.71

(reference 23)

Au+

Average of 3D Splitting of 3D3 Splitting of 3D2 Splitting of 3D1 1D 2

2.63 ... ... ... 3.31

2.47 −0.41 −0.06 1.06 3.86

2.33 −0.49 −0.12 1.34 3.98

2.67 ... ... ... 3.36

2.51 −0.41 −0.06 1.06 3.90

2.39 −0.49 −0.12 1.35 4.04

2.29 −0.42 −0.10 1.15 3.67

(reference 24)

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Wed, 19 Aug 2015 16:05:38

104109-5

K. Krause and W. Klopper

J. Chem. Phys. 142, 104109 (2015)

FIG. 1. Approximated difference density of the averaged 3B1 → 1A1 vertical deexcitation of dithiosuccinimide on the 2c-RI-CC2 level of theory. Yellow corresponds to a gain while red corresponds to loss of electron density (iso-value: 0.002 a 0−3).

FIG. 2. Approximated difference density of the averaged 3A2 → 1A1 vertical deexcitation of 4H -pyran-4-thione on the 2c-RI-CC2 level of theory. Yellow corresponds to a gain while red corresponds to a loss of electron density (iso-value: 0.002 a 0−3).

obtained with the 2c-variant in case of triplet states and lower than the 2c-equivalent in case of the singlet states. In the scalar-relativistic ansatz, spin is a good quantum number, therefore, the states 3P1 and 1P1 or 3D2 and 1D2 do not interact. If spin–orbit coupling is taken into account, those states are of the same symmetry, namely, P1 or D2. Their interaction leads to the stabilization of the state which is lower in energy and destabilization of the higher state.

Here, c denotes the speed of light, me denotes the mass of an electron, and a0 denotes the Bohr radius. The oscillator strength f n is a dimensionless quantity that depends on the excitation frequency and the transition strengths SI0→ n ,  2 f n = ωn S 0→ n . (36) 3 I =x, y, z I Note that unlike the transition moments MI0 ← n and MI0 → n , the transition strength SI0 → n is a real quantity and a physical observable,16 ∗ 1  0← n 0→ n MI MI + MI0← n MI0→ n . SI0→ n = (37) 2 The results of our calculations are compared to oscillator strengths (see Tables II and III) calculated from the transition frequencies and transition moments that Kleinschmidt et al. have obtained with their MRSOCI method,12 where a SOCI matrix is setup starting from spin-free density functional theory/multireference configuration interaction (DFT/MRCI) wavefunctions. In their approach, the spin–orbit interaction is approximated by an atomic mean field derived from the Breit–Pauli Hamiltonian. Since the zero-field splitting is insignificant for the two systems (0.3 meV and 3.4 meV) and the experiments have not been performed at temperatures close to 0 K, it is assumed that the three components of the triplet states are populated equally. Therefore, the high-temperature average

B. Phosphorescence lifetimes

As an example of spin–forbidden transitions, phosphorescence frequencies kT1 of two organic molecules, dithiosuccinimide and 4H-pyran-4-thione, were calculated using the 2cRI-ADC(2) and 2c-RI-CC2 methods in combination with the 2c-X2C decoupling scheme.30 The dhf-TZVPP-2c basis sets with decontracted p shells were used. In the post-Hartree–Fock calculations, the cores [He] and [Ne] were kept frozen for C/N/O and S, respectively. The geometries were optimized for the first excited triplet state using the 1c-RI-ADC(2) and 1c-RI-CC2 methods, respectively, while utilizing the def2TZVPP basis sets and including the spin-free X2C decoupled one-electron operator. The averaged approximated difference densities based on the 2c-RI-CC2 single excitation vectors are shown in Figures 1 and 2. Clearly, both triplet states are of (n → π ∗) character. From the excitation energies of the three components of the first triplet state and the corresponding oscillator strengths f n , the phosphorescence frequencies were calculated according to 2 k n = 2 3 (En − E0)2 f n m e a0 c

(in S.I. units).

1 ⟨τ1⟩ = 3 + + =* ,3 is compared to the experimental values. (

(35)

−1 τ1,1

−1 τ1,2

) −1 −1 τ1,3

−1

 i

kT1, i + -

(38)

TABLE II. Triplet lifetimes in ms and oscillator strengths of the three components A1, A2, and B2 of the lowest triplet state 3B1 of dithiosuccinimide (experimental value: 0.2 ms12). MRSOCI12 State 3B 1 3B 1 3B 1

(A2) (B2) (A1)

⟨3B1⟩

τ 3224 0.20 1.18 0.52

2c-ADC(2) f × 10−9

1.3 2.1 × 10−5 3.6 × 10−6

τ 1149 0.16 1.00 0.38

2c-CC2 f × 10−9

3.5 2.6 × 10−5 4.0 × 10−6

τ

f

760 0.15 1.01

4.6 × 10−9 2.4 × 10−5 3.5 × 10−6

0.38

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Wed, 19 Aug 2015 16:05:38

104109-6

K. Krause and W. Klopper

J. Chem. Phys. 142, 104109 (2015)

TABLE III. Triplet lifetimes in ms and oscillator strengths of the three components B1, B2, and A1 of the lowest triplet state 3A2 of 4H -pyran-4-thione (experimental value: 0.1 ms12). MRSOCI12 State 3A 2 3A 2 3A 2

(B1) (B2) (A1)

⟨3A2⟩

2c-ADC(2)

2c-CC2

τ

f

τ

f

τ

f

140 27 0.11

5.1 × 10−8 2.6 × 10−7 6.5 × 10−5

321 50 0.08

2.2 × 10−8 1.4 × 10−7 8.3 × 10−5

304 39 0.07

1.8 × 10−8 1.4 × 10−7 8.0 × 10−5

0.33

0.25

The optimized structure of dithiosuccinimide is of C2 symmetry (see Fig. 1), however, in accordance with early works of Tatchen et al.,31 the deviations from C2v symmetry are small. Therefore, the nomenclature of the triplet sublevels in Table II assumes C2v symmetry for the sake of clarity. The term symbols of the states in the non-relativistic limit are given as well as the irreducible presentation of the relativistic wavefunctions in parentheses. The first transition, 3 B1 (A2) → 1A1 (A1), is dipole forbidden in C2v symmetry and gains its small transition strength, which is parallel to the z-axis, from the distortion of the structure towards C2. Indeed, in C2 symmetry, this state becomes 3B (A), which is dipole allowed in z-direction. With all three methods, the oscillator strength was computed to be below 5 × 10−9. Thus, the contribution of this transition to the averaged phosphorescence frequency is negligible. The deviations of the three results are mostly due to minor differences in the structure. For the second transition, 3B1 (B2) → 1A1 (A1), an oscillator strength of about 2 × 10−5 was obtained with all of the three methods. This contribution dominates the average phosphorescence time. The oscillator strengths for the third transition, 3B1 (A1) → 1A1 (A1), are in good agreement as well. They result in phosphorescence times of 1.00–1.18 ms. The high temperature averaged phosphorescence lifetimes obtained with 2c-RI-ADC(2) and 2c-RI-CC2 are slightly smaller than the one computed by Kleinschmidt and coworkers. Still, these three values are all of the same order of magnitude as the experimental value. 4H-pyran-4-thione exhibits C2v symmetry in the triplet state 3A2. The oscillator strengths of the first transition, 3 A2 (B1) → 1A1 (A1), are obtained as 2.2 × 10−8 and 1.8 × 10−8 with 2c-RI-CC2 and 2c-RI-ADC(2). The MRSOCI calculation yielded an oscillator strength of 5.1 × 10−8 which is twice as large. Still, since these oscillator strengths are very small, the

0.21

deviations are insignificant. With all three methods, oscillator strengths for the second and third transition were computed to be of the order of magnitude of 10−7 and 10−5, respectively. Hence, all three methods yielded similar phosphorescence lifetimes. Consequently, the computed averaged phosphorescence lifetimes agree well. They are of the same order of magnitude as the experimental result. C. Spin-forbidden contributions to UV/Vis spectra

To demonstrate that our newly implemented 2c-RI methods are applicable to larger molecules as well, the spinforbidden contributions to the low-energy part of the UV/Vis absorption spectrum of the tris(diethylenediamine)cobalt(III) cation were calculated. Since these states mainly reduce to triplet excited states in the SO-free picture, the energy levels are compared to those of triplet excited states obtained from a scalar-relativistic calculation. The geometry of the complex of D3 symmetry was optimized on the scalar-relativistic RI-MP2 level of theory, by applying the 1c-X2C-Hamiltonian. The computations utilized basis sets of polarized valence triple-ζ quality in case of cobalt, nitrogen, and carbon and of split-valence quality in case of hydrogen, all taken from the T repository. The same basis sets were used for the computation of vertical excitation energies of triplet states using the 1c-RI-ADC(2) method. For the two-component computation, i.e., the computation of the vertical excitation energies and oscillator strengths using the 2c-RI-ADC(2) method in combination with the

TABLE IV. Vertical excitation energies (ω in cm−1) and transition strengths ( f ) obtained with the 2c-RI-ADC(2) method and corresponding vertical triplet excitation energies obtained with 1c-RI-ADC(2). State

ω

State

ω

f x+y

fz

1 3E

16 985.2

1E 2E 2 A1 1 A2

16 458.8 16 523.5 16 569.1 17 406.3

3.3 × 10−6 1.9 × 10−6 0 0

0 0 0 7.8 × 10−7

1 3A2

17 150.0

3E 3 A1

17 532.8 18 405.4

7.0 × 10−8 0

0 0

FIG. 3. Low-energy regions of UV/Vis spectra of Co(en)3+ computed with 3 the 2c-RI-ADC(2) method. The arrows indicate the levels of triplet states obtained with the 1c-RI-ADC(2) method. Arbitrarily scaled Gaussians are used to broaden the spectra.

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Wed, 19 Aug 2015 16:05:38

104109-7

K. Krause and W. Klopper

J. Chem. Phys. 142, 104109 (2015)

while utilizing an approximated resolution of the identity in combination with the on-the-fly calculation of doubles. First applications showed no significant difference between the results obtained with CC2 and ADC(2). The computed excitation energies and the resulting zero-field splittings are in reasonable agreement with experiment. The phosphorescence lifetimes computed exemplarily for two organic molecules agree well with both quantum chemical and experimental reference values. The applicability of the method to larger molecules was demonstrated by computing the spin-forbidden contributions to the low-energy region of the UV/Vis spectrum of tris(diethylenediamine)cobalt(III).

ACKNOWLEDGMENTS

FIG. 4. Approximated difference density of the computed 2c-RI-ADC(2) vertical excitations of Co(en)3+ . Yellow corresponds to a loss while red 3 corresponds to a gain of electron density (iso-value: 0.002 a 0−3).

K.K. gratefully acknowledges support by the Fonds der Chemischen Industrie through a Kekulé Mobility Fellowship. We thank Professor Christof Hättig (Bochum) for many helpful discussions. 1O.

2c-X2C-Hamiltonian, the corresponding basis sets optimized for two-component Dirac–Hartree–Fock calculations were utilized. The cores [He] and [Ar] were kept frozen for C/N and Co for all post-Hartree–Fock calculations, resulting in 84 active occupied and 842 virtual spinors. The auxiliary basis set used in the post-Hartree–Fock calculations consisted of 1464 partially contracted basis functions. The scalar-relativistic computation yielded a state of E symmetry and one of A2 symmetry as the lowest triplet states. From group theory, it follows that the 3E state splits into two states of E symmetry, one of A1 and one of A2 symmetry, and that the 3A2 state splits into one of E and one of A1 symmetry if the spin and spatial functions couple. Even though the computation including spin–orbit interaction was performed without considering the symmetry of the complex, the irreducible representations of the final states can be assigned easily based on their degeneracy and the value of the z-component of the transition strength. They agree nicely with the predicted splittings. These excitation energies and transition strengths of the first nine vertical excitations of the Co(en)3+ 3 cation are listed in Table IV. The resulting UV/Vis absorption spectra are visualized in Figure 3. In Figure 4, an approximated difference density is shown, which has been computed as the weighted average of the computed transitions of Table IV as was described in Eq. (33).

V. CONCLUSION

The approximated coupled-cluster singles and doubles method as well as the algebraic diagrammatic construction scheme up to second order has been extended to complex twocomponent reference wavefunctions. These two new methods have been implemented in the T program package

Christiansen, H. Koch, and P. Jørgensen, Chem. Phys. Lett. 243, 409 (1995). 2A. Csehi, G. J. Halász, and Á. Vibók, Int. J. Quantum Chem. 114, 1135 (2014). 3S. Lobsiger, M. A. Trachsel, T. Den, and S. Leutwyler, J. Phys. Chem. B 118, 2973 (2014). 4D. Kánnár and P. G. Szalay, J. Chem. Theory Comput. 10, 3757 (2014). 5C. Hättig and F. Weigend, J. Chem. Phys. 113, 5154 (2000). 6A. Dreuw and M. Wormit, WIREs Comput. Mol. Sci. 5, 82 (2015). 7J. Gao, W. Liu, B. Song, and C. Liu, J. Chem. Phys. 121, 6658 (2004). 8M. Pernpointner, J. Chem. Phys. 140, 084108 (2014). 9D. Peng and M. Reiher, Theor. Chem. Acc. 131, 1081 (2012). 10O. Vahtras, H. Ågren, P. Jørgensen, H.-J. A. Jensen, T. Helgaker, and J. Olsen, J. Chem. Phys. 97, 9178 (1992). 11J. M. Younker and K. D. Dobbs, J. Phys. Chem. C 117, 25714 (2013). 12M. Kleinschmidt, J. Tatchen, and C. M. Marian, J. Chem. Phys. 124, 124101 (2006). 13F. Wang, J. Gauss, and C. van Wüllen, J. Chem. Phys. 129, 064113 (2008). 14Z. Wang, Z. Tu, and F. Wang, J. Chem. Theory Comput. 10, 5567 (2014). 15F. Furche, R. Ahlrichs, C. Hättig, W. Klopper, M. Sierka, and F. Weigend, WIREs Comput. Mol. Sci. 4, 91 (2014). 16C. Hättig, O. Christiansen, and P. Jørgensen, J. Chem. Phys. 108, 8331 (1998). 17J. Schirmer, Phys. Rev. A 26, 2395 (1982). 18C. Hättig, Adv. Quantum Chem. 50, 37 (2005). 19C. Hättig and A. Köhn, J. Chem. Phys. 117, 6939 (2002). 20K. Burns and K. B. Adams, J. Opt. Soc. Am. 46, 94 (1956). 21E. B. Saloman, J. Phys. Chem. Ref. Data 35, 1519 (2006). 22S. Johansson, G. Kalus, T. Brage, D. S. Leckrone, and G. M. Wahlgren, Astrophys. J. 462, 943 (1996). 23A. Kramida, J. Res. Natl. Inst. Stand. Technol. 118, 168 (2013). 24M. Rosberg and J.-F. Wyart, Phys. Scr. 55, 690 (1997). 25D. Figgen, G. Rauhut, M. Dolg, and H. Stoll, Chem. Phys. 311, 227 (2005). 26B. Metz, M. Schweizer, H. Stoll, M. Dolg, and W. Liu, Theor. Chem. Acc. 104, 22 (2000). 27B. O. Roos, R. Lindh, P.-A. Malmqvist, V. Veryazov, and P.-O. Widmark, J. Phys. Chem. A 109, 6575 (2005). 28B. O. Roos, R. Lindh, P.-A. Malmqvist, V. Veryazov, and P.-O. Widmark, J. Phys. Chem. A 108, 2851 (2005). 29K. Krause and W. Klopper, J. Chem. Phys. 139, 191102 (2013). 30D. Peng, N. Middendorf, F. Weigend, and M. Reiher, J. Chem. Phys. 138, 184105 (2013). 31J. Tatchen, M. Kleinschmidt, C. M. Marian, M. Parac, and S. Grimme, Z. Phys. Chem. 217, 205 (2003).

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.170.6.51 On: Wed, 19 Aug 2015 16:05:38

Description of spin-orbit coupling in excited states with two-component methods based on approximate coupled-cluster theory.

A generalization of the approximated coupled-cluster singles and doubles method and the algebraic diagrammatic construction scheme up to second order ...
2MB Sizes 0 Downloads 10 Views