Self-Consistent Constricted Variational Theory RSCF-CV(∞)-DFT

In this paper, the relaxed self-consistent field infinite order constricted variational density functional theory (RSCF-CV(∞)-DFT) for triplet calcu...
1 downloads 10 Views 843KB Size
Subscriber access provided by CORNELL UNIVERSITY LIBRARY

Article

Self-Consistent Constricted Variational Theory RSCF-CV(#)DFT and its Restrictions to Obtain a Numerically Stable #SCFDFT-like Method: Theory and Calculations for Triplet States Young Choon Park, Florian Senn, Mykhaylo Krykunov, and Tom Ziegler J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.6b00333 • Publication Date (Web): 01 Sep 2016 Downloaded from http://pubs.acs.org on September 6, 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 50

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

Self-Consistent Constricted Variational Theory RSCF-CV(∞)-DFT and its Restrictions to Obtain a Numerically Stable ∆SCF-DFT-like Method: Theory and Calculations for Triplet States Young Choon Park,†,¶ Florian Senn,†,¶ Mykhaylo Krykunov,∗,‡ and Tom ˆ

Ziegler†,S †Department of Chemistry, University of Calgary, 2500 University Drive NW, Calgary AB T2N-1N4, Canada ‡Department of Chemistry and Biomolecular Sciences, University of Ottawa, 75 Laurier Ave E, Ottawa, ON K1N-6N5, Canada ¶Y. C. Park and F. Senn contributed equally to this work. ˆ STom Ziegler passed away in March, 2015. E-mail: [email protected]

Abstract In this paper, the relaxed self-consistent field infinite order constricted variational density functional theory (RSCF-CV(∞)-DFT) for triplet calculations is presented. Here, we focus on two main features of our implementation. First, as an extension of our previous work by Krykunov and Ziegler (J. Chem. Theor. Comput. 9, 2761 (2013)), the optimization of the transition matrix representing the orbital transition is implemented and applied for vertical triplet excitations. Second, restricting the

1

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

transition matrix we introduce RSCF-CV(∞)-DFT-based numerically stable ∆SCFDFT-like methods, the most general one of them being SVD-RSCF-CV(∞)-DFT. The reliability of the different methods, RSCF-CV(∞)-DFT and its restricted versions, is examined using the benchmark test set of Silva-Junior et al. (J. Chem. Phys. 129, 104103 (2008)). The obtained excitation energies validate our approach and implementation for RSCF-CV(∞)-DFT and also show that SVD-RSCF-CV(∞)-DFT mimics very well ∆SCF-DFT, as the root mean square deviations between these methods are less than 0.1 eV for all functionals examined.

1

Introduction

For the calculation of electronically excited states, linear response adiabatic Time-Dependent Density Functional Theory 1–7 (TD-DFT) has become a ‘workhorse’ 8 and ‘delivers an excellent compromise between computational efficiency and accuracy’ 9 . But it is known that TD-DFT with standard functionals has difficulties with the description of certain types of excitation, e.g. charge-transfer 10–12 and Rydberg 13–15 excitations. A possible way to circumvent (at least some of) these problems is the use of long-range corrected hybrid functionals by fitting the range separation parameter to the value which equates I and A to −εH and −εL , correspondingly 16 . However, as it was pointed out by Baerends et al. 17 this feature deviates from the exact Kohn-Sham model, since εH = −I holds in exact theory, but there is no similar relation between εL and A. As a result, the exact Kohn-Sham HOMO-LUMO gap is not a good approximation to I −A. Thus, by mixing in the exact Hartree-Fock exchange, the problem of the long-range charge transfer is solved, but “also endow the virtual orbitals with the undesirable property that they are too diffuse, and therefore less suitable to describe excited states” 17 . Another challenge for TD-DFT, namely, single electron excitations with some double transition character, because this problem requires going beyond the adiabatic approximation 18–24 . The frequency dependent (non-adiabatic) TD-DFT improves significantly the energies of such excited states 18–20,25,26 , where a priori knowledge of the doubly excited 2

ACS Paragon Plus Environment

Page 2 of 50

Page 3 of 50

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

states energies is required in the original approach 18,19 but not needed anymore with the recent improvement by Huix-Rotllant et al. 26 . It is therefore not surprising that one looks for alternatives in such cases and there are different DFT-based methods, e.g. Self-Consistent Field DFT (∆SCF-DFT) 27–29 with extensions 30–32 , Ensemble-DFT 33–36 , Constrained Orthogonality Method (COM) 37–39 , Restricted Open-Shell Kohn-Sham (ROKS) 35,40,41 , Constrained DFT (CDFT) 42 , “Taking Orthogonality Constraints Into Account” (TOCIA) 43,44 , Maximum Overlap Method (MOM) 45,46 , Constricted Variational Density Functional Theory (CV-DFT) 47 and extensions 48–51 , Orthogonality Constrained DFT (OCDFT) 52 , Guided SCF 53 , etc., each of them having its advantages and assets. One of these mentioned methods is CV-DFT, reviewed in Ref.

54

; a DFT-based method

where the excited state orbitals are constructed from the ground state orbitals through a perturbation approach over a transition matrix. The correspondence of CV-DFT to second order, CV(2)-DFT, and TD-DFT has been shown theoretically and numerically, see e.g. Ref. 54 , while going to higher orders CV-DFT has shown its strength in theory as well as several applications categorized as e.g. charge-transfer 55–57 and Rydberg excitations 15 . The general RSCF-CV(∞)-DFT method, where the transition matrix is optimized and orbital relaxation is included, has been implemented only for singlet excitations 50 . For triplet excitations the orbital relaxation has been included recently 51 . But the optimization of the transition matrix, whose advantage has been shown for singlets in Ref. 50,57 , has until now not been implemented and therefore never been validated for triplet excitations. We would like to add this missing piece herewith, giving the theoretical expressions for the excitation energy in Section 2.2 as well as showing numerical applications in Section 4.

The oldest and probably simplest of these previously mentioned DFT-based methods to calculate excited states is ∆SCF-DFT 27–29 with its extensions 30–32 , a scheme in which the excited state is obtained by promoting an electron from an occupied to one of the virtual Kohn-Sham orbitals (or a linear combination of them) with a subsequent energy minimiza-

3

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

tion of this new configuration. In spite of its conceptual simplicity it has been successfully applied to challenging tasks as charge-transfer 10–12 and Rydberg 13–15 excitations. Unlike TD-DFT, ∆SCF-DFT includes optimization of the Kohn-Sham orbitals for a particular excited state. Recently the necessity of inclusion of density relaxation effects has been shown by Ronca et al. 58 , where these authors compared the (unrelaxed) response densities from TD-DFT with the corresponding relaxed densities obtained from highly accurate wave function calculations. Unfortunately, the price for the optimization of orbitals in ∆SCF-DFT is a possibility of the variational collapse to the lower states of the same symmetry. Further, the variational principle applied to excited states requires the constraint that the excited state of interest must be orthogonal to all states below it 59 . To solve this problem, various methods use different techniques have been developed to enforce the orthogonality for the excited states Slater determinants, e.g. COM 37–39 introduces orthogonality constraints for the separation of unitary transformations on both the occupied and unoccupied manifold, TOCIA 43,44,60 evaluates of the overlap matrix elements between orbitals of the ground and excited states to project out the lower lying state, OCVDFT 52 uses singular value decomposition technique. Other methods such as guided SCF 53 for excited states do not explicitly enforce the orthogonality. In guided SCF, each SCF cycle is split into several steps in order to maintain the desired orbital occupation arrangement, and thus avoid the variational collapse. Apart from the development of practical algorithms, there is also extensive research devoted to theoretical justification of the variational approach 59,61–66 and ∆SCF-DFT in particular 51,52,67 , since the original Hohenberg-Kohn theorem is valid only for the ground state 68,69 . However, a number of publications 51,52,67 conclude in one way or another that ∆SCF-DFT is as well justified as adiabatic TD-DFT. In a previous publication, Park et al. 51 showed how the ∆SCF-DFT excitation energies can be obtained within the RSCF-CV(∞)-DFT scheme for the single canonical orbital replacement (SOR) type of transition, which has already been claimed by theoretical arguments in Ref. 70 . In this publication we will show how the restriction on the transition type can be un-

4

ACS Paragon Plus Environment

Page 4 of 50

Page 5 of 50

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

stiffened to the much looser restriction of a single natural type orbital 71 replacement. Further we would like to present how one can implement these restrictions leading to SVD-RSCFCV(∞)-DFT (which we will refer to as SVD indicating the singular value decomposition of the transition matrix U) in Section 2.3 to mimic ∆SCF-DFT and show the numerical equivalence of the two methods in Section 4. This SVD method mimics ∆SCF-DFT but without the problem of a possible variational collapse and serves as a bridge, helping to understand the role of different excitations parts but also demonstrating that the electronic structure of the ∆SCF-DFT excitation described by a single Kohn-Sham determinant is much more complex than just a single canonical orbital replacement plus orbital relaxation. The combination of the variational approach of CV-DFT and ∆SCF-DFT also builds on previous investigations of Ref. 47,48,51,70 .

2

Theory

Before we present the theory of our work, the RSCF-CV(∞)-DFT scheme for triplet excitations in section 2.2 and ∆SCF-DFT-like methods within CV-DFT in 2.3, we will first outline the basis of CV-DFT, where we introduce quantities used later.

2.1

The CV-DFT scheme

In this subsection the general CV-DFT scheme will be repeated in a nutshell, for further explanations we refer to the review of CV-DFT, Ref. 54 . We shall consider here only excitations from a closed shell ground state described by a single Slater determinant Ψ0 = |φ1 φ2 . . . φi . . . φnocc |, where nocc is the number of occupied orbitals (nvir the number of virtual orbitals) and {φi } are the ground state orbitals. The general idea is to construct a new set of ‘occupied’ excited state orbitals by mixing the occupied and virtual ground state orbitals with a unitary transformation, Y. Then, the excitation energy is obtained by variational optimization of these orbitals in the excited state

5

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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 50

Slater determinant, Ψ′ = |φ′1 φ′2 ...φ′i ...φ′nocc |, where {φ′i } are the excited state orbitals. However, the matrix elements of Y cannot be used as independent variational parameters, since the unitary relation couples the elements of Y. Therefore, CV-DFT employs an exponential mapping. In the most general form this mapping is given by Y = exp(iH), where H is a Hermitian matrix 72 . Since we deal with real orbitals, only the real traceless antisymmetric contribution (or skew symmetric), U, is retained in the CV-DFT exponential expansion: ∞ X Uk

Y = exp(U) =

k=0

k!

.

(1)

Generally, U would contain ntot (ntot − 1)/2 matrix elements, where ntot = nocc + nvir . However, in the case of single determinant theories such as Hartree-Fock and Kohn-Sham, parameters that involve rotations within the occupied or virtual subspaces are redundant 73 . Therefore, U in Eq. (1) contains only the matrix elements between the occupied and virtual orbitals, and it is termed the transition matrix (which is also referred to as the mixing matrix). With U one can construct a new set of orbitals:   

φ′occ φ′vir



 =

∞ X k=0

k

U k!

!





 φocc   . φvir

(2)

Here, φocc and φvir indicate the column vectors of the sets of occupied {φi ; i = 1, · · · , nocc } and virtual {φa ; a = 1, · · · , nvir } ground state orbitals and φ′occ and φ′vir indicate the resulting sets of occupied and virtual excited state orbitals, respectively. As for other notations, the letters i, j are used to label the indices of occupied orbitals and a, b, c are used for virtual ones. The unitary transformation leads us to have a new set of “occupied” excited state

6

ACS Paragon Plus Environment

Page 7 of 50

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

orbitals

φ′i =

occ X

Yji φj +

vir X

Yai φa ,

(3)

a

j

and the corresponding density,





ρ (1, 1 ) = = +

occ X

i occ X

φ′i (1)φ′i (1′ ) ′

φi (1)φi (1 ) +

i occ occ XX i

vir X occ X a

∆Pai [φi (1)φa (1′ ) + φa (1)φi (1′ )]

i



∆Pij φi (1)φj (1 ) +

vir X vir X a

j

∆Pab φa (1)φb (1′ ) .

(4)

b

The change in the density matrix (∆P) can be written in terms of the matrix elements in the unitary transformation Y, ∆Paj =

occ X

Yai Yji ,

(5)

(Yji Yki − δjk ) ,

(6)

i

∆Pjk =

occ X i

∆Pab =

occ X

Yai Ybi .

(7)

i

The substitution of the excited state density into the Kohn-Sham energy expression yields the CV(n)-DFT excitation energy, ∆E. Note that the order, n, of this theory is defined by the highest power of the Un terms in the energy expression, but not in the exponential representation of Y given by Eq. (1) (see an example with two orbitals rotations in the Appendix). For instance, it is necessary to expand the orbitals up to second order in U in order to obtain the fourth order terms in the excitation energy in the CV(4)-DFT method 74 . The variational optimization of ∆E with respect to U results in a set of non-linear equations if n > 2, which are solved in SCF-CV(n)-DFT scheme 49 . Note that, in CV-DFT, we apply

7

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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 50

the constraint that one electron is fully transferred from occupied to virtual space vir X

occ X

∆Paa = 1 and

a

i

∆Pii = −1 .

(8)

Without this constraint, the energy can collapse to a lower state during the variational procedure. Until here all orbitals not directly participating in the transition are frozen (unchanged). We allow now these orbitals to react on the excitation as well, thus all orbitals are allowed to relax in applying a transformation described by the matrix R (see also Section 2.3). In the framework of the relaxed nth order constricted variational density functional theory, R-CV(n)-DFT (presented in a more general way in Ref. 50 ), the relaxed orbitals are given as

ψi (1) =φi (1) +

vir X c

ψa (1) =φa (1) −

occ X k

vir

Rci φc (1) −

occ

1 XX Rci Rck φk (1) , 2 c k vir

(9)

occ

1 XX Rak φk (1) − Rak Rck φc (1) . 2 c k

(10)

Note that the transposition of the relaxation matrix R (Ria = −Rai ) has been used in Eq. (9) and Eq. (10) and throughout the text in order to avoid confusion with the sign. Also we will refer to the submatrices U(occ × vir) and R(occ × vir) with naming them only U and R, while Uai and Rai stand for single matrix elements. The unitary transformation Y is now applied to the relaxed occupied {ψi ; i = 1 · · · nocc } and virtual orbitals {ψa ; a = 1 · · · nvir }, resulting in relaxed excited state orbitals {ψi′ ; i = 1 · · · nocc } and {ψa′ ; a = 1 · · · nvir }. 







′ ψocc 

ψocc  Y = . ′ ψvir ψvir

(11)

Now, the excited state density is described in terms of the relaxed orbitals and the change

8

ACS Paragon Plus Environment

Page 9 of 50

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

in the density matrix,

ρ′ (1, 1′ ) = +

occ X

ψi (1)ψi (1′ ) +

i occ occ XX i

vir X occ X a

∆Pai [ψi (1)ψa (1′ ) + ψa (1)ψi (1′ )]

i



∆Pij ψi (1)ψj (1 ) +

vir vir X X a

j

∆Pab ψa (1)ψb (1′ ) .

(12)

b

The transition matrix U itself can also be optimized in a self-consistent way, which leads to the self-consistent field nth order constricted variational density functional theory SCFCV(n)-DFT, or RSCF-CV(∞)-DFT if we allow orbital relaxation next to the optimization of the transition matrix. We will describe its framework as well as special cases, in which we have additional restrictions on the transition matrix U, in the next section, 2.2. An overview of the different variations of CV-DFT is given in Table 1. Table 1: Variation of CV-DFT applied. ‘T:’ indicates that it is introduced theoretically. ‘I:’ indicates that it is implemented into the code. ‘*’ indicates that it is introduced first in this work. Order n CV(2)-DFT CV(4)-DFT CV(∞)-DFT SCF-CV(∞)-DFT RSCF-CV(∞)-DFT SOR-R-CV(∞)-DFT COL-RSCF-CV(∞)-DFT SVD-RSCF-CV(∞)-DFT

2nd 4th ∞ ∞ ∞ ∞ ∞ ∞

Transition U Optimization Restrictions No No No Yes Yes No Yes Yes

No No No No No Uai = δab δij Uai = δij γ1 = 1

Relaxation R Order Optimization

Introduced Singlets Triplets

N/A N/A N/A N/A 2nd 2nd 2nd 2nd

T: 47 , I: 74 T: 76 , I: 74 T: 48 , I: 77 T: 48 , I: 49 T: 50 , I: 50 T: 51 T:∗ T:∗

No No No No Yes Yes Yes Yes

T: 49 , I: 75 T: 74 , I: 74 T: 49 , I: 49 T: 49 , I:∗ T: 50 , I:∗ T: 51 , I: 51 T:∗ , I:∗ T:∗ , I:∗

For the infinite order CV-DFT, we introduce the concept of the natural transition orbitals 71 (NTOs), they allow to describe the excitations in a more compact form than with canonical orbitals. We will also make use of this description in Section 2.3. Further, the NTOs will lead us to the constraint on CV-DFT of having a one electron excitation 48 . We will consider here the case of a triplet excitation. The transition matrix U can be decomposed to its singular values Uσα = Vσσ Σ (Wαα )T . 9

ACS Paragon Plus Environment

(13)

Journal of Chemical Theory and Computation

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 50

where Σii = γi and σ = α or β. This leads to the NTO φoi α =

occ X

(Wαα )ji φαj ,

(14)

(Vσσ )ai φσa ,

(15)

j

φvi σ =

vir X a

where j runs over the occupied ground state orbitals and a runs over the virtual ground state orbitals. With these NTO the “occupied” excited state orbital described by the unitary transformation Eq. (3) going to infinite order (n = ∞), 48 , can be rewritten by NTOs, thus φ′j = cos[ηγj ]φoj α + sin[ηγj ]φvj σ ;

j ∈ {occ} ,

(16)

where η is defined by the normalization condition that makes sure that only one electron is transferred from the occupied to the virtual space 48 occ X

sin(ηγi )2 = 1 ,

(17)

i

In the infinite order CV(∞)-DFT this restriction is identical to the restriction in Eq. (8) 48 . If there is only one non-zero singular value γi we have a single NTO transition, thus a transition from a single occupied NTO to a single virtual NTO, see also Ref. 48,54 . In a general case there will be more than one non-zero singular value and an excitation is dominated by single orbital replacement if ηγmax > 1.0 while the transition is best described by several orbital replacements if ηγmax < 1.0 54 .

2.2

RSCF-CV(∞)-DFT for Spin-Flip Transitions

A CV(n)-DFT theory in which the transition matrix U is optimized, the SCF-CV(n)-DFT method, has been previously introduced for the singlet excitations 50 . But for triplet excitations and thus spin-flip transitions, only orbital relaxation has been introduced 51 with a

10

ACS Paragon Plus Environment

Page 11 of 50

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

fixed transition matrix U. We will present here how the transition matrix U can be optimized in terms of resulting energy and derive the necessary gradient expression for spin-flip RSCF-CV(∞)-DFT. As a first we will look how the density is affected by such a transition: For the triplet excitation, the variation that mixes the occupied ground state orbitals of α-spin with the P βα β φa . Therewith the spin virtual ground state orbitals of β-spin is employed: δφαi ∝ a Uai

49 density change for the triplet state, ∆sU T , has to be introduced :

∆sU T =

q

U 2 U 2 2 (∆mU x,T ) + (∆my,T ) + (∆mz,T ) ,

(18)

where

∆mU ρU ρU x,T = ∆˜ T,αβ + ∆˜ T,βα occ(α) vir(β)

X X

=

k

βα β α [φa φk + φαk φβa ] , ∆Pak

∆mU y,T = 0 ,

∆mU ρU ρU z,T = ∆ˆ T,αα − ∆ˆ T,ββ occ(α) occ(α)

=

X X j

(19)

a

αα α α ∆Pjk φj φk

k

(20)

vir(β) vir(β)



X X a

ββ β β φa φb . ∆Pab

(21)

b

For further details we refer to Eq. (S1-3.18c), (S1-3.20b) and (S1-3.21.d) in Ref. 49 . We note here that in this work we follow the Wang-Ziegler definitions 78–80 , consistent with previous CV-DFT publications. A further discussion would go beyond the scope of this work and for a different approach we refer to Ref. 81 . The quantities in Eq. (18)-(21) have an additional superscript U to emphasize their dependence only on U. Now we shall introduce the relaxed ρU densities by substituting expansions (9) and (10) into ∆ˆ ρU T,αα and ∆ˆ T,ββ in Eq. (21) and

11

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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 12 of 50

keeping only the terms up to the second order in R: occ(α) occ(α)

∆ˆ ρU,R T,αα

X X

=

j

k

(2N,0) ∆ˆ ρT,αα

=

αα α α ∆Pjk ψj ψk (0+2N,1)

+ ∆˜ ρT,αα

(0+2N,2)

,

(22)

(0+2N,2)

.

(23)

+ ∆ˆ ρT,αα

vir(β) vir(β)

∆ˆ ρU,R T,ββ

X X

=

ββ β β ψa ψb ∆Pab

a

b (2N,0) ∆ˆ ρT,ββ

=

(0+2N,1)

+ ∆˜ ρT,ββ

+ ∆ˆ ρT,ββ

Here, the first index in parentheses of the superscript, (nU , nR ), means the order in U and (2N,0)

the second one means the order in R. For example ∆ˆ ρT,αα contains contribution from all even order U terms in ∆Pαα starting with 2 (N defines the set of all natural numbers) but (2N,0)

ρU no orbital relaxation contribution, thus ∆ˆ ρT,αα = ∆ˆ T,αα . The first (nR = 1) and second order (nR = 2) relaxed density matrix components in Eq. (22) and (23) are described as, (0+2N,1)

∆˜ ρT,σσ

(2N,1)

(0,1)

ρT,σσ (1, 1′ ) (1, 1′ ) = ∆˜ ρT,σσ (1, 1′ ) + ∆˜ occ(σ) vir(σ)

X X

=

i

(1)σσ

Tai

[φσi (1)φσa (1′ ) + φσa (1)φσi (1′ )] ,

(24)

a

and, (0+2N,2)

∆ˆ ρT,σσ

(2N,2)

(0,2)

ρT,σσ (1, 1′ ) (1, 1′ ) = ∆ˆ ρT,σσ (1, 1′ ) + ∆ˆ occ(σ) occ(σ)

vir(σ) vir(σ)

=

X X a

(2)σσ Tab φσa (1)φσb (1′ )

+

X X i

b

(2)σσ σ φi (1)φσj (1′ ) ,

Tij

(25)

j

where T(1)σσ and T(2)σσ are the changes in density matrix ∆P mixed with first and second order relaxation matrix Rσσ , respectively. Thus, the T(1)σσ has 0+2N order in U and first order Rσσ , whereas the relaxation matrix Rσσ is second order in T(2)σσ . The details of the derivations are given in the supporting material. It also should be emphasized that 12

ACS Paragon Plus Environment

Page 13 of 50

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

the relaxation matrices Rσσ only affect orbitals of the same spin (σ = α or β), there is no ‘crossed-relaxation’ for ∆mU x,T depending on the odd terms in U. We have now defined all components of the excited state density. The substitution of Eq. (19), (22), and (23) into the Kohn-Sham energy expression and the subtraction of the ground state energy yields the following excitation energy expression ∆ETU,R



 i h 1 U,R β 1 U,R = + ∆ρT , ρ0 + ∆ρT − ET ρα0 , ρβ0 2 2   Z 1 U,R β 1 U,R α dν1 ∆ρU,R = FKS ρ0 + ∆ρT , ρ0 + ∆ρT T 2 2 Z 1 1 β U α α β U α = {FKS [ρα0 + ∆mU x,T , ρ0 − ∆mx,T ] − FKS [ρ0 , ρ0 ]}∆mx,T dν1 4 4 Z 1 U,R β 1 U,R U,R α + FKS [ρα0 + ∆ˆ , ρ + ∆ˆ dν1 ρ ρ ]∆ˆ ρT,αα 2 T,αα 0 2 T,ββ Z 1 U,R β 1 U,R U,R β dν1 . ρT,ββ ρT,αα , ρ0 + ∆ˆ ρT,ββ ]∆ˆ [ρα0 + ∆ˆ + FKS 2 2 ETU,R

ρα0

(26)

β α are the α and β components of Here, ρα0 and ρβ0 are the ground state densities, FKS and FKS

the Kohn-Sham Fock operator. The total excitation energy ∆ETU,R can be decomposed into a part with the unrelaxed triplet excitation energy ∆ETU (see, for example, Eq. (S1-3.44c) and (S1-3.44d) in 49 for the details of the derivation) and the relaxation energy ∆ETR . ∆ETU,R = ∆ETU + ∆ETR

(27)

Let us have a look at these two contributions separately. The unrelaxed energy becomes

∆ETU

1 1 β α U α α β U {FKS [ρα0 + ∆mU x,T , ρ0 − ∆mx,T ] − FKS [ρ0 , ρ0 ]}∆mx,T dν1 4 4 Z 1 U 1 U α [ρα0 + ∆ˆ + FKS ]∆ˆ ρU ρT,αα , ρβ0 + ∆ˆ ρ T,αα dν1 2 2 T,ββ Z 1 U 1 U β [ρα0 + ∆ˆ + FKS ρT,αα , ρβ0 + ∆ˆ ρT,ββ ]∆ˆ ρU T,ββ dν1 . 2 2 =

Z

(28)

The relaxation energy is only considered up to second order in the relaxation matrix R. The

13

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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 14 of 50

terms with higher order terms are also derived from the Eq. (26), but their contribution to the excitation energy is negligible (the details of the derivation of higher order terms are given in the supporting material). It is convenient for the numerical analysis to decompose the relaxation energy further into different contributions R(2′ )

∆ETR = ∆ET

R(2′′ )

+ ∆ET

,

(29)

where we have the contributions of first order in R R(2′ ) ∆ET

1 (0+2N,1) β α FKS [ρα0 + ∆ˆ ρU , ρ0 + ∆ˆ ρU ρ T,αα + ∆˜ T,ββ + 2 T,αα Z 1 (0+2N,1) β β [ρα0 + ∆ˆ ρU + FKS ρU ρT,αα , ρ0 + ∆ˆ T,αα + ∆˜ T,ββ + 2

=

Z

1 (0+2N,1) (0+2N,1) ∆˜ ρT,ββ ]∆˜ ρT,αα dν1 2 1 (0+2N,1) (0+2N,1) ρT,ββ dν1 ∆˜ ρT,ββ ]∆˜ 2 (30)

and the second order relaxation contributions, thus second order in R R(2′′ ) ∆ET

=

Z

(0+2N,2) α FKS [ρα0 , ρβ0 ]∆ˆ ρT,αα dν1

+

Z

(0+2N,2)

β [ρα0 , ρβ0 ]∆ˆ ρT,ββ FKS

dν1 .

(31)

Note that U and R are coupled due to the relaxation. The derivation of the energy gradients g follows the procedure of Ref. 50 and the expressions are given in the supplementary material.

2.3

∆SCF-DFT within the RSCF-CV(∞)-DFT Scheme

Previously, it has been shown that the linear response TD-DFT and ∆SCF-DFT can be obtained as special cases of CV-DFT 47,51,82 . Such a connection between two apparently different theories can be established, because CV-DFT has the analytical excitation energy expression for the difference between the ground and excited state. A distinction between two theories within the CV-DFT approach can be briefly sketched as follows. To obtain an equivalence of TD-DFT we restrict the expansion in Eq. (1) to second order (n = 2). 14

ACS Paragon Plus Environment

Page 15 of 50

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

Contrary to this, to mimic ∆SCF-DFT we include all terms in Eq. (1) up to infinity, but we impose additional restrictions on the form of the transition matrix U. Further, the equivalence between adiabatic TD-DFT and CV(2)-DFT has been obtained initially for the Tamm-Dancoff approximation 47 , and only recently it has been shown that these two theories coincide without the Tamm-Dancoff approximation, if the frequency dependent variational principle is employed 82 . The situation with ∆SCF-DFT versus CV-DFT is somewhat similar. Namely, it has been shown earlier that CV-DFT and ∆SCF-DFT are equivalent for transitions of single orbital replacement (SOR) type, which was confirmed numerically for non-hybrid functionals 51 . In this section we shall show that the analogous CV-DFT excitation energy expression can be derived for a much looser restriction of a single NTO transition and how one can obtain the corresponding transition matrix U. Before we proceed with the restrictions on the transition matrix U, we stress that for the second order theory it was enough to demonstrate the equivalence of CV-DFT and adiabatic TD-DFT eigenvalue equations, but a similar proof for ∆SCF-DFT is complicated by additional orbital relaxation. As a consequence, there is no analytical formula for the total ∆SCF-DFT excitation energy with the relaxation. Therefore, we shall confirm our derivations by the results of numerical simulations within the RSCF-CV(∞)-DFT scheme in Section 4.1. We shall now make use of the fact that the ∆SCF-DFT excitation is represented by a single KS-determinant. Therefore, in order to establish a connection to CV-DFT the transition matrix U has to have such a form that excludes the configuration interaction between different states. Let us now analyze the expression of the SCF-CV(n)-DFT triplet excitation energy for the general (or unrestricted) transition matrix U in order to identify such terms. Within the Tamm-Dancoff approximation and by making use of the NTOs, Eq. (28) for the infinite order theory can be written as,

∆ETU = sin2 [ηγi ]

( occ X i

εivβ (ρ0 + 1/2∆ρ′ ) − εioα (ρ0 + 1/2∆ρ′ ) 15

ACS Paragon Plus Environment

)

Journal of Chemical Theory and Computation

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

+

occ occ X X i

Page 16 of 50

sin[ηγi ] cos[ηγi ] sin[ηγj ] cos[ηγj ]Kioα ivβ j oα j vβ

(32)

j

The two electron integral is decomposed into Coulomb and exchange-correlation, C XC Kpq,st = Kpq,st + Kpq,st ,

(33)

with C Kpq,st

=

Z Z

ψp (1)ψq (1)

1 ψs (2)ψt (2)dν1 dν2 r12

(34)

The exchange-correlation integral is given separately for the local (KS) and non-local (HF), XC(KS) Kpq,st

=

Z

ψp (r1 )ψq (r1 )f (r1 )ψs (r1 )ψt (r1 )dr1

(35)

and XC(HF ) Kpq,st

=−

Z Z

ψp (1)ψq (1)

1 ψs (2)ψt (2)dν1 dν2 r12

(36)

where f (r1 ) represents the regular energy kernel. We note that the first term in Eq. (32) includes the orbital energy differences calculated σ as the diagonal matrix elements of the Kohn-Sham Fock operator FKS at the “transition

density”, which includes the density change due to the excitation. The second term of Eq. v

(32) includes the K integrals that describe the interaction between the individual φoi α → φi β v

and φoj α → φj β NTO transitions. It is easy to see that all coefficients in front of the K integrals would be zero, if ηγj = π/2 for j = i whereas γj = 0 for j 6= i. Under these v

conditions the excitation is reduced to the single NTO replacement φoi α → φi β , while the general expression (32) simplifies to ∆ETU = εivβ (ρ0 + 1/2∆ρ′ ) − εioα (ρ0 + 1/2∆ρ′ )

(37)

We note in passing that Eq. (37) the ε are here like in Eq. (32) clearly for NTOs (and thus

16

ACS Paragon Plus Environment

Page 17 of 50

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

a possible linear combination of molecular ground state orbitals). This excitation energy of Eq. (37) corresponds to the following ∆SCF-DFT-type KS-determinant 48 : v

α α α Ψ′I = |φo1α φo2α ...φoi−1 φi β φoi+1 ...φoN/2 φβ1 ...φβN/2 |

(38)

When the orbitals in Ψ′I are optimized without any orthogonality constraint to the lower excited states, like in the ∆SCF-DFT method, such an optimization would lead to the variational collapse to a lower excited states. In CV-DFT the variational collapse to the ground state is prevented due to the constraints of Eq. (8), Eq. (16) respectively and to lower excited states due to the orthogonality of the transition matrices U. It is trivial to v

maintain the orthogonality when ψioα consists of one occupied, and ψi β consists of one virtual canonical orbital (thus the transition matrix U has only a single entry) as in the previously implemented SOR scheme 51 . We shall now describe a procedure for preserving orthogonality in CV-DFT between the excitations stemming from a more general single NTO to a single NTO transition. In order to restrict the general SCF-CV(n)-DFT excited state KS-determinant to the ∆SCF-DFT type given by Eq. (38), the transition has to be enforced to a single NTO to a single NTO transition, which corresponds to a transition where the transition matrix U has a single singular value being non-zero. This can be obtained in enforcing the rank one approximation, thus retaining only the largest singular value in Eq. (13). It results in the following structure of the transition matrix U: U = v1 w1T

(39)

with v1 and w1 for the vectors of the largest singular value of Vββ and Wαα out of Eq. (13). But such a transition matrix U does in general not fulfill the orthogonality constraint and neither does the transition matrix Ui = Ui−1 + ∆U, thus the transition matrix during an SCF optimization in cycle i where ∆U is calculated according to the general formulas Eq. 17

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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 18 of 50

(SI58 ∼ SI60) in the supplementary material. Thus, it is necessary to maintain the orthogonality between the transition matrices U and restrict them to the rank one approximation simultaneously. Orthogonality to lower excitations in form of trace orthogonality 50 simplifies for transitions matrices UI and UII with the rank one approximation as following: XX a

i

IT II Uia Uai =

XX i

vaI wiI vaII wiII = v1I , v1II

a



 w1I , w1II = 0

(40)

  As we can see, Eq. (40) is equivalent to v1I , v1II = 0 or w1I , w1II = 0.

We shall now describe how we obtain such a transition matrix, starting from a transition

matrix UI (independent if UI = U0 or UI = Ui + ∆U during the SCF optimization) of excitation I. We shall use the symbol ‘tilde’ to mark orthogonalized matrices and vectors. For the excitation I we orthogonalize the vector xI1 (x being either v or w) to the corresponding ˜ II vector x 1 for each lower excitation II (starting with the lowest one), depending on which projection leads to the smaller spectral norm of the difference in the transition matrix of ′ ˜ I ||2 ). We note hereby that x ˜ II excitation I (||U I − U 1 is orthogonalized to possible correspond-

ing vectors of lower lying excitations III if (but only if) for excitation I the excitation II is ˜ 1I and projected out in the same space (being V or W) as excitation III. With these vectors v ˜ 1I we construct the transition matrix being trace orthogonal to all lower lying transitions w and of rank one, ˜I = v ˜ 1IT ˜ 1I w U

(41)

˜ I }, the KS-determinants It can be shown that with this orthonormal dual basis set {˜ vI , w of type (38) are orthogonal to each other and to the ground state. So far we have considered only the optimization of the restricted transition matrix U without any orbital relaxation. In the general RSCF-CV(∞)-DFT method the relaxation is introduced according to the formulas given in Eq. (9) and (10). Note that in a procedure described in this section the transition matrix U can include all possible canonical occupied to canonical virtual ground state orbital transitions. As a result, the same canonical virtual 18

ACS Paragon Plus Environment

Page 19 of 50

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

orbitals might contribute both to the excitation and relaxation. In order to see how the relaxation of orbitals that participate in the transition is avoided in RSCF-CV(∞)-DFT, let us consider first the relaxation in a single canonical orbital replacement ψk → ψc , in which the matrix elements of the transition matrix are given by Uai = δca δki . The relaxed excited state in this case is given by |ΨR k→c i =

vir occ X X a

i

˜ ai |Φca R ki i .

(42)

Here, |Φca ki i is the doubly excited Slater determinant that includes ψi → ψa and ψk → ψc transitions, while the matrix elements of the relaxation vector exclude Rck component, i.e. ˜ ai = Rai (1 − δca )(1 − δki ). The state in Eq. (42) is orthogonal to the ground state, because R it consists of the orthonormal set of occupied orbitals, in which two of them were replaced by the virtual ones. ˜ in the SVD method is orthogonalized to the On the other hand, the relaxation vector R corresponding transition vector: ˜ = R − (UTo v R)Uio →iv . R i →i

(43)

˜ to the NTO basis. We can do it by constructing a set of It is convenient to transform R ψj o → ψj v transitions that are orthogonal to ψio → ψiv , i.e. UTjo →j v Uio →iv = 0. Based on these auxiliary NTO transitions we can construct the following projections: ˜ , Rj v j o = (UTjo →j v R)

(44)

˜ is given by Eq. (43). As a result, the relaxed excited state in the SVD method will where R have the form: |ΨR io →iv i

=

occ X

j o 6=io

19

v v

Rj v j o |Φiio jjo i .

ACS Paragon Plus Environment

(45)

Journal of Chemical Theory and Computation

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 50

v v

The formula (45) is a generalization of Eq. (42), in which |Φiio jjo i consists of the orthonormal set of the occupied orbitals {ψio ; i = 1, ..., nocc }, in which two of them were replaced by the virtual ones. Now it is easy to see that the relaxation and excitation are separated in the NTO basis. In other words, the relaxation from the ‘active’ (participating in a transition) and ‘passive’ orbitals maintains the orthogonality to the ground state. Note that we don’t distinguish in Eqs. (42) and (45) between the spin flip and spin conserving transitions, so Eqs. (42) and (45) are fulfilled in both cases. Note that currently a simpler approach is implemented in the CV-DFT code given by Eq. (43), but the transformation of the relaxation matrix to the NTOs basis as in Eq. (45) would give us a tool to estimate the weight of a particular doubly excited configuration present in the relaxed state. In this regard, we would like to mention here one important aspect of the orbital relaxation in RSCF-CV(∞)-DFT, namely, mixing in the double excitations to the single ones 83 . This interaction is required, for example, for an adequate description of photochemical processes 24 , and CV-DFT allows us to demonstrate that this feature is naturally incorporated in the ∆SCF-DFT method. As it was shown before 83 , it is enough to consider the relaxation energy up to first order in Rσσ . It is given for the considered in this section single state excitations by the following formula: R(1′ ) ∆ET

=

occ vir X X a

i

Rai Kai,ivβ ivβ − Kai,ioα ioα



(46)

If we now expand the matrix elements expressed in terms of the NTOs into the canonical orbitals, we obtain the molecular integrals with four different indices that connect doubles with singles. Note that similar integrals but with the exact two-electron operator r12 instead of the exchange-correlation kernel, fXC , are introduced ad hoc in the dressed TD-DFT 25 . However, the explicit incorporation of the double excitations into the eigenvalue equations represent a challenge for the efficient implementation. To overcome this problem only certain

20

ACS Paragon Plus Environment

Page 21 of 50

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

double excitations are included. Contrary to this, RSCF-CV(∞)-DFT as well as ∆SCF-DFT v

includes implicitly all possible double excitations through the response to the φoi α → φi β transition. We should stress here that double excitations are for many chemists understood as an excitation of two electrons, while in CV(∞)-DFT excitations are restricted to be of one electron excitations. However, the total density change (excitation plus relaxation, either occ-occ or vir-vir) in RSCF-CV(∞)-DFT is not normalized to one anymore, because the excited density change is normalized to one, but the relaxation density is not 83 . Therefore, we have a mixing of the doubly excited configurations. Another aspect of the doubly excited configurations is the interaction with the ground state. Although the overlap between |Ψ0 i and |Φca ki i is zero, these states can mix indirectly, ˆ ca i = ˆ is the standard molecular Hamiltonian. This is true for i.e. hΨ0 |H|Φ 6 0, where H ki the configuration interaction method which includes at least doubles and singles (SDCI). Although the RSCF-CV(∞)-DFT is not equivalent to the SDCI, we can make a connection to the wave function theory by analyzing the RSCF-CV(∞)-DFT scheme within an exchangeonly Hartree-Fock formalism. For this purpose we have to analyze the matrix elements that might connect the doubles to the ground state in a similar way as we have done in the analysis of Eq. (46). First of all, we notice that the matrix elements either in Eq. (42) or (45) include the first order change in R as well as in U. Strictly speaking, it is the change in the density due to odd order terms with regard to the infinite order CV-DFT. Anyways, the relaxed excitation energy expression of RSCF-CV(∞)-DFT does not include such terms (see Eq. (26) and Eq. S38 of Supporting Information): ˆ ca Uck Rai hΨ0 |H|Φ ki i = Uck Rai Kck,ai

(47)

Such terms could be included in principle, but they were omitted in the formalism, be-

21

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

cause the relaxation in RSCF-CV(∞)-DFT was introduced in a way similar to ∆SCF-DFT. Namely, the relaxation energy terms include cross terms between the second order density change in U (or even order contribution, i.e. occ-occ and vir-vir blocks of the excited state density matrix) and the first and second order relaxaion density change. At the moment when the RSCF-CV(∞)-DFT was formulated 50 , the role of the omitted terms given by Eq. (47) was not clear. However, now we can see that these terms connect the doubly excited configurations with the ground state (at least for the spin conserving transitions). This feature is important for the calculation of the conical intersections 84 . Currently, only a few DFT-based methods are considered as appropriate tools for such studies, where we would like to mention: ensemble DFT 84,85 , spin-flip TDDFT (with the high-spin reference state) 84,86 , Brillouin corrected dressed DFT 25,87 , and the configuration interaction-corrected TDA method 88 . In the last two methods the interaction between the ground and excited states is achieved by the Brilloiun correction which can be considered as ad hoc approach unlike ensemble DFT. Thus, the inclusion of the above mentioned cross terms will enable the conical intersections modeling within the RSCF-CV(∞)-DFT framework. However, this will be the topic of our future studies. So far we have outlined the key points of the modified RSCF-CV(∞)-DFT method which coincides with ∆SCF-DFT in a sense that it is restricted to single state excitations. We stress here once again that numerical simulations are necessary to corroborate this theoretical finding. In order to investigate the importance of the transition matrix U optimization, the calculations with three different ∆SCF-DFT-like schemes constructed by gradually loosening the restriction on the U matrix have been performed. A graphical illustration of these different restricted transition matrices is given in Figure 1 (the illustration shows submatrices of the vertical triplet 1 A’ excitation in propanamide discussed in Section 4.2). The first and the most simple one is the single orbital replacement, which has been implemented before in our previous publication 51 . Since it has already been described in

22

ACS Paragon Plus Environment

Page 22 of 50

Page 23 of 50

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

1 0.5 0 !0.5 !1

Figure 1: Schematic illustration of the three different restricted transition matrices U: on SOR COL the left the SOR, where Uai = δab δij . In the middle COL, where Uai = δij and therefore SVD an entire column of U can have values 6= 0. On the right SVD, where all Uai can have values 6= 0 as long γ1 = 1 and γi = 0 ∀i ≥ 2. The matrices shown are submatrices of the vertical triplet 1 A’ excitation in propanamide discussed in Section 4.2.

details, here we just mention that in this approximation the transition is from one canonical occupied to one canonical virtual ground state orbital φαk → φβc , which are frozen, and, thus, only one element Uck is non zero. We will refer to the RSCF-CV(∞)-DFT with this special restriction as SOR-R-CV(∞)-DFT or simply SOR within this manuscript. The second approximation involves the excitation of an electron from a frozen occupied v

orbital to the linear combination of virtual orbitals, φαk → φkβ . During optimization all v

virtual ground state orbitals can be included in φkβ . We shall term this approximation COL-RSCF-CV(∞)-DFT and refer to it within this manuscript as COL, since the transition matrix U has only one column of non zero matrix elements (see Figure 1). Of course, there is another generalization of the SOR transition in which the occupied orbital, φoi α , is optimized, while the virtual one is frozen. However, we think that the COL approximation is more practical, since this type of a transition might take place, for example, in the excitations from the core orbitals 45 . The third single state approximation is the most general, and it involves the optimization of both the occupied and virtual orbitals in the excited state KS-determinant (38) according to the description in this section, so the transition matrix U has the form given by Eq.

23

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

(39). We shall term this approximation SVD-RSCF-CV(∞)-DFT and refer to it within this manuscript as SVD. Note that all three approximations involve full orbital relaxation as described in Section 2.2. We also note here that these variants of RSCF-CV(∞)-DFT (SOR, COL, SVD) are in general restricted to abelian groups, then for non-abelian groups some excitations result from the direct product of two or more transitions (e.g. in Td an A1 resulting from t1x × t1x + t1y × t1y + t1z × t1z ), nicely visible for SOR where these excitations need more than a single matrix element. In the original ∆SCF-DFT method such intractable cases can often be resolved by evaluation of certain two-electron integrals 89 . However, this is not the case for the general RSCF-CV(∞)-DFT method 83 . First of all, because the RSCF-CV(∞)-DFT can describe excitations that consists of a linear combination of several individual NTO transitions. Second of all, it contains additional terms in the energy expression that are absent in the ∆SCF-DFT-like methods, but present in the adiabatic TDDFT (i.e. configuration interaction).

3

Computational Details

All calculations were preformed using an implementation in a developer’s version of the ADF 2012 program 90–92 . The standard triple-ζ Slater type orbital basis with one set of polarization function (TZP) 93 was used and all electrons were treated without the use of frozen core approximation 90 . For the reliability of the calculations including exact exchange, more diffuse functions 92 were added in all calculations. The functionals adapted in this calculation are the local density approximation in the VWN parametrization (LDA) 94 , the generalized gradient approximation (GGA) of Becke’s 1988 exchange with Perdew’s 1986 correlation (BP86) 95,96 , and the two hybrid functionals including Becke three-parameter hybrid functional combined with Lee, Yang, and Parr (LYP) correlation (B3LYP) 95,97,98 and Becke’s half and half HF/DFT hybrid exchange with the LYP correlation (BHLYP) 97,99 .

24

ACS Paragon Plus Environment

Page 24 of 50

Page 25 of 50

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

These two hybrid functionals, B3LYP and BHLYP, have different portions of Hartree-Fock exchange (20% and 50%, respectively). In order to compare with the results reported in the publication of Silva-Junior et al. 100 , we used the same geometries which are reported by Schreiber et al. 101 .

4

Results and Discussion

We have the general method RSCF-CV(∞)-DFT and its restricted versions mimicking ∆SCF-DFT (SOR, COL SVD) to evaluate and compare. For this we calculate the different vertical excitation energies obtained from the triplet states of the test set in Ref. 100 (only benzene is excluded as it belongs to a non-abelian point-group, see section 2.3) using different functionals. The data is presented and discussed as following: As a first issue we compare in Section 4.1 how well the different ∆SCF-DFT-like versions of RSCF-CV(∞)-DFT (SOR, COL, SVD) reproduce the excitation energy obtained by ∆SCF-DFT. Then we compare in Section 4.2 the different ∆SCF-DFT-like versions of RSCF-CV(∞)-DFT with each other and thereby have a look at the different restrictions and contributions. As a last part, Section 4.3, we compare these ∆SCF-DFT-like energies and the general RSCF-CV(∞)-DFT energies with the values obtained by Silva-Junior et al. 100 with highly accurate ab initio methods.

4.1

Comparison of ∆SCF-DFT-like RSCF-CV(∞)-DFT Methods with ∆SCF-DFT

For SOR and COL we restrict the nature of the single natural type orbital transition and for all three methods (SOR, COL, SVD) orbital relaxation is included only up to second order, therefore we expect for all methods to obtain energies higher than the corresponding ∆SCFDFT transition energies. This is the case for the large majority of all obtained energies and we found a maximum negative deviation of 0.04 eV (obtained for pyridine 1 A2 with BHLYP). 25

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

A comparison between the three restricted RSCF-CV(∞)-DFT methods and ∆SCF-DFT is given in Table 2. We note here that the ∆SCF-DFT calculations for imidazole 3 A’ when using LDA and BP86 as well as 2 A’ and 3 A’ when using B3LYP and BHLYP did suffer from convergence problems, therefore these excitations are not included in the comparison. Let us have now a look at these values regarding the influence of the functional. Using the local density approximation (LDA), the root mean square deviation (RMSD) of the excitation energies obtained by a ∆SCF-DFT-like method (SOR, COL, SVD) compared to the ∆SCF-DFT excitation energies is for each method smaller than 0.10 eV. Very similar excitation energy differences are obtained when using the generalized gradient approximation BP86 (see Table SV in the supplementary material). Adding Hartree-Fock exchange by using the hybrid functionals B3LYP and BHLYP, we observe for SOR an increase of the energy differences, which is especially large with BHLYP, where more Hartree-Fock exchange is included. For the other two methods, COL and SVD, where the transition matrix U is less restricted, the ∆SCF-DFT excitation energies are reproduced with a similar quality when using B3LYP as when using one of the non-hybrid functionals. This changes when using BHLYP, but the energy differences still stay significantly lower than the ones obtained by SOR (RMSD=0.09 eV for COL and 0.06 eV for SVD compared to 0.44 eV for SOR). Thus, an increase of the exact Hartree-Fock exchange in the functional requires more flexibility in the transition matrix U in order to reproduce ∆SCF-DFT results, while the SOR approximation is in general quite effective for non-hybrids, as already found in Ref. 51 .

4.2

Comparison of the Different Restrictions and Contributions

Let us have a look at the differences in the total RMSD between the different ∆SCF-DFTlike methods. As visible in Table 2, the RMSD and mean deviations confirm throughout the trend to mimic ∆SCF-DFT more successfully when imposing less restrictions on the transition matrix U. This is also shown in Figure 2, where one can see that the ∆SCF-DFT excitation energies are with only one exception (2 B1u state of s-tetrazine) better reproduced 26

ACS Paragon Plus Environment

Page 26 of 50

Page 27 of 50

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

Table 2: Comparison of excitation energies[a] obtained by ∆SCF-DFT-like RSCF-CV(∞)DFT methods (SOR, COL, SVD) with respect to those obtained with ∆SCF-DFT.[b] [a]: Vertical excitation excitations of the test set in Ref. 100 . The values for each excitation are given in the supplementary material. [b]: Values are given in eV. Functional LDA

B3LYP

BHLYP

RMSD Mean dev. Max. (+) dev. RMSD Mean dev. Max. (+) dev. RMSD Mean dev. Max. (+) dev.

SOR

COL

SVD

0.09 0.05 0.29 0.13 0.08 0.45 0.44 0.25 1.55

0.06 0.04 0.20 0.05 0.03 0.18 0.09 0.05 0.35

0.02 0.01 0.08 0.02 0.01 0.08 0.06 0.02 0.32

using the COL approximation than with SOR and in general even better reproduced using the SVD method. A remarkable exception is the imidazole 3 A’ excitation, where due to orthogonality to lower lying states SVD results in a clearly higher excitation energy than COL (0.42 eV with BHLYP, see Table SIV of the supplementary material) The 2 B1u state of s-tetrazine is indicative of another special case, for which the SOR type of the transition matrix U accidentally is the optimal one for all three methods, resulting in the same excitation energy. It is now logical to look at the differences in the transition matrix U that cause such changes in the excitation energy. For the comparison of the transition matrices U of the three ∆SCF-DFT-like methods we choose as an example the 1 A’ transition in propanamide using BHLYP, where we find the largest excitation energy deviation (the excitation energies obtained by SOR and ∆SCF-DFT differ by 1.55 eV).  ” ” SOR 5 a ,4 a = The transition matrix U has for SOR only one element unequal zero; Uai  SOR 5 a” ,4 a” = 0. For COL the excitation is still restricted to stem from only one 1, Uai

occupied orbital, here 4 a” , but we allow to have all virtual orbitals contributing to the excitation and thus instead of the single virtual orbital entry 5 a” we find in the transition

27

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

Figure 2: Excitation energy differences from ∆SCF-DFT using the B3LYP functional for SOR (black squares), COL (red circles) and SVD (blue triangles). matrix U a linear combination of several virtual orbitals contributing to the excitation: 



 0.837      UCOL [β , α ] = −0.515 vir occ ai     0.150

αocc = {ψ4 a” } βvir

= {ψ5 a” , ψ6 a” , ψ7 a” }

where only the matrix elements with an absolute weight larger than 0.100 are given. As we can see, the weight of the major transition 4 a” → 5 a” is reduced from 1.000 to 0.837, while two other transitions (4 a” → 6 a” and 4 a” → 7 a” ) contribute significantly. As a result of this weaker restriction the excitation energy of COL drops from 6.99 eV (in SOR) to 5.61 eV. The transition matrix U of the same excitation but for the SVD method incorporates

28

ACS Paragon Plus Environment

Page 28 of 50

Page 29 of 50

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

more “degrees of freedom” and we obtain 



0.141 0.803   0.153     [β , α ] = USVD −0.096 −0.088 −0.503 vir occ ai     0.027 0.025 0.141

αocc = {ψ2 a” , ψ3 a” , ψ4 a” } βvir

= {ψ5 a” , ψ6 a” , ψ7 a” }

where we only show matrix row and column combination leading to an element with absolute weight larger than 0.100. Although, three contributions (U5 a” ,4 a” , U6 a” ,4 a” and U7 a” ,4 a” ) stemming from an excitation out of the occupied orbital 4 a” are very similar to the corresponding ones obtained by COL, but there are also transitions from 2 a” and 3 a” occupied orbitals to the same virtual ones (5 a” , 6 a” and 7 a” ). The resulting excitation energy is further stabilized in comparison to COL method by 0.21 eV, and it differs only by 0.04 eV from the corresponding energy obtained with ∆SCF-DFT. It is easy to see that the submatrix of the transition matrix U given above, USVD ai [βvir , αocc ], can be represented in the form of Eq. (39), with v1ββ [βvir ] = NvC (0.803, −0.503, 0.141) and w1αα [αocc ] = NwC (0.190, 0.176, 1.00), where NvC and NwC are the normalization constants. As a further aspect we look at the importance of the orbital relaxation. Indeed, it was shown before 51 that the relaxation energy contribution can be very substantial for the SOR method. Is it so important for COL and SVD, where the U matrix is allowed to vary? In order to answer this question we have plotted different contributions to the triplet excitation energy for the SVD method obtained with the B3LYP functional in Figure 3. These different contributions are the unrelaxed energy according to Eq. (37) and the two parts of the contribution due to the orbital relaxation according to Eq. (30) and Eq. (31). These contributions of the orbital relaxation are not necessarily large but can be large. While we obtain for the 1 B2 state of norbornadiene only a total relaxation energy contribution of −0.11 eV this contribution sums up to −1.71 eV for the 1 A2 state of pyridine. In the second case the stabilizing part of the relaxation energy given by Eq. (30) contributes even with as much as −2.45 eV. By contrast the contribution given by Eq. (31) is positive and numerically 29

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

Figure 3: Contributions to the total excitation energy of the triplet states obtained with SVD with B3LYP: Excitation energies without orbital relaxation (black squares, Eq. (37)) and the two different second order orbital relaxation contributions (dark blue circles for the contribution stemming from Eq. (30) and light blue triangles for the part out of (31)).

αα

Further shown are the norms of the relaxation matrices R (maroon crosses) and Rββ (orange plus signs). For the visibility, the values of norms are multiplied by ten.

30

ACS Paragon Plus Environment

Page 30 of 50

Page 31 of 50

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

smaller (with a maximum of 0.73 eV for the pyridine 1 A2 transition). We obtained total relaxation energies contributing more than 0.80 eV only for n → π ∗ excitations, but out of these all except the ones from s-tetrazine and p-benzoquinone show this large total relaxation energy. This is because the change of density stemming from the n → π ∗ transition can be substantial in some region of space, whereas in s-tetrazine and p-benzoquinone the n orbitals are spread over several atomic centers in molecules resulting in minor relaxations 51 . Therewith we get also for optimized transition matrices U a similar picture to the one for the SOR method results based on non-hybrid functionals published previously 51 . But despite the possible large energy contributions, the norm of the relaxation matrix R is in general small:

over this test set with the B3LYP functional the mean and maximum values for Rαα are



0.08 and 0.16, while for Rββ 0.11 and 0.21, see Figure 3. We note here that first Rββ

usually takes on larger values but does not always do (over the test set using the B3LYP



functional we obtain Rββ < Rαα in 9 cases). Second, while there is a trend between

the sum of the relaxation matrix norms and the relaxation energy, they are not necessarily

correlated which we can be seen from the two following cases: For the pyridine 1A2 state

we obtain a relaxation energy of −1.71 eV with relaxation matrix norms Rαα = 0.15 and

ββ

R = 0.19, while for the p-benzoquinone 1B1u state we have relaxation matrix norms of

αα

R = 0.16 and Rββ = 0.19 resulting in a relaxation energy of −0.27 eV (in both cases values are obtained using B3LYP as functional).

4.3

Comparison with Benchmark Values

To this point we compared the excitation energies of the different ∆SCF-DFT-like versions of RSCF-CV(∞)-DFT versus ∆SCF-DFT to see, how good we can mimic ∆SCF-DFT. But we did neither compare them with the more general RSCF-CV(∞)-DFT nor did we have a look what these excitation energies worth, why we will now compare the results from this work with the “best” results for the same benchmark test set obtained with highly accurate ab initio methods obtained by Silva-Junior et al. 100 . An overview of this comparison is given 31

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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 32 of 50

in Table 3 for the general RSCF-CV(∞)-DFT, its restricted version SVD (which is, as we showed in Section 4.1 numerically equivalent to ∆SCF-DFT) and TD-DFT as an indication. Table 3: Comparison of excitation energies[a] obtained by TD-DFT, RSCF-CV(∞)-DFT and SVD-RSCF-CV(∞)-DFT (SVD) with respect to the “best” estimate values out of Ref. 100 .[b] [a]:Vertical excitation excitations of the test set in Ref. 100 . The values for each excitation are given in the supplementary material. [b]: Values are given in eV and only consider excitations with successful ∆SCF-DFT-calculations. Functional

LDA

B3LYP

BHLYP

RMSD Mean dev. Max. (+) dev. Max. (-) dev. RMSD Mean dev. Max. (+) dev. Max. (-) dev. RMSD Mean dev. Max. (+) dev. Max. (-) dev.

TD-DFT

RSCF-CV(∞)-DFT

SVD

0.56 -0.32 0.51 -1.51 0.32 -0.22 0.33 -0.95 0.31 -0.14 1.03 -0.61

0.51 -0.23 0.58 -1.76 0.37 -0.24 0.33 -1.57 0.42 -0.16 0.82 -1.10

0.48 -0.13 0.73 -1.64 0.39 -0.15 0.76 -1.45 0.47 -0.05 1.08 -1.09

The general RSCF-CV(∞)-DFT method as well as the SVD method come closest to the “best” values when using the B3LYP functional, where the general RSCF-CV(∞)-DFT yields an RMSD of 0.37 eV, which is only 0.05 eV larger than the corresponding TD-DFT value and only marginally better than SVD (0.39 eV). The general RSCF-CV(∞)-DFT method and its restricted version SVD yield very similar RMSD for all functionals, where the largest difference is found for BHLYP with a difference of less than 0.1 eV. This is not quite obvious, then there are excitations in the test set that require at least two transitions for the description of the excited state and thus can not be properly represented as a single NTO transition, for example the 1 Ag excitation in E-butadiene, E-hexatriene and E-octatetraene (with γmax values of 0.83, 0.90 and 0.91 for RSCF-CV(∞)-DFT with B3LYP). Indeed, RSCF-CV(∞)-DFT performs better for these cases, for example for the three 1 Ag excitations mentioned before, its energies are 5.35 eV, 32

ACS Paragon Plus Environment

Page 33 of 50

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

4.35 eV and 3.62 eV when using LDA (see Table SI of the Supporting Information), while the SVD method yields 5.80 eV, 4.61 eV and 3.82 eV, which is further from the “best” values (5.08 eV, 4.15 eV and 3.55 eV). So it is clear that RSCF-CV(∞)-DFT correctly predicts the lowering of the energy due to the variational optimization of the transition matrix U. However, the difference for these three excitations between the SVD method and “best” values is not that large for the results obtained with LDA and so these excitations do not significantly change the RMSD value. Therefore, on average both methods (SVD and RSCFCV(∞)-DFT) yield very similar RMSD. Such a performance is similar for BP86 and B3LYP. In contrast the situation when using BHLYP, where these three excitation energies obtained with RSCF-CV(∞)-DFT (4.64 eV, 3.97 eV and 3.40 eV) are lower than the “best” values but the values obtained with the SVD method are approximately 1 eV higher than the “best” values (6.20 eV, 5.29 eV and 4.59 eV, see Table SIV in the Supporting Information). Such differences between RSCF-CV(∞)-DFT and SVD are reflected in the results of statistical averaging (see Table 3) over all excitations (0.42 eV against 0.47 eV). When we look now at the mean deviations, we see the average energies obtained with SVD coming closer to the “best” value than the ones obtained with RSCF-CV(∞)-DFT. This is neither a contradiction nor does it mean SVD to be superior to the more general RSCF-CV(∞)-DFT and can be understood as following. The DFT-based mean deviations are negative, meaning the excitation energies being in average lower than the “best” ones. SVD is a restricted version of RSCF-CV(∞)-DFT and limiting the transition as being a one NTO to one NTO transition, the corresponding energy is higher than the energy obtained without this restriction, thus less ‘too low’ when compared to the in average again higher “best” values. This “systematic underestimation” of triplet excitation energies has also been found for TD-DFT by Silva-Junior et al. 100 . But we would like to note that even with a “systematic underestimation” of the triplet excitation energies in meaning of a negative mean deviation we clearly also have excitations where the excitation energy compared to the “best” results is overestimated, as visible in Table 3, although the maximal overestimation

33

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

is always smaller than the maximal underestimation.

5

Conclusion

We have implemented and tested in this work several variational methods for the calculation of the triplet excited states. The most general one, RSCF-CV(∞)-DFT, completes the CVDFT-theory in providing the necessary optimization of the transition matrix U for triplet excitations. We validated our implementation by calculating the triplet excitation energies of the benchmark set out of Ref. 100 and comparison with its “best” estimate values as well as with corresponding TD-DFT values. With this possibility of variational optimization of the transition matrix U and orbital relaxation we introduced possible restrictions imposed on the transition matrix U. The goal of these restrictions is to develop a theoretical method that reproduces the ∆SCF-DFT results as close as possible but does not suffer from the variational collapse. A comparison of these ∆SCF-DFT-like methods (SOR, COL, and SVD) to the ∆SCF-DFT method has demonstrated that one of them, namely SVD, is numerically equivalent (within 0.1 eV) to ∆SCF-DFT for all types of functionals employed in this work. These small differences in excitation energies are attributed to the cut off in the relaxation energy after the second order terms. Thus, the implementation of the SVD method completes our previous work on theoretical justification of the ∆SCF-DFT method 47,48,51,70 . A comparison of the SVD method to the general RSCF-CV(∞)-DFT scheme has demonstrated that these two methods perform very similarly for transitions satisfactorily described with a single Slater determinant. But when there are at least two Slater determinants needed for an acceptable description of the excitation as we have seen for the 1 Ag transitions of E-butadiene, E-hexatriene and E-octatetraene, the restriction of the SVD method is too strong and the use of RSCF-CV(∞)-DFT is advantageous. In spite of this deficiency the stronger restricted methods (SOR and COL) are compu-

34

ACS Paragon Plus Environment

Page 34 of 50

Page 35 of 50

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

tationally cheaper than the general one, because they include fewer variational parameters in the transition matrix U. The excitation energy expression is also in all cases simpler than the general one. Moreover, SOR, COL, and SVD yield the orthogonal excited state KS-determinants, and additional configuration interaction between these states can be implemented. However, the last approach is beyond the scope of this publication.

6

Appendix

We shall consider here the two orbitals rotation in CV-DFT. In this case the transition matrix U has only one independent variational parameter: 



 0 γ  U= , −γ 0

(48)

so the unitary matrix Y is given by 



 cos(ηγ) sin(ηγ)  Y= , − sin(ηγ) cos(ηγ)

(49)

where the scaling parameter in this case is chosen from sin2 (ηγ) = 1

(50)

so that ηγ = π/2. This yields the excited state orbital φ′occ = φvir , and the density change is given by ∆ρ = φvir φvir − φocc φocc . Let us now analyze what happens if we truncate the Tailor expansion for the sine and cosine functions in Eq. (49). For the first order expansion, i.e. sin(x) = x and cos(x) = 1, the normalization condition yields x = 1 and det(Y) = 2. For the second order expansion, i.e. sin(x) = x and cos(x) = 1 − x2 /2, det(Y) = 1.25. For the third order expansion, i.e. 35

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

sin(x) = x − x3 /6 and cos(x) = 1 − x2 /2, it is impossible to normalize to one the truncated sin(x) function (see Figure 4). At least, this function has the maximum value in the range √ 2. In this case, det(Y) ≈ 0.889. For the fourth order expansion, we obtain [0:π/2]: det(Y) ≈ 0.917. For the fifth and sixth order expansions, the normalization condition for sin(x) = x − x3 /6 + x5 /120 yields x ≈ 1.491 and det(Y) ≈ 1.01. The results for other orders of the expansion are shown in Table 4. Eventually, for the ninth order expansion we obtain det(Y) ≈ 1.0001 and x ≈ 1.568, which is quite close to π/2. As we can see, the exponential

Figure 4: The truncated Tailor series expansion for the sine function in the range [0:π/2]. mapping is not quite accurate, when it is truncated, because Y† Y ≈ I for n < ∞.

Acknowledgement We would like to thank Dr. Yuriy Zinchenko for his help and valuable discussions. This work is supported by Natural Sciences and Engineering Research Council of Canada (NSERC). F.S. thanks the Eyes High program for the financial support. 36

ACS Paragon Plus Environment

Page 36 of 50

Page 37 of 50

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

Table 4: The results of truncating a Tailor series expansion for sine and cosine functions in Eq. (49) Order of expanstion 1 2 3 4 5 6 7 8 9 ∞

ηγ 1.0 1.0 1.414 1.414 1.491 1.491 1.568 1.568 1.568 1.571

det(Y) 2.0 1.25 0.8889 0.9167 1.0088 1.0062 0.9997 0.9997 1.0001 1.0

Supporting Information Available Full tables with the excitation energies, statistical values and detailed derivations of the key equations for RSCF-CV(∞)-DFT for the spin-flip transition. This material is available free of charge via the Internet at http://pubs.acs.org/.

References (1) Runge, E.; Gross, E. K. U. Density-Functional Theory for Time-Dependent Systems. Phys. Rev. Lett. 1984, 52, 997 – 1000. (2) Casida, M. E. In Recent Advances in Density Functional Methods. Part I ; Chong, D. P., Ed.; World Scientific, 1995; Chapter 5, pp 155 – 192. (3) van Gisbergen, S. J. A.; Snijders, J. G.; Baerends, E. J. A Density Functional Theory Study of Frequency-Dependent Polarizabilities and Van der Waals Dispersion Coefficients for Polyatomic Molecules. J. Chem. Phys. 1995, 103, 9347 – 9354. (4) Petersilka, M.; Gossmann, U.; Gross, E. K. U. Excitation Energies from TimeDependent Density-Functional Theory. Phys. Rev. Lett. 1996, 76, 1212 – 1215.

37

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

(5) Bauernschmitt, R.; Ahlrichs, R. Reatment of Electronic Excitations Within the Adiabatic Approximation of Time Dependent Density Functional Theory. Chem. Phys. Lett. 1996, 256, 454 – 464. (6) Furche, F. On the Density Matrix Based Approach to Time-Dependent Density Functional Response Theory. J. Chem. Phys. 2001, 114, 5982 – 5992. (7) Furche, F.; Ahlrichs, R. Adiabatic Time-Dependent Density Functional Methods for Excited State Properties. J. Chem. Phys. 2002, 117, 7433 – 7447. (8) Autschbach, J. Charge-Transfer Excitations and Time-Dependent Density Functional Theory: Problems and Some Proposed Solutions. ChemPhysChem 2009, 10, 1757 – 1760. (9) Ullrich, C. A. Time-Dependent Density-Functional Theory; Oxford Graduate Texts, 2012. (10) Dreuw, A.; Weisman, J. L.; Head-Gordon, M. Long-Range Charge-Transfer Excited States in Time-Dependent Density Functional Theory Require Non-Local Exchange. J. Chem. Phys. 2003, 119, 2943 – 2946. (11) Tozer, D. J. Relationship Between Long-Range Charge-Transfer Excitation Energy Error and Integer Discontinuity in Kohn-Sham Theory. J. Chem. Phys. 2003, 119, 12697 – 12699. (12) Peach, M. J. G.; Benfield, P.; Helgaker, T.; Tozer, D. J. Excitation Energies in Density Functional Theory: An Evaluation and a Diagnostic Test. J. Chem. Phys. 2008, 128, 044118. (13) Casida, M. E.; Jamorski, C.; Casida, K. C.; Salahub, D. R. Molecular Excitation Energies to High-Lying Bound States from Time-Dependent Density-Functional Response

38

ACS Paragon Plus Environment

Page 38 of 50

Page 39 of 50

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

Theory: Characterization and Correction of the Time-Dependent Local Density Approximation Ionization Threshold. J. Chem. Phys. 1998, 108, 4439 – 4449. (14) Tozer, D. J.; Handy, N. C. On the Determination of Excitation Energies Using Density Functional Theory. Phys. Chem. Chem. Phys. 2000, 2, 2117 – 2121. (15) Seidu, I.; Krykunov, M.; Ziegler, T. Applications of Time-Dependent and TimeIndependent Density Functional Theory to Rydberg Transitions. J. Phys. Chem. A 2015, 119, 5107 – 5116. (16) Stein, T.; Kronik, L.; Baer, R. Reliable Prediction of Charge Transfer Excitations in Molecular Complexes Using Time-Dependent Density Functional Theory. J. Am. Chem. Soc. 2009, 131, 2818 – 2820. (17) Baerends, E. J.; Gritsenko, O. V.; van Meer, R. The Kohn-Sham Gap, The Fundamental Gap and the Optical Gap: The Physical Meaning of Occupied and Virtual Kohn-Sham Orbital Energies. Phys. Chem. Chem. Phys. 2013, 15, 16408 – 16425. (18) Cave, R. J.; Zhang, F.; Maitra, N. T.; Burke, K. A Dressed TDDFT Treatment of the 2 1Ag States of Butadiene and Hexatriene. Chem. Phys. Lett. 2004, 389, 39 – 42. (19) Mazur, G.; Wlodarczyk, R. Application of the dressed time-dependent density functional theory for the excited states of linear polyenes. J. Comput. Chem. 2009, 30, 811 – 817. (20) Casida, M. E. Propagator Corrections to Adiabatic Time-Dependent DensityFunctional Theory Linear Response Theory. J. Chem. Phys. 2005, 122, 054111. (21) Romaniello, P.; Sangalli, D.; Berger, J. A.; Sottile, F.; Molinari, L. G.; Reining, L.; Onida, G. Double Excitations in Finite Systems. J. Chem. Phys. 2009, 130, 044108. (22) Gritsenko, O.; Baerends, E. J. Asymptotic Correction of the Exchange-Correlation

39

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

Kernel of Time-Dependent Density Functional Theory for Long-Range ChargeTransfer Excitations. J. Chem. Phys. 2004, 121, 655 – 660. (23) Hellgren, M.; Gross, E. K. U. Discontinuities of the Exchange-Correlation Kernel and Charge-Transfer Excitations in Time-Dependent Density-Functional Theory. Phys. Rev. A 2012, 85, 022514. (24) Levine, B. G.; Ko, C.; Quenneville, J.; Mart´ınez, T. J. Conical Intersections and Double Excitations in Time-Dependent Density Functional Theory. Mol. Phys. 2006, 104, 1039 – 1051. (25) Casida, M.; Huix-Rotllant, M. In Density-Functional Methods for Excited States; Ferr´e, N., Filatov, M., Huix-Rotllant, M., Eds.; Top. Curr. Chem.; Springer International Publishing, 2016; Vol. 368; pp 1 – 60. (26) Huix-Rotllant, M.; Ipatov, A.; Rubio, A.; Casida, M. E. Assessment of dressed timedependent density-functional theory for the low-lying valence states of 28 organic chromophores. Chem. Phys. 2011, 391, 120 – 129. (27) Slater, J. C. A Simplification of the Hartree-Fock Method. Phys. Rev. 1951, 81, 385 – 390. (28) Slater, J. C.; Wood, J. H. Statistical Exchange and the Total Energy of a Crystal. Int. J. Quantum Chem. 1970, 5, 3 – 34. (29) Slater, J. C. In Statistical Exchange-Correlation in the Self-Consistent Field ; L¨owdin, P.-O., Ed.; Adv. Quantum Chem.; Academic Press, 1972; Vol. 6; pp 1 – 92. ´ Transition Functional Method in the Density-Functional Theory. Phys. Rev. (30) Nagy, A. A 1996, 53, 3660 – 3663.

40

ACS Paragon Plus Environment

Page 40 of 50

Page 41 of 50

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

(31) Liu, T.; Han, W.-G.; Himo, F.; Ullmann, G. M.; Bashford, D.; Toutchkine, A.; Hahn, K. M.; Noodleman, L. Density Functional Vertical Self-Consistent Reaction Field Theory for Solvatochromism Studies of Solvent-Sensitive Dyes. J. Phys. Chem. A 2004, 108, 3545 – 3555. (32) Gavnholt, J.; Olsen, T.; Engelund, M.; Schiotz, J. Delta self-consistent field method to obtain potential energy surfaces of excited molecules on surfaces. Phys. Rev. B 2008, 78, 075441. (33) Oliveira, L.; Gross, E. K. U.; Kohn, W. Density-Functional Theory for Ensembles of Fractionally Occupied States. II. Application to the He Atom. Phys. Rev. A 1988, 37, 2821 – 2833. (34) Filatov, M.; Shaik, S. Spin-Restricted Density Functional Approach to the Open-Shell Problem. Chem. Phys. Lett. 1998, 288, 689 – 697. (35) Filatov, M.; Shaik, S. A spin-Restricted Ensemble-Referenced Kohn–Sham Method and its Application to Diradicaloid Situations. Chem. Phys. Lett. 1999, 304, 429 – 437. (36) Gidopoulos, N.; Papaconstantinou, P.; Gross, E. K. U. Spurious Interactions, and Their Correction, in the Ensemble-Kohn-Sham Scheme for Excited States. Phys. Rev. Lett. 2002, 88, 033003. (37) Pederson, M. R.; Klein, B. M. Improved theoretical methods for studies of defects in insulators: Application to the F center in LiF. Phys. Rev. B 1988, 37, 10319 – 10331. (38) Baruah, T.; Pederson, M. R. Density functional study on a light-harvesting carotenoidporphyrin-C60 molecular triad. J. Chem. Phys. 2006, 125, 164706. (39) Baruah, T.; Pederson, M. R. DFT Calculations on Charge-Transfer States of a

41

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

Carotenoid-Porphyrin-C60 Molecular Triad. J. Chem. Theory Comput. 2009, 5, 834 – 843. (40) Frank, I.; Hutter, J.; Marx, D.; Parrinello, M. Molecular Dynamics in Low-Spin Excited States. J. Chem. Phys. 1998, 108, 4060 – 4069. (41) Okazaki, I.; Sato, F.; Yoshihiro, T.; Ueno, T.; Kashiwagi, H. Development of a Restricted Open Shell Kohn-Sham Program and its Application to a Model Heme Complex. J. Mol. Struct.:THEOCHEM 1998, 451, 109 – 119. (42) Wu, Q.; Van Voorhis, T. Direct Optimization Method to Study Constrained Systems Within Density-Functional Theory. Phys. Rev. A 2005, 72, 024502. (43) Glushkov, V. N.; Levy, M. Optimized Effective Potential Method for Individual LowLying Excited States. J. Chem. Phys. 2007, 126, 174106. (44) Glushkov, V. N.; Assfeld, X. Doubly, Triply, and Multiply Excited States from a Constrained Optimized Effective Potential Method. J. Chem. Phys. 2010, 132, 204106. (45) Besley, N. A.; Gilbert, A. T. B.; Gill, P. M. W. Self-Consistent-Field Calculations of Core Excited States. J. Chem. Phys. 2009, 130, 124308. (46) Gilbert, A. T. B.; Besley, N. A.; Gill, P. M. W. Self-Consistent Field Calculations of Excited States Using the Maximum Overlap Method (MOM). J. Phys. Chem. A 2008, 112, 13164 – 13171. (47) Ziegler, T.; Seth, M.; Krykunov, M.; Autschbach, J.; Wang, F. On the Relation Between Time-Dependent and Variational Density Functional Theory Approaches for the Determination of Excitation Energies and Transition Moments. J. Chem. Phys. 2009, 130, 154102. (48) Cullen, J.; Krykunov, M.; Ziegler, T. The formulation of a self-consistent constricted

42

ACS Paragon Plus Environment

Page 42 of 50

Page 43 of 50

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

variational density functional theory for the description of excited states. Chem. Phys. 2011, 391, 11 – 18. (49) Ziegler, T.; Krykunov, M.; Cullen, J. The Implementation of a Self-Consistent Constricted Variational Density Functional Theory for the Description of Excited States. J. Chem. Phys. 2012, 136, 124107. (50) Krykunov, M.; Ziegler, T. Self-Consistent Formulation of Constricted Variational Density Functional Theory with Orbital Relaxation. Implementation and Applications. J. Chem. Theory Comput. 2013, 9, 2761 – 2773. (51) Park, Y. C.; Krykunov, M.; Ziegler, T. On the Relation Between Adiabatic Time Dependent Density Functional Theory (TDDFT) and the ∆SCF-DFT Method. Introducing a Numerically Stable ∆SCF-DFT Scheme for Local Functionals Based on Constricted Variational DFT. Mol. Phys. 2015, 113, 1636 – 1647. (52) Evangelista, F. A.; Shushkov, P.; Tully, J. C. Orthogonality Constrained Density Functional Theory for Electronic Excited States. J. Phys. Chem. A 2013, 117, 7378 – 7392. (53) Peng, B.; Van Kuiken, B. E.; Ding, F.; Li, X. A Guided Self-Consistent-Field Method for Excited-State Wave Function Optimization: Applications to Ligand-Field Transitions in Transition-Metal Complexes. J. Chem. Theory Comput. 2013, 9, 3933 – 3938. (54) Ziegler, T.; Krykunov, M.; Seidu, I.; Park, Y. C. In Density-Functional Methods for Excited States; Ferr´e, N., Filatov, M., Huix-Rotllant, M., Eds.; Top. Curr. Chem.; Springer International Publishing, 2016; Vol. 368; pp 61 – 95. (55) Ziegler, T.; Seth, M.; Krykunov, M.; Autschbach, J.; Wang, F. Is Charge Transfer Transitions Really too Difficult for Standard Density Functionals or are they just a Problem for Time-Dependent Density Functional Theory Based on a Linear Response Approach. J. Mol. Struct.:THEOCHEM 2009, 914, 106 – 109. 43

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

(56) Krykunov, M.; Grimme, S.; Ziegler, T. Accurate Theoretical Description of the 1La and 1Lb Excited States in Acenes Using the All Order Constricted Variational Density Functional Theory Method and the Local Density Approximation. J. Chem. Theory Comput. 2012, 8, 4434 – 4440. (57) Krykunov, M.; Seth, M.; Ziegler, T. Introducing Constricted Variational Density Functional Theory in its Relaxed Self-Consistent Formulation (RSCF-CV-DFT) as an Alternative to Adiabatic Time Dependent Density Functional Theory for Studies of Charge Transfer Transitions. J. Chem. Phys. 2014, 140, 18A502. (58) Ronca, E.; Angeli, C.; Belpassi, L.; De Angelis, F.; Tarantelli, F.; Pastore, M. Density Relaxation in Time-Dependent Density Functional Theory: Combining Relaxed Density Natural Orbitals and Multireference Perturbation Theories for an Improved Description of Excited States. J. Chem. Theory Comput. 2014, 10, 4014 – 4024. ´ Variational Density-Functional Theory for an Individual Excited (59) Levy, M.; Nagy, A. State. Phys. Rev. Lett. 1999, 83, 4361 – 4364. (60) Tsaune, A.;

Glushkov, V.;

Aprasyukhin, A. Towards an Optimal Zeroth-

Approximation for Excited State Energy. J. Mol. Struct.:THEOCHEM 1994, 312, 289 – 295. (61) Theophilou, A. K. The Energy Density Functional Formalism for Excited States. J. Phys. C: Solid State Physics 1979, 12, 5419 – 5430. (62) Fritsche, L. Generalized Kohn-Sham Theory for Electronic Excitations in Realistic Systems. Phys. Rev. B 1986, 33, 3976 – 3989. (63) G¨orling, A. Density-Functional Theory for Excited States. Phys. Rev. A 1996, 54, 3912 – 3915.

44

ACS Paragon Plus Environment

Page 44 of 50

Page 45 of 50

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

(64) G¨orling, A. Density-Functional Theory Beyond the Hohenberg-Kohn Theorem. Phys. Rev. A 1999, 59, 3359 – 3374. (65) Ayers, P. W.; Levy, M. Time-Independent (Static) Density-Functional Theories for Pure Excited States: Extensions and Unification. Phys. Rev. A 2009, 80, 012508. ´ Time-Independent Density-Functional Theory for (66) Ayers, P. W.; Levy, M.; Nagy, A. Excited States of Coulomb Systems. Phys. Rev. A 2012, 85, 042518. (67) Kowalczyk, T.; Yost, S. R.; Voorhis, T. V. Assessment of the ∆SCF Density Functional Theory Approach for Electronic Excitations in Organic Dyes. J. Chem. Phys. 2011, 134, 054128. (68) Gaudoin, R.; Burke, K. Lack of Hohenberg-Kohn Theorem for Excited States. Phys. Rev. Lett. 2004, 93, 173001. (69) Hohenberg, P.; Kohn, W. Inhomogeneous Electron Gas. Phys. Rev. 1964, 136, B864 – B871. (70) Zhekova, H.; Krykunov, M.; Autschbach, J.; Ziegler, T. Applications of Time Dependent and Time Independent Density Functional Theory to the First π to π ⋆ Transition in Cyanine Dyes. J. Chem. Theory Comput. 2014, 10, 3299 – 3307. (71) Martin, R. L. Natural Transition Orbitals. J. Chem. Phys. 2003, 118, 4775 – 4777. (72) Olsen, J. In Lecture Notes in Quantum Chemistry: European Summer School in Quantum Chemistry; Roos, B. O., Ed.; Springer Berlin Heidelberg: Berlin, Heidelberg, 1992; pp 37 – 87. (73) Helgaker, T.; Larsen, H.; Olsen, J.; Jrgensen, P. Direct optimization of the AO density matrix in Hartree-Fock and Kohn-Sham theories. Chem. Phys. Lett. 2000, 327, 397 – 403.

45

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

(74) Ziegler, T.; Krykunov, M. On the Calculation of Charge Transfer Transitions with Standard Density Functionals Using Constrained Variational Density Functional Theory. J. Chem. Phys. 2010, 133, 074104. (75) Seidu, I.; Zhekova, H. R.; Seth, M.; Ziegler, T. Calculation of Exchange Coupling Constants in Triply-Bridged Dinuclear Cu(II) Compounds Based on Spin-Flip Constricted Variational Density Functional Theory. J. Phys. Chem. A 2012, 116, 2268 – 2277. (76) Ziegler, T.; Seth, M.; Krykunov, M.; Autschbach, J. A Revised Electronic Hessian for Approximate Time-Dependent Density Functional Theory. J. Chem. Phys. 2008, 129, 184114. (77) Ziegler, T.; Krykunov, M.; Cullen, J. The Application of Constricted Variational Density Functional Theory to Excitations Involving Electron Transitions from Occupied Lone-Pair Orbitals to Virtual π ⋆ Orbitals. J. Chem. Theory Comput. 2011, 7, 2485 – 2491. (78) Wang, F.; Ziegler, T. Time-dependent density functional theory based on a noncollinear formulation of the exchange-correlation potential. J. Chem. Phys. 2004, 121, 12191 – 12196. (79) Wang, F.; Ziegler, T. The performance of time-dependent density functional theory based on a noncollinear exchange-correlation potential in the calculations of excitation energies. J. Chem. Phys. 2005, 122, 074109. (80) Wang, F.; Ziegler, T. Use of noncollinear exchange-correlation potentials in multiplet resolutions by time-dependent density functional theory. Int. J. Quantum Chem. 2006, 106, 2545 – 2550. (81) Peralta, J. E.; Scuseria, G. E.; Frisch, M. J. Noncollinear magnetism in density functional calculations. Phys. Rev. B 2007, 75, 125119.

46

ACS Paragon Plus Environment

Page 46 of 50

Page 47 of 50

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

(82) Ziegler, T.; Krykunov, M.; Autschbach, J. Derivation of the RPA (Random Phase Approximation) Equation of ATDDFT (Adiabatic Time Dependent Density Functional Ground State Response Theory) from an Excited State Variational Approach Based on the Ground State Functional. J. Chem. Theory Comput. 2014, 10, 3980 – 3986. (83) Seidu, I.; Krykunov, M.; Ziegler, T. Applications of Time-Dependent and TimeIndependent Density Functional Theory to Electronic Transitions in Tetrahedral d0 Metal Oxides. J. Chem. Theory Comput. 2015, 11, 4041 – 4053. (84) Huix-Rotllant, M.; Nikiforov, A.; Thiel, W.; Filatov, M. In Density-Functional Methods for Excited States; Ferr´e, N., Filatov, M., Huix-Rotllant, M., Eds.; Top. Curr. Chem.; Springer International Publishing, 2016; pp 445 – 476. (85) Filatov, M. Assessment of Density Functional Methods for Obtaining Geometries at Conical Intersections in Organic Molecules. J. Chem. Theory Comput. 2013, 9, 4526 – 4541. (86) Huix-Rotllant, M.; Natarajan, B.; Ipatov, A.; Muhavini Wawire, C.; Deutsch, T.; Casida, M. E. Assessment of noncollinear spin-flip Tamm-Dancoff approximation timedependent density-functional theory for the photochemical ring-opening of oxirane. Phys. Chem. Chem. Phys. 2010, 12, 12811 – 12825. (87) Huix-Rotllant, M. Improved correlation kernels for linear-response time-dependent density-functional theory. Theses, Universit´e de Grenoble, 2011. (88) Li, S. L.; Marenich, A. V.; Xu, X.; Truhlar, D. G. Configuration Interaction-Corrected TammDancoff Approximation: A Time-Dependent Density Functional Method with the Correct Dimensionality of Conical Intersections. J. Phys. Chem. Lett. 2014, 5, 322 – 328. (89) Dickson, R. M.; Ziegler, T. A density functional study of the electronic spectrum of permanganate. Int. J. Quantum Chem. 1996, 58, 681 – 687. 47

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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

(90) te Velde, G.; Bickelhaupt, F. M.; Baerends, E. J.; Fonseca Guerra, C.; van Gisbergen, S. J. A.; Snijders, J. G.; Ziegler, T. Chemistry with ADF. J. Comput. Chem. 2001, 22, 931 – 967. (91) Fonseca Guerra, C.; Snijders, J. G.; te Velde, G.; Baerends, E. J. Towards an Order-N DFT Method. Theor. Chem. Acc. 1998, 99, 391 – 403. (92) Baerends, E.; Ziegler, T.; Autschbach, J.; Bashford, D.; B´erces, A.; Bickelhaupt, F.; Bo, C.; Boerrigter, P.; Cavallo, L.; Chong, D.; Deng, L.; Dickson, R.; Ellis, D.; van Faassen, M.; Fan, L.; Fischer, T.; Fonseca Guerra, C.; Franchini, M.; Ghysels, A.; Giammona, A.; van Gisbergen, S.; G¨otz, A.; Groeneveld, J.; Gritsenko, O.; Gr¨ uning, M.; Gusarov, S.; Harris, F.; van den Hoek, P.; Jacob, C.; Jacobsen, H.; Jensen, L.; Kaminski, J.; van Kessel, G.; Kootstra, F.; Kovalenko, A.; Krykunov, M.; van Lenthe, E.; McCormack, D.; Michalak, A.; Mitoraj, M.; Morton, S.; Neugebauer, J.; Nicu, V.; Noodleman, L.; Osinga, V.; Patchkovskii, S.; Pavanello, M.; Philipsen, P.; Post, D.; Pye, C.; Ravenek, W.; Rodr´ıguez, J.; Ros, P.; Schipper, P.; van Schoot, H.; Schreckenbach, G.; Seldenthuis, J.; Seth, M.; Snijders, J.; Sol´a, M.; Swart, M.; Swerhone, D.; te Velde, G.; Vernooijs, P.; Versluis, L.; Visscher, L.; Visser, O.; Wang, F.; Wesolowski, T.; van Wezenbeek, E.; Wiesenekker, G.; Wolff, S.; Woo, T.; Yakovlev, A. ADF. 2012. (93) Van Lenthe, E.; Baerends, E. J. Optimized Slater-type Basis Sets for the Elements 1-118. J. Comput. Chem. 2003, 24, 1142 – 1156. (94) Vosko, S. H.; Wilk, L.; Nusair, M. Accurate Spin-Dependent Electron Liquid Correlation Energies for Local Spin Density Calculations: A Critical Analysis. Can. J. Phys. 1980, 58, 1200 – 1211. (95) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A 1988, 38, 3098 – 3100.

48

ACS Paragon Plus Environment

Page 48 of 50

Page 49 of 50

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

(96) Perdew, J. P. Density-Functional Approximation for the Correlation Energy of the Inhomogeneous Electron Gas. Phys. Rev. B 1986, 33, 8822 – 8824. (97) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B 1988, 37, 785 – 789. (98) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648 – 5652. (99) Becke, A. D. A New Mixing of Hartree-Fock and Local Density-Functional Theories. J. Chem. Phys. 1993, 98, 1372 – 1377. (100) Silva-Junior, M. R.; Schreiber, M.; Sauer, S. P. A.; Thiel, W. Benchmarks for Electronically Excited States: Time-Dependent Density Functional Theory and Density Functional Theory Based Multireference Configuration Interaction. J. Chem. Phys. 2008, 129, 104103. (101) Schreiber, M.; Silva-Junior, M. R.; Sauer, S. P. A.; Thiel, W. Benchmarks for Electronically Excited States: CASPT2, CC2, CCSD, and CC3. J. Chem. Phys. 2008, 128, 134110.

49

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

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 ACS Paragon Plus Environment

Page 50 of 50