Two-component hybrid time-dependent density functional theory within the TammDancoff approximation Michael Kühn and Florian Weigend Citation: The Journal of Chemical Physics 142, 034116 (2015); doi: 10.1063/1.4905829 View online: http://dx.doi.org/10.1063/1.4905829 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/142/3?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Phosphorescence lifetimes of organic light-emitting diodes from two-component time-dependent density functional theory J. Chem. Phys. 141, 224302 (2014); 10.1063/1.4902013 Performance of Tamm-Dancoff approximation on nonadiabatic couplings by time-dependent density functional theory J. Chem. Phys. 140, 054106 (2014); 10.1063/1.4862904 Valence excitation energies of alkenes, carbonyl compounds, and azabenzenes by time-dependent density functional theory: Linear response of the ground state compared to collinear and noncollinear spin-flip TDDFT with the Tamm-Dancoff approximation J. Chem. Phys. 138, 134111 (2013); 10.1063/1.4798402 Analytical Hessian of electronic excited states in time-dependent density functional theory with Tamm-Dancoff approximation J. Chem. Phys. 135, 014113 (2011); 10.1063/1.3605504 Nonadiabatic coupling vectors for excited states within time-dependent density functional theory in the Tamm–Dancoff approximation and beyond J. Chem. Phys. 133, 194104 (2010); 10.1063/1.3503765
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
THE JOURNAL OF CHEMICAL PHYSICS 142, 034116 (2015)
Two-component hybrid time-dependent density functional theory within the Tamm-Dancoff approximation Michael Kühn1 and Florian Weigend1,2,a) 1
Institut für Physikalische Chemie, Karlsruher Institut für Technologie, Kaiserstraße 12, 76131 Karlsruhe, Germany 2 Institut für Nanotechnologie, Karlsruher Institut für Technologie, Postfach 3640, 76021 Karlsruhe, Germany
(Received 11 November 2014; accepted 29 December 2014; published online 20 January 2015) We report the implementation of a two-component variant of time-dependent density functional theory (TDDFT) for hybrid functionals that accounts for spin-orbit effects within the Tamm-Dancoff approximation (TDA) for closed-shell systems. The influence of the admixture of Hartree-Fock exchange on excitation energies is investigated for several atoms and diatomic molecules by comparison to numbers for pure density functionals obtained previously [M. Kühn and F. Weigend, J. Chem. Theory Comput. 9, 5341 (2013)]. It is further related to changes upon switching to the local density approximation or using the full TDDFT formalism instead of TDA. Efficiency is demonstrated for a comparably large system, Ir(ppy)3 (61 atoms, 1501 basis functions, lowest 10 excited states), which is a prototype molecule for organic light-emitting diodes, due to its “spin-forbidden” triplet-singlet transition. C 2015 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4905829]
I. INTRODUCTION
Time-dependent density functional theory (TDDFT)1,2 has become a popular method for the calculation of electronic excitations of large systems. For compounds containing heavy elements, relativistic effects on spectra become important and need to be incorporated. Mostly this is limited to scalar relativistic effects since their incorporation does not change the structure of common quantum-chemical programs. However, in one-component scalar relativistic (1c) TDDFT, the influence of spin-orbit coupling (SOC) is neglected. Therefore, it is not possible to properly describe phenomena based on “spinforbidden” transitions, e.g., phosphorescence, which is the mechanism of interest in molecules applied in organic lightemitting diodes (OLEDs).3,4 When going beyond 1c TDDFT and including SOC, the most fundamental way is fully relativistic TDDFT based on the Dirac-Coulomb Hamiltonian, yielding four-component (4c) complex spinors instead of 1c real molecular orbitals (MOs).5–7 A more economic alternative of very similar accuracy for chemical problems is quasirelativistic TDDFT using twocomponent (2c) complex spinors.8,9 In the past few years, several groups reported 2c (Kramers-restricted) closed-shell TDDFT implementations using a noncollinear exchange-correlation (XC) kernel:8,9 Wang et al.8 in the Amsterdam density functional program package (ADF),10 Peng et al.9 in the Beijing density functional program suite (BDF),11–13 and Nakata et al.14 in the GAMESS (US) program package15 within the Tamm-Dancoff approximation (TDA).16 All three implementations are based on the all-electron 2c zeroth-order regular approximation (ZORA).17–19 Additionally, the all-electron exact 2c (X2C) a)Author to whom correspondence should be addressed. Electronic mail:
[email protected]. Fax: +49 721 608-47225.
0021-9606/2015/142(3)/034116/8/$30.00
approach20–29 has been applied within the implementation in the BDF program.30–32 Furthermore, Bast et al.7 claim that their 4c TDDFT implementation in the DIRAC program package33 based on the Dirac-Coulomb Hamiltonian is also fully operational for the X2C approach. Recently, we reported an implementation34 in the TURBOMOLE program suite,35 that can be used in combination with the X2C approach as well as 2c Dirac-HartreeFock effective core potentials (dhf-ECPs), which are available for elements beyond Kr36–42 and represent the most efficient choice to include relativistic effects, since only valence electrons have to be treated explicitly within the calculations, in contrast to all-electron methods. We demonstrated the high efficiency of our code by calculating the electronic spectrum of Au20 (2620 basis functions, polarized quadruple-ζ valence basis sets, 2c dhf-ECPs, lowest 57 excited states)34 and Au25(SCH3)−18 (1695 basis functions, polarized double-ζ valence basis sets, 2c dhf-ECPs, lowest 192 excited states),43 which are the largest systems studied so far at the 2c TDDFT level. Our implementation also facilitated a detailed investigation of the “spin-forbidden,” i.e., phosphorescent, transitions in an eight-membered set of iridium-containing organometallic complexes that are potential candidates for OLEDs.44 Since the emitting state in such compounds usually is of chargetransfer character (to be more precise: spin-orbit coupled triplet metal-to-ligand charge-transfer character, 3MLCT), TDDFT using pure density XC functionals largely underestimates the experimental phosphorescence energies. This is a widely known problem of TDDFT, which can—at least partly—be circumvented by inclusion of Hartree-Fock (HF) exchange. To the best of our knowledge, so far, only two 2c (and 4c) TDDFT/TDA implementations have been reported, that are able to treat HF exchange and therefore the important class of hybrid XC functionals using proper kernels (i.e., full derivatives of XC functionals and inclusion of HF exchange
142, 034116-1
© 2015 AIP Publishing LLC
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-2
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
in the coupling matrix): the implementations by Bast et al. in the DIRAC program package7 and Nakata et al. in the GAMESS (US) program suite,14 which both additionally work for the class of long-range corrected hybrid functionals. However, the systems investigated using the aforementioned implementations are of rather moderate size. To the best of our knowledge, the largest system studied so far is [UO2Cl4]2− (allelectron polarized triple-ζ valence basis sets, lowest 15 excited states),45 using the implementation of Bast et al.7 This indicates that the calculation of systems with up to about 150 atoms (a few thousand basis functions), e.g., the OLED molecules studied in Ref. 44, using the aforementioned 2c (and 4c) hybrid TDDFT/TDA codes is a challenge. The goal of this work is to extend our recently implemented efficient 2c (Kramers-restricted) closed-shell TDDFT code34 in the TURBOMOLE program suite to hybrid XC functionals using proper kernels, which is most conveniently and efficiently done within the 2c TDA making optimal use of already existing routines. Our code allows for treatments of systems with a few thousand basis functions on a single processor as well as reasonable scaling for moderately parallelized tasks. The paper is organized as follows. In Sec. II, we summarize the theory of 2c (Kramers-restricted) closed-shell TDDFT and TDA. In Sec. III, the accuracy of our 2c hybrid TDA implementation is assessed by comparing the excitation energies of the s → p transitions in Cd and Hg to the 2c TDA variant of Nakata et al.14 and the 4c TDDFT approach of Bast et al.7 as well as to experimental data. The excitation energies of the d → s transitions in Au+ and the lowest transitions in I2 and TlH are compared to experimental results. Finally, we demonstrate the efficiency of our code by calculating the phosphorescence energies of the popular OLED model compound Ir(ppy)3 using 1501 basis functions (polarized triple-ζ valence basis sets, 2c dhf-ECP, lowest 10 excited states), which is the largest system calculated so far using a 2c or 4c hybrid TDDFT/TDA algorithm.
The complex (Hermitian) orbital rotation Hessians are defined as Ai σ˜ aτ˜ j σ˜ ′ bτ˜ ′ = (ϵ aτ˜ − ϵ i σ˜ )δ i j δ ab δσ˜ σ˜ ′δτ˜ τ˜ ′ +Ci σ˜ aτ˜ j σ˜ ′ bτ˜ ′, Bi σ˜ aτ˜ j σ˜ ′ bτ˜ ′ = Ci σ˜ aτ˜ bτ˜ ′ j σ˜ ′,
(2) (3)
with the coupling matrix Ci σ˜ aτ˜ j σ˜ ′ bτ˜ ′ = ⟨Φi σ˜ Φbτ˜ ′|Φaτ˜ Φ j σ˜ ′⟩ + ⟨Φi σ˜ Φbτ˜ ′|fXC|Φaτ˜ Φ j σ˜ ′⟩ − cX⟨Φi σ˜ Φbτ˜ ′|Φ j σ˜ ′Φaτ˜ ⟩.
(4)
δ pq is the Kronecker delta and ϵ p σ˜ describes the energy of 2c spinor Φ p σ˜ (⃗r ) =
NBF α *c µβ p σ˜ + ξ µ (⃗r ) c ˜µ=1 , µ p σ
(5)
that is not necessarily an eigenfunction of Sˆ z = 21 σ z . ξ µ represent real basis functions (µ = 1,2, ..., NBF) and cαµ p σ˜ and c µβ p σ˜ are complex expansion coefficients. The first term in Eq. (4) is the Coulomb contribution in Dirac notation. The second term contains the matrix elements of the noncollinear XC kernel.5,6,8,9,47,48 The last term is the HF exchange contribution scaled by a hybrid mixing parameter cX. Note that there is no additional SOC contribution to the coupling matrix. SOC effects on excitations arise via the reference-state spinors and their energies only, as only one-electron relativistic potentials are included in the Hamiltonian (via ECPs, X2C, ZORA, . . . ). In the following, we only consider systems possessing a (Kramers-restricted) closed-shell ground state, for which Nα˜ = Nβ˜ = N2 and ϵ p α˜ = ϵ p β˜ since the Kramers partners Φ p α˜ and Φ p β˜ are either both occupied or unoccupied. The noncollinear XC kernel within the adiabatic approximation for (Kramersrestricted) closed-shell systems reads8,9 fXC(⃗r ,⃗r ′) =
δ2 EXC 1⊗1 δ ρ(⃗r )δ ρ(⃗r ′) δ2 EXC + δ(⃗r −⃗r ′)σ k ⊗ σ k . 2(⃗ δs r ) k=x, y, z
(6)
⃗ (⃗r )| is the noncollinear spin density with the magnes(⃗r ) = | m tization vector
II. TWO-COMPONENT TIME-DEPENDENT DENSITY FUNCTIONAL THEORY
The 2c non-Hermitian time-dependent Kohn-Sham (TDKS) eigenvalue problem that accounts for SOC effects on excitation energies ω n (n = 1,2,3, ...)2,7–9,34,46 reads
⃗ (⃗r ) = m
Nσ˜ σ ˜ i=1
Φ†i σ˜ (⃗r )⃗ σ Φi σ˜ (⃗r ).
(7)
The total density is given as ′ ′ * Ai σ˜ aτ˜ j σ˜ bτ˜ ∗ ′ ′ , Bi σ˜ aτ˜ j σ˜ bτ˜
1 = ωn * ,0
Bi σ˜ aτ˜ j σ˜ ′ bτ˜ ′ + * X j σ˜ ′ bτ˜ ′, n + ∗ Ai σ˜ aτ˜ j σ˜ ′ bτ˜ ′ - , Yj σ˜ ′ bτ˜ ′, n -
0 + * X i σ˜ aτ,˜ n + , −1- , Yi σ˜ aτ,˜ n -
ρ(⃗r ) =
Nσ˜ σ ˜ i=1
(1)
where the Einstein summation convention is utilized and indices i, j, ... are used for occupied, a, b, ... for virtual, and p, q, ... for general molecular spinors. Greek indices σ, ˜ τ, ˜ ˜ refer to 2c time-reversal symmetry-adapted σ ˜ ′, τ˜ ′ ∈ { α, ˜ β} Kramers partners describing Nα˜ moment-up (α) ˜ and Nβ˜ ˜ electrons. It is noted that Eq. (1) also holds moment-down ( β) for 4c TDDFT.5–7,9
Φ†i σ˜ (⃗r )Φi σ˜ (⃗r ).
(8)
In Eq. (6), EXC is the nonrelativistic XC energy depending on both ρ and s. The left side of the tensor product refers to electron 1 (indices i and a in Eq. (4)), whereas the right side refers to electron 2 (indices b and j). To arrive at a convenient formulation of 2c TDDFT, that is able to make optimal use of the existing 1c TDDFT code49,50 in the TURBOMOLE program suite,35 the goal is to reduce the two-superdimensional non-Hermitian eigenvalue problem in Eq. (1) to a one-superdimensional Hermitian eigenvalue problem. In the following, we will discuss two possibilities: 2c
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-3
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
TDDFT within Ziegler’s approximation8,9,34 in Sec. II A and the 2c TDA in Sec. II B.16 A. Two-component time-dependent density functional theory within Ziegler’s approximation
In the case of pure density functionals, i.e., cX = 0 in the coupling matrix (Eq. (4)), an approximate 2c TDKS Hermitian eigenvalue problem for the square of the excitation energies of (Kramers-restricted) closed-shell systems was derived by Wang, Ziegler, van Lenthe, van Gisbergen, and Baerends8 as well as Peng, Zou, and Liu.9 It reads Ωi σ˜ aτ˜ j σ˜ ′ bτ˜ ′ (X +Y )′j σ˜ ′ bτ˜ ′, n = ω2n (X +Y )i′σ˜ aτ,˜ n ,
(9)
Ωi σ˜ aτ˜ j σ˜ ′ bτ˜ ′ = (ϵ aτ˜ − ϵ i σ˜ )2δ i j δ ab δσ˜ σ˜ ′δτ˜ τ˜ ′ √ + 2 ϵ aτ˜ − ϵ i σ˜ Ci σ˜ aτ˜ j σ˜ ′ bτ˜ ′ ϵ bτ˜ ′ − ϵ j σ˜ ′.
(10)
with
Note that the transformation from Eq. (1) to Eq. (9) is only approximate if SOC is present since the (polarized) XC kernel (Eq. (6)) breaks time-reversal (or Kramers) symmetry, which has to be exploited in the transformation. However, Wang et al. claim that test calculations show errors always less than 0.05 eV compared to a calculation according to Eq. (1).8 Peng et al. state that this approximation may not work if either the excitation energy is not close to the orbital energy difference or if there are highly inhomogeneous induced fields (e.g., localized non-fully filled f spinors of a formally closedshell system).9 In case the XC contribution to the coupling matrix in Eq. (4) is neglected, i.e., only the Coulomb contribution is included, it has been shown that Eq. (9) is obtained by an exact transformation also in the presence of SOC.51 In all cases, time-reversal symmetry has to be exploited so that the transformation is only applicable when (Kramers-restricted) closed-shell systems are considered. For open-shell systems, it is not possible to apply an analogous transformation that arrives at Eq. (9). At the 1c level in contrast, the transformation to the 1c analog of Eq. (9) is exact for all cases, as long as HF exchange is excluded, i.e., cX = 0. Recently, we reported the implementation of the formalism presented in this section in the ESCF excited state module49,50 of the TURBOMOLE program suite,35 see Ref. 34 for details.
of TDA is that transition moments (e.g., the transition dipole moment) and therefore oscillator strengths and radiative lifetimes are ill-defined due to the lack of gauge invariance.50 In practice, this often leads to poor results for TDA oscillator strengths.2,52–55 It is noted in passing that in contrast to Eq. (9), Eq. (11) is also applicable when open-shell systems are considered, however, the expression for the noncollinear XC kernel (Eq. (6)) becomes quite involved.9,48 Both 2c TDA (Eq. (11)) and 2c TDDFT within Ziegler’s approximation (Eq. (9)), that are based on a noncollinear XC kernel, have the correct 1c limit, where moment-up/moment˜ electrons become ordinary spin-up/spin-down down (α/ ˜ β) (α/ β) electrons. As a result, the coupling matrix is real, the MOs Φ pσ are eigenfunctions of Sˆ z , and the XC kernel is collinear. Besides this correct 1c limit, all of the above equations can be easily extended to the fully relativistic 4c formalism based on the Dirac-Coulomb Hamiltonian, simply by replacing all 2c quantities by their 4c analogs.5–7,9 2c TDA (and also 2c TDDFT within Ziegler’s approximation) is computationally more demanding than 1c TDA. For the 2c (Kramers-restricted) closed-shell treatment, we expect an additional factor of eight in both computation time and required memory compared to the 1c closed-shell formalism, due to the larger size of the excitation vectors X i σ˜ aτ,˜ n (n = const): a factor of four arises since the number of both occupied and virtual molecular orbitals/spinors is doubled and an additional factor of two is due to the introduction of imaginary parts (complex spinors). It is noted that by exploitation of double group symmetry together with time-reversal symmetry, the computational cost for the 2c (and also 4c) treatment can be reduced,56 as done in Refs. 7, 8, and 30–32. The 2c TDA presented in this section was implemented in the ESCF excited-state module34,49,50 of the TURBOMOLE program suite,35 mainly by extending the existing 1c routines to the 2c formalism and making use of our recently implemented 2c TDDFT code.34 Details of the implementation are given in Figure 1 of the supplementary material.57 Our code is parallelized employing the technique by van Wüllen58 in the same way as for the 1c case.
III. APPLICATION A. Computational details
B. Two-component Tamm-Dancoff approximation
To arrive at a 2c excited state theory that is able to treat HF exchange (and therefore hybrid functionals), i.e., cX , 0 in the coupling matrix (Eq. (4)), most conveniently the TDA,16 that neglects Y in Eq. (1), is used. The 2c Hermitian eigenvalue problem for excitation energies of (Kramers-restricted) closedshell systems within the TDA reads14 Ai σ˜ aτ˜ j σ˜ ′ bτ˜ ′ X j σ˜ ′ bτ˜ ′, n = ωTDA X i σ˜ aτ,˜ n . n
(11)
In the case of cX = 1 and no XC contribution to the coupling matrix (Eq. (4)), 2c TDA is identical to the 2c configuration interaction singles (CIS) method. Compared to the full twosuperdimensional formalism in Eq. (1), usage of TDA leads to systematically higher excitation energies. A disadvantage
All ground-state closed-shell structures of diatomic molecules were optimized using 1c DFT with the generalized gradient approximation (GGA) functional of Becke59 and Perdew60 (BP86) in combination with polarized quadruple-ζ valence basis sets (dhf-QZVP-2c).61 The dhf-bases are optimized for the usage in connection with dhf-ECPs, which account for scalar relativistic effects in 1c calculations and additionally for SOC effects in 2c calculations. For I a pseudopotential (dhf-ECP-28)62 covering the inner 28 electrons was used, for Tl a pseudopotential (dhf-ECP-60)39 covering the inner 60 electrons. The structure of the lowest excited triplet state of Ir(ppy)3 was optimized without symmetry constraints using 1c ground-state DFT in combination with Becke’s threeparameter hybrid functional with Lee-Yang-Parr correlation
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-4
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
(B3LYP)63 and polarized triple-ζ valence basis sets (def2TZVP).64 See also Ref. 44 for further information on the choice of the level of theory for the structure optimization. For Ir, a Wood-Boring-type pseudopotential (mwb-ECP-60)65 was used, which covers the inner 60 electrons and accounts for scalar relativistic effects only. Coordinates of all optimized structures are available in Tables I-III of the supplementary material.57 Vertical excitation energies (VEEs) were computed starting from a (Kramers-restricted) singlet ground state at both the 1c and 2c TDA16 as well as the TDDFT level34,49,50, using the Perdew-Wang parameterization of the local spin-density approximation (LSDA),66 the GGA functional BP86,59,60 and the hybrid functional B3LYP.63 Additionally, for Ir(ppy)3 at the 1c level, phosphorescence energies were calculated directly by vertical de-excitations from the lowest excited triplet state (T1) to the singlet ground state (S0) using the spin-flip TammDancoff approximation (SF-TDA).67 For atoms and diatomic molecules, dhf-QZVP-2c basis sets were used, for Ir(ppy)3 basis sets of polarized triple-ζ valence quality (dhf-TZVP2c).61 For Cd, a pseudopotential (dhf-ECP-28)38 covering the inner 28 electrons, for Au, Hg, and Ir pseudopotentials (dhfECP-60)37,38 covering the inner 60 electrons were used. The iterative Davidson procedure in the TDDFT/TDA calculations is stopped when the Euclidean norm of the residual vector is smaller than 10−5 a.u. In all calculations using pure density functionals, the resolution of the identity (RI) approximation68,69 was employed. The applied auxiliary basis sets70 are suited for bases optimized for dhf-ECPs.71 Reference state energies and density matrices were converged to 10−7 a.u. Fine quadrature grids of size m4 (Ref. 72) were employed. 1. Assessment of accuracy: Excitation energies and fine-structure splittings of selected closed-shell atoms, ions, and diatomic molecules
Accuracy of the method was studied for the s → p transitions in the atoms Cd and Hg, the d → s transitions in Au+ as well as for the lowest transitions in I2 and TlH. For all cases, we consider the influence of the admixture of HF exchange by
comparison with the (pure) GGA functional BP86 and relate it to both the effect of the gradient correction by comparison with the LDA functional LSDA and to the influence of the TDA by comparison to the full TDDFT formalism (for the LSDA functional). Further, results are compared to other implementations of hybrid TDDFT7,14 and to the experiment, as far as available. First, we investigate the s → p transitions in the Cd atom as well as the Hg atom in detail, see Table I. The quantities of interest are the average vertical excitation energies for the triplet excitation (3 P), the respective splittings due to SOC, and the energy difference to the singlet (1 P1) excitation (“shift”). For both atoms and all quantities, the differences between B3LYP and BP86 are very small. The maximum difference is observed for the shift of the 1 P1 transition in case of Cd, ca. 0.1 eV. Energies for the triplet splittings are essentially the same for these two as well as for the other methods. For LSDA, comparably large differences to the BP86 treatment for the average 3 P VEE are observed, typically 0.3 eV. Employing the full TDDFT formalism instead of the TDA variant (at LSDA level) leads to VEEs that are slightly smaller (by less than 0.1 eV for the average 3 P VEE and by somewhat more than 0.1 eV for the 1 P1 VEE). Further, 2c (average) VEEs and 1c VEEs (given in parentheses in Table I) are very similar for all cases with maximum differences amounting to 0.06 eV. Our B3LYP results very well agree with those obtained within a 2c ZORA approach, Ref. 14 and at the 4c Dirac-Coulomb level, Ref. 7. In case of Cd, our results are more similar to those of Ref. 14; in case of Hg, they are always between those from the two references. This is not unexpected as in Ref. 7 a 4c (Dirac) code is used, in Ref. 14 a 2c ZORA approach. Our method, a 2c approach with ECPs fitted to 4c Dirac results, could be regarded as a mixture of both. Also, agreement with experimental data is very good, deviations for splittings are below 0.03 eV, for excitation energies up to ca. 0.15 eV, but it has to be noted that this is not better than the numbers obtained with the significantly more economic BP86 functional. The d → s excitation in Au+ reveals a similar picture concerning the sublevel splitting of the 3 D transition, see Table II; almost identical numbers are obtained with the four types of calculations that moreover agree well with the experiment. A marked difference to Cd or Hg is observed for the
TABLE I. Two-component (2c) VEEs for s → p transitions, fine-structure splittings, and shifts (given relatively to average of triplet state) of Cd and Hg. Unless noted otherwise the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with dhf-QZVP-2c basis sets and dhf-ecp-2c spin-orbit effective core potentials (ECPs). Values in parentheses are scalar relativistic (1c) results for the 3 P and 1 P states. Degeneracies g J = (2J + 1) are also given. All values are in eV. TDDFT results are from Ref. 34. 2c TDA VEEs using ZORA obtained by Nakata et al.14 are listed in column “B3LYP,14” 4c TDDFT VEEs based on the Dirac-Coulomb Hamiltonian by Bast et al.7 in column “B3LYP.7” Experimental results are from Ref. 76. gJ
B3LYP
BP86
LSDA
LSDA (TDDFT)
B3LYP14
B3LYP7
Expt.
Average of triplet state Splitting of 3 P0 Splitting of 3 P1 Splitting of 3 P2 Shift of 1 P1
1 3 5 3
3.921 (3.924) −0.134 −0.069 0.068 1.514 (1.513)
3.909 (3.912) −0.136 −0.070 0.069 1.620 (1.619)
4.236 (4.240) −0.139 −0.072 0.071 1.425 (1.424)
4.195 (4.199) −0.141 −0.074 0.072 1.271 (1.269)
3.91 −0.14 −0.07 0.07 1.63
3.81 −0.15 −0.08 0.07 1.33
3.874 −0.140 −0.073 0.072 1.543
Average of triplet state Splitting of 3 P0 Splitting of 3 P1 Splitting of 3 P2 Shift of 1 P1
1 3 5 3
5.073 (5.119) −0.482 −0.268 0.257 1.358 (1.341)
5.151 (5.201) −0.483 −0.266 0.256 1.436 (1.422)
5.452 (5.506) −0.491 −0.277 0.264 1.251 (1.236)
5.411 (5.470) −0.492 −0.287 0.270 1.156 (1.131)
5.17 −0.41 −0.23 0.22 1.45
5.03 −0.52 −0.30 0.28 1.24
5.181 −0.514 −0.295 0.280 1.523
Atom Cd
Hg
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-5
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
TABLE II. Two-component (2c) VEEs for d → s transitions, fine-structure splittings, and shifts (given relatively to average of triplet state) of Au+. Unless noted otherwise, the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with the dhf-QZVP-2c basis set and dhf-ecp-2c spin-orbit effective core potential (ECP). Values in parentheses are scalar relativistic (1c) results for the 3 D and 1 D states. Degeneracies g J = (2J + 1) are also given. All values are in eV. Experimental results are from Ref. 76.
Average of triplet state Splitting of 3 D3 Splitting of 3 D2 Splitting of 3 D1 Shift of 1 D 2
gJ
B3LYP
BP86
LSDA
LSDA (TDDFT)
Expt.
7 5 3 5
1.766 (1.903) −0.451 −0.015 1.077 1.581 (1.078)
1.546 (1.677) −0.459 −0.007 1.083 1.621 (1.120)
1.554 (1.687) −0.442 −0.020 1.065 1.550 (1.039)
1.498 (1.642) −0.454 −0.020 1.092 1.531 (1.002)
2.288 −0.423 −0.101 1.155 1.385
average triplet excitation energy: for both Cd and Hg, values for B3LYP and BP86 were almost identical and LSDA significantly differed from those, but for Au+ BP86 and LSDA are very similar and B3LYP is larger by ca. 0.2 eV, thus being closer to the experimental value; nevertheless, the difference to the experiment, ca. 0.5 eV, is significant. This error can be traced back to the different self-interaction experienced by the 5 d and the 6 s levels, as noted in the work of Wang et al.8 The 5 d levels are more compact than the 6 s level and therefore are more prone to unphysical self-interaction. As a result, the energies of the former are overestimated compared to the energy of the latter, consequently leading to an underestimation of the 5 d → 6 s excitation energies, which in zeroth order are given by the differences of the respective orbital/spinor energies. Concerning functional changes, see above, the same trends are obtained for the calculation at the 1c level, but here the 3 D VEEs are larger by ca. 0.15 eV throughout. This difference may be explained by the splitting of the fully occupied 5 d orbitals by SOC, which is not present in the 1c calculation. The 3 D excitation in the 2c case mainly involves the d 5/2 level which is higher in energy than the d shell in the 1c case.
For the diatomic molecules I2 and TlH, see Tables III and IV, B3LYP energies for excitations to the Π states (or their 2c analogs, respectively) are slightly (ca. 0.1–0.2 eV) higher than those obtained with BP86; this is somewhat more than changes from BP86 to LSDA/TDA or from the latter to LSDA/TDDFT (both below 0.1 eV). Excitations to the Σ states of I2 are very similar for B3LYP and BP86, here larger differences (0.2 eV) are observed between the functionals BP86 and LSDA. This holds for both 1c and 2c treatments. Differences to the few available experimental data are smallest for B3LYP in all cases, excitations to the Π states are underestimated by up to 0.4 eV, the experimental excitation energy to the 1u (3Σu+) state of I2 is almost exactly met (−0.03 eV). Overall, for the systems discussed so far, the gain from the comparable costly (see below) B3LYP is rather moderate; the accuracy of level splittings is very satisfying also for the cheaper non-hybrid (pure) density functionals, the accuracy of averaged 2c excitation energies (and 1c excitation energies) is systematically improved by B3LYP only for the transitions in Au+ as well as the Π states of the diatomic molecules.
TABLE III. Scalar relativistic (1c) (states expressed in the Λ − Σ notation) as well as two-component (2c) VEEs (states expressed in the ω − ω notation) of I2. Unless noted otherwise the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with the dhf-QZVP-2c basis set and dhf-ecp-2c spin-orbit effective core potential (ECP). Degeneracies are also given. All values are in eV. TDDFT results are from Ref. 34. Excited state Scalar relativistic (1c) 3Π u 1Π u 3Π g 1Π g 3Σ + u Two-component (2c) 2u 1u 0−u 0+u 1u 2g 1g 0−g 0+g 1g 0−u 1u
Composition/contribution
Degeneracy
0+u , 0−u , 1u , 2u 1u 0+g , 0−g , 1g , 2g 1g 0−u , 1u 3Π
u
3Π
u
3Π
u
3Π
u
1Π
u 3Π g 3Π g 3Π g 3Π g 1Π g 3Σ + u 3Σ + u
B3LYP
BP86
LSDA
LSDA (TDDFT)
6 2 6 2 3
1.60 2.21 3.21 3.80 4.35
1.49 2.11 3.03 3.63 4.39
1.52 2.08 3.08 3.61 4.59
1.50 2.00 3.07 3.56 4.55
2 2 1 1 2 2 2 1 1 2 1 2
1.32 1.49 1.90 2.03 2.30 2.97 3.14 3.63 3.65 4.03 4.51 4.54
1.22 1.39 1.81 1.93 2.22 2.80 2.97 3.44 3.46 3.85 4.53 4.56
1.25 1.41 1.86 1.97 2.21 2.86 3.01 3.51 3.53 3.86 4.73 4.75
1.22 1.38 1.84 1.94 2.14 2.84 2.99 3.49 3.51 3.81 4.69 4.71
Expt.
1.6977 1.8477 2.1377 2.3778 2.4977
4.5779
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-6
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
TABLE IV. Scalar relativistic (1c) (states expressed in the Λ − Σ notation) as well as two-component (2c) VEEs (states expressed in the ω − ω notation) of TlH. Unless noted otherwise, the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with dhf-QZVP-2c basis sets and the dhf-ecp-2c spin-orbit effective core potential (ECP) for Tl. Degeneracies are also given. All values are in eV. TDDFT results are from Ref. 34. Experimental results are from Ref. 79. Excited state Scalar relativistic (1c) 3Π 1Π Two-component (2c) 0− 1 0+ 1 2
Composition/contribution
Degeneracy
B3LYP
BP86
LSDA
LSDA (TDDFT)
0+, 0−, 1, 2 1
6 2
2.24 2.75
2.17 2.69
2.15 2.59
2.14 2.54
1 2 1 2 2
1.97 2.15 2.17 2.69 2.99
1.91 2.08 2.10 2.63 2.95
1.92 2.08 2.14 2.65 2.90
1.89 2.05 2.09 2.64 2.86
3Π 3Π 3Π 1Π 3Π
2. Demonstration of efficiency: Phosphorescence energies of Ir(ppy)3
Ir(ppy)3, see Figure 1, is a prototype molecule for OLEDs; the lowest transition of this molecule is “spin-forbidden” with pronounced charge-transfer (CT) character. Thorough studies of this system with pure density functionals were presented recently.44 Calculations for this molecule at the B3LYP level were carried out in the present work both for demonstration of efficiency and in order to investigate the effect of admixing HF exchange. For this transition, B3LYP significantly improves the accuracy of the average triplet excitation energy (phospho-
FIG. 1. Structure of Ir(ppy)3 (fac-tris[2-(2-pyridinyl-κ N )phenyl-κC]iridium) in Lewis notation.
Expt.
2.20 3.02 3.00
rescence energy), see Table V. With B3LYP, the difference to the experimental value of 2.43 eV (Refs. 73 and 74) amounts to 0.3 eV, for BP86, LSDA/TDA and LSDA/TDDFT more than 0.6 eV are obtained. This improvement is expected, as HF copes much better with CT states than DFT, thus leading to significantly higher excitation energies. Consequently, the same trend is observed for the 1c treatments, which show excitation energies higher by ca. 0.1 eV throughout, thus being even closer to the experiment. It is noted that B3LYP despite being superior to LDA and GGA functionals still has limitations concerning the calculation of CT transitions. It is commonly believed that within TDDFT, the latter are described most accurately by long-range corrected hybrid functionals. Using those, the 1c triplet excitation energy for Ir(ppy)3 was calculated to 2.3 eV, which is slightly closer to the experimental value than the B3LYP result.75 A gain in accuracy similar to that of B3LYP can be obtained in a more economic way, namely, by correcting the LSDA numbers with excitation energies obtained from SF-TDA,48,67 as described previously,44,67 see second column from the right in Table V. The size of the triplet splitting is overestimated by all methods. The difference between the highest and the lowest sublevel amounts to 0.010 eV for the experiment, to 0.024 eV for B3LYP, and to 0.036 eV for BP86. LSDA yields 0.032 and LSDA/TDDFT 0.028 eV. Apparently, B3LYP is closest to the experiment, but it has to be noted that the triplet sublevel splitting is very sensitive to the structure parameters of the molecule, as demonstrated previously,44 and thus a reliable comparison is difficult. The increase of computational costs compared to, e.g., BP86 for the improvement of B3LYP is comparably large,
TABLE V. Two-component (2c) phosphorescence energies and fine-structure splittings of triplet substates TSO 1, m (m = 1, 2, 3) (given relatively to average of triplet state) of Ir(ppy)3. Unless noted otherwise, the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with dhf-TZVP-2c basis sets and the dhf-ecp-2c spin-orbit effective core potential (ECP) for Ir. Values in parentheses are scalar relativistic (1c) results. All values are in eV. TDDFT results are from Ref. 34, SF-TDA results from Ref. 44. In case of the latter, the 1c phosphorescence energies were calculated directly as de-excitations from T1 to S0 using the 1c SF-TDA. The absolute energies of the triplet substates TSO 1, m (obtained with 2c TDDFT) were shifted by the difference between the 1c SF-TDA and the 1c TDDFT results. Experimental results are from Refs. 73 and 74.
Average of triplet state Splitting of TSO 1,1 Splitting of TSO 1,2 Splitting of TSO 1,3
B3LYP
BP86
LSDA
LSDA (TDDFT)
LSDA (SF-TDA)
Expt.
2.1488 (2.2114) −0.0090 −0.0064 0.0154
1.7993 (1.8842) −0.0137 −0.0084 0.0221
1.7831 (1.8675) −0.0125 −0.0073 0.0198
1.7770 (1.8618) −0.0110 −0.0063 0.0173
2.1858 (2.2706) −0.0110 −0.0063 0.0173
2.4311 −0.0040 −0.0023 0.0063
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-7
M. Kühn and F. Weigend
at least in case of TURBOMOLE, where non-hybrid DFT is implemented in a very efficient way by making use of efficient numerical grids as well as the RI approximation for the Coulomb part. This makes HF exchange by far the dominant part in hybrid DFT treatments, for 1c and for 2c calculations and for the ground state as well as for excited states. The CPU time for the determination of the lowest ten excitations at the 2c B3LYP level using 1501 basis functions amounts to 30 days and 16 h on a single Intel Xeon X5650 2.67 GHz CPU; for convergence of the excitation vectors, 18 iterations were needed. This effort is much higher than at the 2c BP86 level, 18 h and 25 min for 11 iterations. Thus, the time per iteration is higher for B3LYP than for BP86 by about a factor of 24.4. Almost the same factor between B3LYP and BP86 for the CPU times per iteration is obtained for 1c calculations, 23.2. Here, nine iterations are needed for BP86 and 14 for B3LYP. The ratio between CPU times per iteration for 2c and 1c treatments amounts to a factor of 8.2 for B3LYP and to 7.8 for BP86, which is close to the expected value of eight, see Sec. II B.
IV. CONCLUSION
We presented the implementation of a two-component variant of hybrid TDDFT. It allows for routine treatments of spin-orbit effects in the calculation of excitations of closedshell molecules with about 50 atoms. Concerning efficiency and treatable system size, this is a remarkable advance over already existing codes.7,14 In contrast to the latter, twocomponent effective core potentials are used here, which at least in part is the reason for the higher efficiency. Results for excitation energies and sublevel splittings of atoms and small molecules obtained with this approach usually are between those published for the four-component hybrid TDDFT based on the Dirac-Coulomb Hamiltonian7 and the two-component hybrid TDDFT within the ZORA approach.14 Significant improvement over GGA or LDA functionals is observed mainly for the average values of excitation energies in case of charge-transfer transitions, which comes along with an increase of CPU times by a factor of about 25, very similar to one-component treatments. Excitation energies for other cases and sublevel splittings in general are similar for B3LYP and the other functionals. ACKNOWLEDGMENTS
M.K. is funded by the Carl Zeiss Foundation. The author thanks TURBOMOLE GmbH for financial support. 1Time-Dependent
Density Functional Theory, Lecture Notes in Physics, edited by M. A. L. Marques, C. A. Ullrich, F. Nogueira, A. Rubio, K. Burke, and E. K. U. Gross (Springer, Heidelberg, 2006). 2F. Furche, J. Chem. Phys. 114, 5982 (2001). 3M. A. Baldo, D. F. O’Brien, Y. You, A. Shoustikov, S. Sibley, M. E. Thompson, and S. R. Forrest, Nature 395, 151 (1998). 4R. C. Evans, P. Douglas, and C. J. Winscom, Coord. Chem. Rev. 250, 2093 (2006). 5J. Gao, W. Liu, B. Song, and C. Liu, J. Chem. Phys. 121, 6658 (2004). 6J. Gao, W. Zou, W. Liu, Y. Xiao, D. Peng, B. Song, and C. Liu, J. Chem. Phys. 123, 054102 (2005). 7R. Bast, H. J. A. Jensen, and T. Saue, Int. J. Quantum Chem. 109, 2091 (2009).
J. Chem. Phys. 142, 034116 (2015) 8F.
Wang, T. Ziegler, E. van Lenthe, S. van Gisbergen, and E. J. Baerends, J. Chem. Phys. 122, 204103 (2005). Peng, W. Zou, and W. Liu, J. Chem. Phys. 123, 144101 (2005). 10See http://www.scm.com for Amsterdam Density Functional (ADF) version 2013.01a. 11W. Liu, G. Hong, D. Dai, L. Li, and M. Dolg, Theor. Chem. Acc. 96, 75 (1997). 12W. Liu, F. Wang, and L. Li, J. Theor. Comput. Chem. 2, 257 (2003). 13W. Liu, F. Wang, and L. Li, Recent Advances in Relativistic Molecular Theory, Recent Advances in Computational Chemistry edited by K. Hirao and Y. Ishikawa (World Scientific, Singapore, 2004), Vol. 5, p. 257. 14A. Nakata, T. Tsuneda, and K. Hirao, J. Chem. Phys. 135, 224106 (2011). 15M. W. Schmidt, K. K. Baldridge, J. A. Boatz, S. T. Elbert, M. S. Gordon, J. H. Jensen, S. Koseki, N. Matsunaga, K. A. Nguyen, S. J. Su, T. L. Windus, M. Dupuis, and J. A. Montgomery, J. Comput. Chem. 14, 1347 (1993). 16S. Hirata and M. Head-Gordon, Chem. Phys. Lett. 314, 291 (1999). 17E. van Lenthe, E. J. Baerends, and J. G. Snijders, J. Chem. Phys. 99, 4597 (1993). 18E. van Lenthe, E. J. Baerends, and J. G. Snijders, J. Chem. Phys. 101, 9783 (1994). 19E. van Lenthe, A. W. Ehlers, and E. J. Baerends, J. Chem. Phys. 110, 8943 (1999). 20W. Kutzelnigg and W. Liu, J. Chem. Phys. 123, 241102 (2005). 21W. Kutzelnigg and W. Liu, Mol. Phys. 104, 2225 (2006). 22W. Liu and W. Kutzelnigg, J. Chem. Phys. 126, 114107 (2007). 23W. Liu and D. Peng, J. Chem. Phys. 125, 044102 (2006). 24D. Peng, D. Liu, Y. Xiao, and L. Cheng, J. Chem. Phys. 127, 104106 (2007). 25M. Iliaš and T. Saue, J. Chem. Phys. 126, 064102 (2007). 26D. Peng and M. Reiher, Theor. Chem. Acc. 131, 1081 (2012). 27D. Peng and M. Reiher, J. Chem. Phys. 136, 244108 (2012). 28D. Peng, N. Middendorf, F. Weigend, and M. Reiher, J. Chem. Phys. 138, 184105 (2013). 29W. Liu and D. Peng, J. Chem. Phys. 131, 031104 (2009). 30W. Xu, J. Ma, D. Peng, W. Zou, W. Liu, and V. Staemmler, Chem. Phys. 356, 219 (2009). 31W. Xu, Y. Zhang, and W. Liu, Sci. China, Ser. B: Chem. 52, 1945 (2009). 32Y. Zhang, W. Xu, Q. Sun, W. Zou, and W. Liu, J. Comput. Chem. 31, 532 (2010). 33H. J. Aa. Jensen, R. Bast, T. Saue, L. Visscher, V. Bakken, K. G. Dyall, S. Dubillard, U. Ekström, E. Eliav, T. Enevoldsen, T. Fleig, O. Fossgaard, A. S. P. Gomes, T. Helgaker, J. K. Lærdahl, Y. S. Lee, J. Henriksson, M. Iliaš, Ch. R. Jacob, S. Knecht, S. Komorovský, O. Kullie, C. V. Larsen, H. S. Nataraj, P. Norman, G. Olejniczak, J. Olsen, Y. C. Park, J. K. Pedersen, M. Pernpointner, K. Ruud, P. Sałek, B. Schimmelpfennig, J. Sikkema, A. J. Thorvaldsen, J. Thyssen, J. van Stralen, S. Villaume, O. Visser, T. Winther, and S. Yamamoto, DIRAC, a relativistic ab initio electronic structure program, Release DIRAC12 (2012), see http://www.diracprogram.org. 34M. Kühn and F. Weigend, J. Chem. Theory Comput. 9, 5341 (2013). 35Local version of TURBOMOLE V6.6 2014, a development of University of Karlsruhe and Forschungszentrum Karlsruhe GmbH, 1989-2007, TURBOMOLE GmbH since 2007. 36K. A. Peterson, D. Figgen, M. Dolg, and H. Stoll, J. Chem. Phys. 126, 124101 (2007). 37D. Figgen, K. A. Peterson, M. Dolg, and H. Stoll, J. Chem. Phys. 130, 164108 (2009). 38D. Figgen, G. Rauhut, M. Dolg, and H. Stoll, Chem. Phys. 311, 227 (2005). 39B. Metz, H. Stoll, and M. Dolg, J. Chem. Phys. 113, 2563 (2000). 40K. A. Peterson, D. Figgen, E. Goll, H. Stoll, and M. Dolg, J. Chem. Phys. 119, 11113 (2003). 41I. S. Lim, P. Schwerdtfeger, B. Metz, and H. Stoll, J. Chem. Phys. 122, 104103 (2005). 42I. S. Lim, H. Stoll, and P. Schwerdtfeger, J. Chem. Phys. 124, 034107 (2006). 43D. E. Jiang, M. Kühn, Q. Tang, and F. Weigend, J. Phys. Chem. Lett. 5, 3286 (2014). 44M. Kühn and F. Weigend, J. Chem. Phys. 141, 224302 (2014). 45P. Tecmer, R. Bast, K. Ruud, and L. Visscher, J. Phys. Chem. A 116, 7397 (2012). 46M. E. Casida, in Recent Advances in Density Functional Methods, Part I, edited by D. P. Chong (World Scientific, Singapore, 1995), pp. 155–192. 47F. Wang and W. Liu, J. Chin. Chem. Soc. 50, 597 (2003). 48F. Wang and T. Ziegler, J. Chem. Phys. 121, 12191 (2004). 49R. Bauernschmitt and R. Ahlrichs, Chem. Phys. Lett. 256, 454 (1996). 50F. Furche and D. Rappoport, Computational Photochemistry, Theoretical and Computational Chemistry edited by M. Olivucci (Elsevier, Amsterdam, 2005), Vol. 16, p. 93. 9D.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-8 51M.
M. Kühn and F. Weigend
Kühn, J. Chem. Theory Comput. 10, 623 (2014). Oddershede and P. Jørgensen, J. Chem. Phys. 66, 1541 (1977). 53T. D. Bouman and A. E. Hansen, J. Chem. Phys. 66, 3460 (1977). 54M. E. Casida and M. Huix-Rotllant, Annu. Rev. Phys. Chem. 63, 287 (2012). 55M. Grüning, A. Marini, and X. Gonze, Nano Lett. 9, 2820 (2009). 56D. Peng, J. Ma, and W. Liu, Int. J. Quantum Chem. 109, 2149 (2009). 57See supplementary material at http://dx.doi.org/10.1063/1.4905829 for detailed description of the implemented 2c TDA algorithm and coordinates of all structures. 58C. van Wüllen, J. Comput. Chem. 32, 1195 (2010). 59A. D. Becke, Phys. Rev. A 38, 3098 (1988). 60J. P. Perdew, Phys. Rev. B 33, 8822 (1986). 61F. Weigend and A. Baldes, J. Chem. Phys. 133, 174102 (2010). 62K. A. Peterson, B. C. Shepler, D. Figgen, and H. Stoll, J. Phys. Chem. A 110, 13877 (2006). 63A. D. Becke, J. Chem. Phys. 98, 5648 (1993). 64F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys. 7, 3297 (2005). 65D. Andrae, U. Häußermann, M. Dolg, H. Stoll, and H. Preuß, Theor. Chim. Acta 77, 123 (1990). 66J. P. Perdew and Y. Wang, Phys. Rev. B 45, 13244 (1992). 52J.
J. Chem. Phys. 142, 034116 (2015) 67M.
Kühn and F. Weigend, ChemPhysChem 12, 3331 (2011). Eichkorn, F. Weigend, O. Treutler, and R. Ahlrichs, Theor. Chem. Acc. 97, 119 (1997). 69R. Bauernschmitt, M. Häser, O. Treutler, and R. Ahlrichs, Chem. Phys. Lett. 264, 573 (1997). 70F. Weigend, Phys. Chem. Chem. Phys. 8, 1057 (2006). 71M. K. Armbruster, W. Klopper, and F. Weigend, Phys. Chem. Chem. Phys. 8, 4862 (2006). 72O. Treutler and R. Ahlrichs, J. Chem. Phys. 102, 346 (1995). 73M. A. Baldo, S. Lamansky, P. E. Burrows, M. E. Thompson, and S. R. Forrest, Appl. Phys. Lett. 75, 4 (1999). 74W. J. Finkenzeller and H. Yersin, Chem. Phys. Lett. 377, 299 (2003). 75J. Fine, K. Diri, A. I. Krylov, C. Nemirow, Z. Lu, and C. Wittig, Mol. Phys. 110, 1849 (2012). 76C. E. Moore, Atomic Energy Levels Vol. 1-3, Natl. Bur. Stand. (U.S.) Circ. No. 467 (U.S. Government Printing Office, Washington, DC, 1949, 1952, 1958). 77J. Tellinghuisen, J. Chem. Phys. 58, 2821 (1973). 78J. Tellinghuisen, J. Chem. Phys. 76, 4736 (1982). 79K. P. Huber and G. Herzberg, Molecular Spectra and Molecular Structure (Van Nostrand and Reinhold, New York, 1979), Vol. IV. 68K.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
THE JOURNAL OF CHEMICAL PHYSICS 142, 034116 (2015)
Two-component hybrid time-dependent density functional theory within the Tamm-Dancoff approximation Michael Kühn1 and Florian Weigend1,2,a) 1
Institut für Physikalische Chemie, Karlsruher Institut für Technologie, Kaiserstraße 12, 76131 Karlsruhe, Germany 2 Institut für Nanotechnologie, Karlsruher Institut für Technologie, Postfach 3640, 76021 Karlsruhe, Germany
(Received 11 November 2014; accepted 29 December 2014; published online 20 January 2015) We report the implementation of a two-component variant of time-dependent density functional theory (TDDFT) for hybrid functionals that accounts for spin-orbit effects within the Tamm-Dancoff approximation (TDA) for closed-shell systems. The influence of the admixture of Hartree-Fock exchange on excitation energies is investigated for several atoms and diatomic molecules by comparison to numbers for pure density functionals obtained previously [M. Kühn and F. Weigend, J. Chem. Theory Comput. 9, 5341 (2013)]. It is further related to changes upon switching to the local density approximation or using the full TDDFT formalism instead of TDA. Efficiency is demonstrated for a comparably large system, Ir(ppy)3 (61 atoms, 1501 basis functions, lowest 10 excited states), which is a prototype molecule for organic light-emitting diodes, due to its “spin-forbidden” triplet-singlet transition. C 2015 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4905829]
I. INTRODUCTION
Time-dependent density functional theory (TDDFT)1,2 has become a popular method for the calculation of electronic excitations of large systems. For compounds containing heavy elements, relativistic effects on spectra become important and need to be incorporated. Mostly this is limited to scalar relativistic effects since their incorporation does not change the structure of common quantum-chemical programs. However, in one-component scalar relativistic (1c) TDDFT, the influence of spin-orbit coupling (SOC) is neglected. Therefore, it is not possible to properly describe phenomena based on “spinforbidden” transitions, e.g., phosphorescence, which is the mechanism of interest in molecules applied in organic lightemitting diodes (OLEDs).3,4 When going beyond 1c TDDFT and including SOC, the most fundamental way is fully relativistic TDDFT based on the Dirac-Coulomb Hamiltonian, yielding four-component (4c) complex spinors instead of 1c real molecular orbitals (MOs).5–7 A more economic alternative of very similar accuracy for chemical problems is quasirelativistic TDDFT using twocomponent (2c) complex spinors.8,9 In the past few years, several groups reported 2c (Kramers-restricted) closed-shell TDDFT implementations using a noncollinear exchange-correlation (XC) kernel:8,9 Wang et al.8 in the Amsterdam density functional program package (ADF),10 Peng et al.9 in the Beijing density functional program suite (BDF),11–13 and Nakata et al.14 in the GAMESS (US) program package15 within the Tamm-Dancoff approximation (TDA).16 All three implementations are based on the all-electron 2c zeroth-order regular approximation (ZORA).17–19 Additionally, the all-electron exact 2c (X2C) a)Author to whom correspondence should be addressed. Electronic mail:
[email protected]. Fax: +49 721 608-47225.
0021-9606/2015/142(3)/034116/8/$30.00
approach20–29 has been applied within the implementation in the BDF program.30–32 Furthermore, Bast et al.7 claim that their 4c TDDFT implementation in the DIRAC program package33 based on the Dirac-Coulomb Hamiltonian is also fully operational for the X2C approach. Recently, we reported an implementation34 in the TURBOMOLE program suite,35 that can be used in combination with the X2C approach as well as 2c Dirac-HartreeFock effective core potentials (dhf-ECPs), which are available for elements beyond Kr36–42 and represent the most efficient choice to include relativistic effects, since only valence electrons have to be treated explicitly within the calculations, in contrast to all-electron methods. We demonstrated the high efficiency of our code by calculating the electronic spectrum of Au20 (2620 basis functions, polarized quadruple-ζ valence basis sets, 2c dhf-ECPs, lowest 57 excited states)34 and Au25(SCH3)−18 (1695 basis functions, polarized double-ζ valence basis sets, 2c dhf-ECPs, lowest 192 excited states),43 which are the largest systems studied so far at the 2c TDDFT level. Our implementation also facilitated a detailed investigation of the “spin-forbidden,” i.e., phosphorescent, transitions in an eight-membered set of iridium-containing organometallic complexes that are potential candidates for OLEDs.44 Since the emitting state in such compounds usually is of chargetransfer character (to be more precise: spin-orbit coupled triplet metal-to-ligand charge-transfer character, 3MLCT), TDDFT using pure density XC functionals largely underestimates the experimental phosphorescence energies. This is a widely known problem of TDDFT, which can—at least partly—be circumvented by inclusion of Hartree-Fock (HF) exchange. To the best of our knowledge, so far, only two 2c (and 4c) TDDFT/TDA implementations have been reported, that are able to treat HF exchange and therefore the important class of hybrid XC functionals using proper kernels (i.e., full derivatives of XC functionals and inclusion of HF exchange
142, 034116-1
© 2015 AIP Publishing LLC
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-2
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
in the coupling matrix): the implementations by Bast et al. in the DIRAC program package7 and Nakata et al. in the GAMESS (US) program suite,14 which both additionally work for the class of long-range corrected hybrid functionals. However, the systems investigated using the aforementioned implementations are of rather moderate size. To the best of our knowledge, the largest system studied so far is [UO2Cl4]2− (allelectron polarized triple-ζ valence basis sets, lowest 15 excited states),45 using the implementation of Bast et al.7 This indicates that the calculation of systems with up to about 150 atoms (a few thousand basis functions), e.g., the OLED molecules studied in Ref. 44, using the aforementioned 2c (and 4c) hybrid TDDFT/TDA codes is a challenge. The goal of this work is to extend our recently implemented efficient 2c (Kramers-restricted) closed-shell TDDFT code34 in the TURBOMOLE program suite to hybrid XC functionals using proper kernels, which is most conveniently and efficiently done within the 2c TDA making optimal use of already existing routines. Our code allows for treatments of systems with a few thousand basis functions on a single processor as well as reasonable scaling for moderately parallelized tasks. The paper is organized as follows. In Sec. II, we summarize the theory of 2c (Kramers-restricted) closed-shell TDDFT and TDA. In Sec. III, the accuracy of our 2c hybrid TDA implementation is assessed by comparing the excitation energies of the s → p transitions in Cd and Hg to the 2c TDA variant of Nakata et al.14 and the 4c TDDFT approach of Bast et al.7 as well as to experimental data. The excitation energies of the d → s transitions in Au+ and the lowest transitions in I2 and TlH are compared to experimental results. Finally, we demonstrate the efficiency of our code by calculating the phosphorescence energies of the popular OLED model compound Ir(ppy)3 using 1501 basis functions (polarized triple-ζ valence basis sets, 2c dhf-ECP, lowest 10 excited states), which is the largest system calculated so far using a 2c or 4c hybrid TDDFT/TDA algorithm.
The complex (Hermitian) orbital rotation Hessians are defined as Ai σ˜ aτ˜ j σ˜ ′ bτ˜ ′ = (ϵ aτ˜ − ϵ i σ˜ )δ i j δ ab δσ˜ σ˜ ′δτ˜ τ˜ ′ +Ci σ˜ aτ˜ j σ˜ ′ bτ˜ ′, Bi σ˜ aτ˜ j σ˜ ′ bτ˜ ′ = Ci σ˜ aτ˜ bτ˜ ′ j σ˜ ′,
(2) (3)
with the coupling matrix Ci σ˜ aτ˜ j σ˜ ′ bτ˜ ′ = ⟨Φi σ˜ Φbτ˜ ′|Φaτ˜ Φ j σ˜ ′⟩ + ⟨Φi σ˜ Φbτ˜ ′|fXC|Φaτ˜ Φ j σ˜ ′⟩ − cX⟨Φi σ˜ Φbτ˜ ′|Φ j σ˜ ′Φaτ˜ ⟩.
(4)
δ pq is the Kronecker delta and ϵ p σ˜ describes the energy of 2c spinor Φ p σ˜ (⃗r ) =
NBF α *c µβ p σ˜ + ξ µ (⃗r ) c ˜µ=1 , µ p σ
(5)
that is not necessarily an eigenfunction of Sˆ z = 21 σ z . ξ µ represent real basis functions (µ = 1,2, ..., NBF) and cαµ p σ˜ and c µβ p σ˜ are complex expansion coefficients. The first term in Eq. (4) is the Coulomb contribution in Dirac notation. The second term contains the matrix elements of the noncollinear XC kernel.5,6,8,9,47,48 The last term is the HF exchange contribution scaled by a hybrid mixing parameter cX. Note that there is no additional SOC contribution to the coupling matrix. SOC effects on excitations arise via the reference-state spinors and their energies only, as only one-electron relativistic potentials are included in the Hamiltonian (via ECPs, X2C, ZORA, . . . ). In the following, we only consider systems possessing a (Kramers-restricted) closed-shell ground state, for which Nα˜ = Nβ˜ = N2 and ϵ p α˜ = ϵ p β˜ since the Kramers partners Φ p α˜ and Φ p β˜ are either both occupied or unoccupied. The noncollinear XC kernel within the adiabatic approximation for (Kramersrestricted) closed-shell systems reads8,9 fXC(⃗r ,⃗r ′) =
δ2 EXC 1⊗1 δ ρ(⃗r )δ ρ(⃗r ′) δ2 EXC + δ(⃗r −⃗r ′)σ k ⊗ σ k . 2(⃗ δs r ) k=x, y, z
(6)
⃗ (⃗r )| is the noncollinear spin density with the magnes(⃗r ) = | m tization vector
II. TWO-COMPONENT TIME-DEPENDENT DENSITY FUNCTIONAL THEORY
The 2c non-Hermitian time-dependent Kohn-Sham (TDKS) eigenvalue problem that accounts for SOC effects on excitation energies ω n (n = 1,2,3, ...)2,7–9,34,46 reads
⃗ (⃗r ) = m
Nσ˜ σ ˜ i=1
Φ†i σ˜ (⃗r )⃗ σ Φi σ˜ (⃗r ).
(7)
The total density is given as ′ ′ * Ai σ˜ aτ˜ j σ˜ bτ˜ ∗ ′ ′ , Bi σ˜ aτ˜ j σ˜ bτ˜
1 = ωn * ,0
Bi σ˜ aτ˜ j σ˜ ′ bτ˜ ′ + * X j σ˜ ′ bτ˜ ′, n + ∗ Ai σ˜ aτ˜ j σ˜ ′ bτ˜ ′ - , Yj σ˜ ′ bτ˜ ′, n -
0 + * X i σ˜ aτ,˜ n + , −1- , Yi σ˜ aτ,˜ n -
ρ(⃗r ) =
Nσ˜ σ ˜ i=1
(1)
where the Einstein summation convention is utilized and indices i, j, ... are used for occupied, a, b, ... for virtual, and p, q, ... for general molecular spinors. Greek indices σ, ˜ τ, ˜ ˜ refer to 2c time-reversal symmetry-adapted σ ˜ ′, τ˜ ′ ∈ { α, ˜ β} Kramers partners describing Nα˜ moment-up (α) ˜ and Nβ˜ ˜ electrons. It is noted that Eq. (1) also holds moment-down ( β) for 4c TDDFT.5–7,9
Φ†i σ˜ (⃗r )Φi σ˜ (⃗r ).
(8)
In Eq. (6), EXC is the nonrelativistic XC energy depending on both ρ and s. The left side of the tensor product refers to electron 1 (indices i and a in Eq. (4)), whereas the right side refers to electron 2 (indices b and j). To arrive at a convenient formulation of 2c TDDFT, that is able to make optimal use of the existing 1c TDDFT code49,50 in the TURBOMOLE program suite,35 the goal is to reduce the two-superdimensional non-Hermitian eigenvalue problem in Eq. (1) to a one-superdimensional Hermitian eigenvalue problem. In the following, we will discuss two possibilities: 2c
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-3
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
TDDFT within Ziegler’s approximation8,9,34 in Sec. II A and the 2c TDA in Sec. II B.16 A. Two-component time-dependent density functional theory within Ziegler’s approximation
In the case of pure density functionals, i.e., cX = 0 in the coupling matrix (Eq. (4)), an approximate 2c TDKS Hermitian eigenvalue problem for the square of the excitation energies of (Kramers-restricted) closed-shell systems was derived by Wang, Ziegler, van Lenthe, van Gisbergen, and Baerends8 as well as Peng, Zou, and Liu.9 It reads Ωi σ˜ aτ˜ j σ˜ ′ bτ˜ ′ (X +Y )′j σ˜ ′ bτ˜ ′, n = ω2n (X +Y )i′σ˜ aτ,˜ n ,
(9)
Ωi σ˜ aτ˜ j σ˜ ′ bτ˜ ′ = (ϵ aτ˜ − ϵ i σ˜ )2δ i j δ ab δσ˜ σ˜ ′δτ˜ τ˜ ′ √ + 2 ϵ aτ˜ − ϵ i σ˜ Ci σ˜ aτ˜ j σ˜ ′ bτ˜ ′ ϵ bτ˜ ′ − ϵ j σ˜ ′.
(10)
with
Note that the transformation from Eq. (1) to Eq. (9) is only approximate if SOC is present since the (polarized) XC kernel (Eq. (6)) breaks time-reversal (or Kramers) symmetry, which has to be exploited in the transformation. However, Wang et al. claim that test calculations show errors always less than 0.05 eV compared to a calculation according to Eq. (1).8 Peng et al. state that this approximation may not work if either the excitation energy is not close to the orbital energy difference or if there are highly inhomogeneous induced fields (e.g., localized non-fully filled f spinors of a formally closedshell system).9 In case the XC contribution to the coupling matrix in Eq. (4) is neglected, i.e., only the Coulomb contribution is included, it has been shown that Eq. (9) is obtained by an exact transformation also in the presence of SOC.51 In all cases, time-reversal symmetry has to be exploited so that the transformation is only applicable when (Kramers-restricted) closed-shell systems are considered. For open-shell systems, it is not possible to apply an analogous transformation that arrives at Eq. (9). At the 1c level in contrast, the transformation to the 1c analog of Eq. (9) is exact for all cases, as long as HF exchange is excluded, i.e., cX = 0. Recently, we reported the implementation of the formalism presented in this section in the ESCF excited state module49,50 of the TURBOMOLE program suite,35 see Ref. 34 for details.
of TDA is that transition moments (e.g., the transition dipole moment) and therefore oscillator strengths and radiative lifetimes are ill-defined due to the lack of gauge invariance.50 In practice, this often leads to poor results for TDA oscillator strengths.2,52–55 It is noted in passing that in contrast to Eq. (9), Eq. (11) is also applicable when open-shell systems are considered, however, the expression for the noncollinear XC kernel (Eq. (6)) becomes quite involved.9,48 Both 2c TDA (Eq. (11)) and 2c TDDFT within Ziegler’s approximation (Eq. (9)), that are based on a noncollinear XC kernel, have the correct 1c limit, where moment-up/moment˜ electrons become ordinary spin-up/spin-down down (α/ ˜ β) (α/ β) electrons. As a result, the coupling matrix is real, the MOs Φ pσ are eigenfunctions of Sˆ z , and the XC kernel is collinear. Besides this correct 1c limit, all of the above equations can be easily extended to the fully relativistic 4c formalism based on the Dirac-Coulomb Hamiltonian, simply by replacing all 2c quantities by their 4c analogs.5–7,9 2c TDA (and also 2c TDDFT within Ziegler’s approximation) is computationally more demanding than 1c TDA. For the 2c (Kramers-restricted) closed-shell treatment, we expect an additional factor of eight in both computation time and required memory compared to the 1c closed-shell formalism, due to the larger size of the excitation vectors X i σ˜ aτ,˜ n (n = const): a factor of four arises since the number of both occupied and virtual molecular orbitals/spinors is doubled and an additional factor of two is due to the introduction of imaginary parts (complex spinors). It is noted that by exploitation of double group symmetry together with time-reversal symmetry, the computational cost for the 2c (and also 4c) treatment can be reduced,56 as done in Refs. 7, 8, and 30–32. The 2c TDA presented in this section was implemented in the ESCF excited-state module34,49,50 of the TURBOMOLE program suite,35 mainly by extending the existing 1c routines to the 2c formalism and making use of our recently implemented 2c TDDFT code.34 Details of the implementation are given in Figure 1 of the supplementary material.57 Our code is parallelized employing the technique by van Wüllen58 in the same way as for the 1c case.
III. APPLICATION A. Computational details
B. Two-component Tamm-Dancoff approximation
To arrive at a 2c excited state theory that is able to treat HF exchange (and therefore hybrid functionals), i.e., cX , 0 in the coupling matrix (Eq. (4)), most conveniently the TDA,16 that neglects Y in Eq. (1), is used. The 2c Hermitian eigenvalue problem for excitation energies of (Kramers-restricted) closedshell systems within the TDA reads14 Ai σ˜ aτ˜ j σ˜ ′ bτ˜ ′ X j σ˜ ′ bτ˜ ′, n = ωTDA X i σ˜ aτ,˜ n . n
(11)
In the case of cX = 1 and no XC contribution to the coupling matrix (Eq. (4)), 2c TDA is identical to the 2c configuration interaction singles (CIS) method. Compared to the full twosuperdimensional formalism in Eq. (1), usage of TDA leads to systematically higher excitation energies. A disadvantage
All ground-state closed-shell structures of diatomic molecules were optimized using 1c DFT with the generalized gradient approximation (GGA) functional of Becke59 and Perdew60 (BP86) in combination with polarized quadruple-ζ valence basis sets (dhf-QZVP-2c).61 The dhf-bases are optimized for the usage in connection with dhf-ECPs, which account for scalar relativistic effects in 1c calculations and additionally for SOC effects in 2c calculations. For I a pseudopotential (dhf-ECP-28)62 covering the inner 28 electrons was used, for Tl a pseudopotential (dhf-ECP-60)39 covering the inner 60 electrons. The structure of the lowest excited triplet state of Ir(ppy)3 was optimized without symmetry constraints using 1c ground-state DFT in combination with Becke’s threeparameter hybrid functional with Lee-Yang-Parr correlation
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-4
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
(B3LYP)63 and polarized triple-ζ valence basis sets (def2TZVP).64 See also Ref. 44 for further information on the choice of the level of theory for the structure optimization. For Ir, a Wood-Boring-type pseudopotential (mwb-ECP-60)65 was used, which covers the inner 60 electrons and accounts for scalar relativistic effects only. Coordinates of all optimized structures are available in Tables I-III of the supplementary material.57 Vertical excitation energies (VEEs) were computed starting from a (Kramers-restricted) singlet ground state at both the 1c and 2c TDA16 as well as the TDDFT level34,49,50, using the Perdew-Wang parameterization of the local spin-density approximation (LSDA),66 the GGA functional BP86,59,60 and the hybrid functional B3LYP.63 Additionally, for Ir(ppy)3 at the 1c level, phosphorescence energies were calculated directly by vertical de-excitations from the lowest excited triplet state (T1) to the singlet ground state (S0) using the spin-flip TammDancoff approximation (SF-TDA).67 For atoms and diatomic molecules, dhf-QZVP-2c basis sets were used, for Ir(ppy)3 basis sets of polarized triple-ζ valence quality (dhf-TZVP2c).61 For Cd, a pseudopotential (dhf-ECP-28)38 covering the inner 28 electrons, for Au, Hg, and Ir pseudopotentials (dhfECP-60)37,38 covering the inner 60 electrons were used. The iterative Davidson procedure in the TDDFT/TDA calculations is stopped when the Euclidean norm of the residual vector is smaller than 10−5 a.u. In all calculations using pure density functionals, the resolution of the identity (RI) approximation68,69 was employed. The applied auxiliary basis sets70 are suited for bases optimized for dhf-ECPs.71 Reference state energies and density matrices were converged to 10−7 a.u. Fine quadrature grids of size m4 (Ref. 72) were employed. 1. Assessment of accuracy: Excitation energies and fine-structure splittings of selected closed-shell atoms, ions, and diatomic molecules
Accuracy of the method was studied for the s → p transitions in the atoms Cd and Hg, the d → s transitions in Au+ as well as for the lowest transitions in I2 and TlH. For all cases, we consider the influence of the admixture of HF exchange by
comparison with the (pure) GGA functional BP86 and relate it to both the effect of the gradient correction by comparison with the LDA functional LSDA and to the influence of the TDA by comparison to the full TDDFT formalism (for the LSDA functional). Further, results are compared to other implementations of hybrid TDDFT7,14 and to the experiment, as far as available. First, we investigate the s → p transitions in the Cd atom as well as the Hg atom in detail, see Table I. The quantities of interest are the average vertical excitation energies for the triplet excitation (3 P), the respective splittings due to SOC, and the energy difference to the singlet (1 P1) excitation (“shift”). For both atoms and all quantities, the differences between B3LYP and BP86 are very small. The maximum difference is observed for the shift of the 1 P1 transition in case of Cd, ca. 0.1 eV. Energies for the triplet splittings are essentially the same for these two as well as for the other methods. For LSDA, comparably large differences to the BP86 treatment for the average 3 P VEE are observed, typically 0.3 eV. Employing the full TDDFT formalism instead of the TDA variant (at LSDA level) leads to VEEs that are slightly smaller (by less than 0.1 eV for the average 3 P VEE and by somewhat more than 0.1 eV for the 1 P1 VEE). Further, 2c (average) VEEs and 1c VEEs (given in parentheses in Table I) are very similar for all cases with maximum differences amounting to 0.06 eV. Our B3LYP results very well agree with those obtained within a 2c ZORA approach, Ref. 14 and at the 4c Dirac-Coulomb level, Ref. 7. In case of Cd, our results are more similar to those of Ref. 14; in case of Hg, they are always between those from the two references. This is not unexpected as in Ref. 7 a 4c (Dirac) code is used, in Ref. 14 a 2c ZORA approach. Our method, a 2c approach with ECPs fitted to 4c Dirac results, could be regarded as a mixture of both. Also, agreement with experimental data is very good, deviations for splittings are below 0.03 eV, for excitation energies up to ca. 0.15 eV, but it has to be noted that this is not better than the numbers obtained with the significantly more economic BP86 functional. The d → s excitation in Au+ reveals a similar picture concerning the sublevel splitting of the 3 D transition, see Table II; almost identical numbers are obtained with the four types of calculations that moreover agree well with the experiment. A marked difference to Cd or Hg is observed for the
TABLE I. Two-component (2c) VEEs for s → p transitions, fine-structure splittings, and shifts (given relatively to average of triplet state) of Cd and Hg. Unless noted otherwise the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with dhf-QZVP-2c basis sets and dhf-ecp-2c spin-orbit effective core potentials (ECPs). Values in parentheses are scalar relativistic (1c) results for the 3 P and 1 P states. Degeneracies g J = (2J + 1) are also given. All values are in eV. TDDFT results are from Ref. 34. 2c TDA VEEs using ZORA obtained by Nakata et al.14 are listed in column “B3LYP,14” 4c TDDFT VEEs based on the Dirac-Coulomb Hamiltonian by Bast et al.7 in column “B3LYP.7” Experimental results are from Ref. 76. gJ
B3LYP
BP86
LSDA
LSDA (TDDFT)
B3LYP14
B3LYP7
Expt.
Average of triplet state Splitting of 3 P0 Splitting of 3 P1 Splitting of 3 P2 Shift of 1 P1
1 3 5 3
3.921 (3.924) −0.134 −0.069 0.068 1.514 (1.513)
3.909 (3.912) −0.136 −0.070 0.069 1.620 (1.619)
4.236 (4.240) −0.139 −0.072 0.071 1.425 (1.424)
4.195 (4.199) −0.141 −0.074 0.072 1.271 (1.269)
3.91 −0.14 −0.07 0.07 1.63
3.81 −0.15 −0.08 0.07 1.33
3.874 −0.140 −0.073 0.072 1.543
Average of triplet state Splitting of 3 P0 Splitting of 3 P1 Splitting of 3 P2 Shift of 1 P1
1 3 5 3
5.073 (5.119) −0.482 −0.268 0.257 1.358 (1.341)
5.151 (5.201) −0.483 −0.266 0.256 1.436 (1.422)
5.452 (5.506) −0.491 −0.277 0.264 1.251 (1.236)
5.411 (5.470) −0.492 −0.287 0.270 1.156 (1.131)
5.17 −0.41 −0.23 0.22 1.45
5.03 −0.52 −0.30 0.28 1.24
5.181 −0.514 −0.295 0.280 1.523
Atom Cd
Hg
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-5
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
TABLE II. Two-component (2c) VEEs for d → s transitions, fine-structure splittings, and shifts (given relatively to average of triplet state) of Au+. Unless noted otherwise, the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with the dhf-QZVP-2c basis set and dhf-ecp-2c spin-orbit effective core potential (ECP). Values in parentheses are scalar relativistic (1c) results for the 3 D and 1 D states. Degeneracies g J = (2J + 1) are also given. All values are in eV. Experimental results are from Ref. 76.
Average of triplet state Splitting of 3 D3 Splitting of 3 D2 Splitting of 3 D1 Shift of 1 D 2
gJ
B3LYP
BP86
LSDA
LSDA (TDDFT)
Expt.
7 5 3 5
1.766 (1.903) −0.451 −0.015 1.077 1.581 (1.078)
1.546 (1.677) −0.459 −0.007 1.083 1.621 (1.120)
1.554 (1.687) −0.442 −0.020 1.065 1.550 (1.039)
1.498 (1.642) −0.454 −0.020 1.092 1.531 (1.002)
2.288 −0.423 −0.101 1.155 1.385
average triplet excitation energy: for both Cd and Hg, values for B3LYP and BP86 were almost identical and LSDA significantly differed from those, but for Au+ BP86 and LSDA are very similar and B3LYP is larger by ca. 0.2 eV, thus being closer to the experimental value; nevertheless, the difference to the experiment, ca. 0.5 eV, is significant. This error can be traced back to the different self-interaction experienced by the 5 d and the 6 s levels, as noted in the work of Wang et al.8 The 5 d levels are more compact than the 6 s level and therefore are more prone to unphysical self-interaction. As a result, the energies of the former are overestimated compared to the energy of the latter, consequently leading to an underestimation of the 5 d → 6 s excitation energies, which in zeroth order are given by the differences of the respective orbital/spinor energies. Concerning functional changes, see above, the same trends are obtained for the calculation at the 1c level, but here the 3 D VEEs are larger by ca. 0.15 eV throughout. This difference may be explained by the splitting of the fully occupied 5 d orbitals by SOC, which is not present in the 1c calculation. The 3 D excitation in the 2c case mainly involves the d 5/2 level which is higher in energy than the d shell in the 1c case.
For the diatomic molecules I2 and TlH, see Tables III and IV, B3LYP energies for excitations to the Π states (or their 2c analogs, respectively) are slightly (ca. 0.1–0.2 eV) higher than those obtained with BP86; this is somewhat more than changes from BP86 to LSDA/TDA or from the latter to LSDA/TDDFT (both below 0.1 eV). Excitations to the Σ states of I2 are very similar for B3LYP and BP86, here larger differences (0.2 eV) are observed between the functionals BP86 and LSDA. This holds for both 1c and 2c treatments. Differences to the few available experimental data are smallest for B3LYP in all cases, excitations to the Π states are underestimated by up to 0.4 eV, the experimental excitation energy to the 1u (3Σu+) state of I2 is almost exactly met (−0.03 eV). Overall, for the systems discussed so far, the gain from the comparable costly (see below) B3LYP is rather moderate; the accuracy of level splittings is very satisfying also for the cheaper non-hybrid (pure) density functionals, the accuracy of averaged 2c excitation energies (and 1c excitation energies) is systematically improved by B3LYP only for the transitions in Au+ as well as the Π states of the diatomic molecules.
TABLE III. Scalar relativistic (1c) (states expressed in the Λ − Σ notation) as well as two-component (2c) VEEs (states expressed in the ω − ω notation) of I2. Unless noted otherwise the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with the dhf-QZVP-2c basis set and dhf-ecp-2c spin-orbit effective core potential (ECP). Degeneracies are also given. All values are in eV. TDDFT results are from Ref. 34. Excited state Scalar relativistic (1c) 3Π u 1Π u 3Π g 1Π g 3Σ + u Two-component (2c) 2u 1u 0−u 0+u 1u 2g 1g 0−g 0+g 1g 0−u 1u
Composition/contribution
Degeneracy
0+u , 0−u , 1u , 2u 1u 0+g , 0−g , 1g , 2g 1g 0−u , 1u 3Π
u
3Π
u
3Π
u
3Π
u
1Π
u 3Π g 3Π g 3Π g 3Π g 1Π g 3Σ + u 3Σ + u
B3LYP
BP86
LSDA
LSDA (TDDFT)
6 2 6 2 3
1.60 2.21 3.21 3.80 4.35
1.49 2.11 3.03 3.63 4.39
1.52 2.08 3.08 3.61 4.59
1.50 2.00 3.07 3.56 4.55
2 2 1 1 2 2 2 1 1 2 1 2
1.32 1.49 1.90 2.03 2.30 2.97 3.14 3.63 3.65 4.03 4.51 4.54
1.22 1.39 1.81 1.93 2.22 2.80 2.97 3.44 3.46 3.85 4.53 4.56
1.25 1.41 1.86 1.97 2.21 2.86 3.01 3.51 3.53 3.86 4.73 4.75
1.22 1.38 1.84 1.94 2.14 2.84 2.99 3.49 3.51 3.81 4.69 4.71
Expt.
1.6977 1.8477 2.1377 2.3778 2.4977
4.5779
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-6
M. Kühn and F. Weigend
J. Chem. Phys. 142, 034116 (2015)
TABLE IV. Scalar relativistic (1c) (states expressed in the Λ − Σ notation) as well as two-component (2c) VEEs (states expressed in the ω − ω notation) of TlH. Unless noted otherwise, the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with dhf-QZVP-2c basis sets and the dhf-ecp-2c spin-orbit effective core potential (ECP) for Tl. Degeneracies are also given. All values are in eV. TDDFT results are from Ref. 34. Experimental results are from Ref. 79. Excited state Scalar relativistic (1c) 3Π 1Π Two-component (2c) 0− 1 0+ 1 2
Composition/contribution
Degeneracy
B3LYP
BP86
LSDA
LSDA (TDDFT)
0+, 0−, 1, 2 1
6 2
2.24 2.75
2.17 2.69
2.15 2.59
2.14 2.54
1 2 1 2 2
1.97 2.15 2.17 2.69 2.99
1.91 2.08 2.10 2.63 2.95
1.92 2.08 2.14 2.65 2.90
1.89 2.05 2.09 2.64 2.86
3Π 3Π 3Π 1Π 3Π
2. Demonstration of efficiency: Phosphorescence energies of Ir(ppy)3
Ir(ppy)3, see Figure 1, is a prototype molecule for OLEDs; the lowest transition of this molecule is “spin-forbidden” with pronounced charge-transfer (CT) character. Thorough studies of this system with pure density functionals were presented recently.44 Calculations for this molecule at the B3LYP level were carried out in the present work both for demonstration of efficiency and in order to investigate the effect of admixing HF exchange. For this transition, B3LYP significantly improves the accuracy of the average triplet excitation energy (phospho-
FIG. 1. Structure of Ir(ppy)3 (fac-tris[2-(2-pyridinyl-κ N )phenyl-κC]iridium) in Lewis notation.
Expt.
2.20 3.02 3.00
rescence energy), see Table V. With B3LYP, the difference to the experimental value of 2.43 eV (Refs. 73 and 74) amounts to 0.3 eV, for BP86, LSDA/TDA and LSDA/TDDFT more than 0.6 eV are obtained. This improvement is expected, as HF copes much better with CT states than DFT, thus leading to significantly higher excitation energies. Consequently, the same trend is observed for the 1c treatments, which show excitation energies higher by ca. 0.1 eV throughout, thus being even closer to the experiment. It is noted that B3LYP despite being superior to LDA and GGA functionals still has limitations concerning the calculation of CT transitions. It is commonly believed that within TDDFT, the latter are described most accurately by long-range corrected hybrid functionals. Using those, the 1c triplet excitation energy for Ir(ppy)3 was calculated to 2.3 eV, which is slightly closer to the experimental value than the B3LYP result.75 A gain in accuracy similar to that of B3LYP can be obtained in a more economic way, namely, by correcting the LSDA numbers with excitation energies obtained from SF-TDA,48,67 as described previously,44,67 see second column from the right in Table V. The size of the triplet splitting is overestimated by all methods. The difference between the highest and the lowest sublevel amounts to 0.010 eV for the experiment, to 0.024 eV for B3LYP, and to 0.036 eV for BP86. LSDA yields 0.032 and LSDA/TDDFT 0.028 eV. Apparently, B3LYP is closest to the experiment, but it has to be noted that the triplet sublevel splitting is very sensitive to the structure parameters of the molecule, as demonstrated previously,44 and thus a reliable comparison is difficult. The increase of computational costs compared to, e.g., BP86 for the improvement of B3LYP is comparably large,
TABLE V. Two-component (2c) phosphorescence energies and fine-structure splittings of triplet substates TSO 1, m (m = 1, 2, 3) (given relatively to average of triplet state) of Ir(ppy)3. Unless noted otherwise, the TDA was used. All calculated values were obtained using the B3LYP, BP86, and LSDA functionals in combination with dhf-TZVP-2c basis sets and the dhf-ecp-2c spin-orbit effective core potential (ECP) for Ir. Values in parentheses are scalar relativistic (1c) results. All values are in eV. TDDFT results are from Ref. 34, SF-TDA results from Ref. 44. In case of the latter, the 1c phosphorescence energies were calculated directly as de-excitations from T1 to S0 using the 1c SF-TDA. The absolute energies of the triplet substates TSO 1, m (obtained with 2c TDDFT) were shifted by the difference between the 1c SF-TDA and the 1c TDDFT results. Experimental results are from Refs. 73 and 74.
Average of triplet state Splitting of TSO 1,1 Splitting of TSO 1,2 Splitting of TSO 1,3
B3LYP
BP86
LSDA
LSDA (TDDFT)
LSDA (SF-TDA)
Expt.
2.1488 (2.2114) −0.0090 −0.0064 0.0154
1.7993 (1.8842) −0.0137 −0.0084 0.0221
1.7831 (1.8675) −0.0125 −0.0073 0.0198
1.7770 (1.8618) −0.0110 −0.0063 0.0173
2.1858 (2.2706) −0.0110 −0.0063 0.0173
2.4311 −0.0040 −0.0023 0.0063
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-7
M. Kühn and F. Weigend
at least in case of TURBOMOLE, where non-hybrid DFT is implemented in a very efficient way by making use of efficient numerical grids as well as the RI approximation for the Coulomb part. This makes HF exchange by far the dominant part in hybrid DFT treatments, for 1c and for 2c calculations and for the ground state as well as for excited states. The CPU time for the determination of the lowest ten excitations at the 2c B3LYP level using 1501 basis functions amounts to 30 days and 16 h on a single Intel Xeon X5650 2.67 GHz CPU; for convergence of the excitation vectors, 18 iterations were needed. This effort is much higher than at the 2c BP86 level, 18 h and 25 min for 11 iterations. Thus, the time per iteration is higher for B3LYP than for BP86 by about a factor of 24.4. Almost the same factor between B3LYP and BP86 for the CPU times per iteration is obtained for 1c calculations, 23.2. Here, nine iterations are needed for BP86 and 14 for B3LYP. The ratio between CPU times per iteration for 2c and 1c treatments amounts to a factor of 8.2 for B3LYP and to 7.8 for BP86, which is close to the expected value of eight, see Sec. II B.
IV. CONCLUSION
We presented the implementation of a two-component variant of hybrid TDDFT. It allows for routine treatments of spin-orbit effects in the calculation of excitations of closedshell molecules with about 50 atoms. Concerning efficiency and treatable system size, this is a remarkable advance over already existing codes.7,14 In contrast to the latter, twocomponent effective core potentials are used here, which at least in part is the reason for the higher efficiency. Results for excitation energies and sublevel splittings of atoms and small molecules obtained with this approach usually are between those published for the four-component hybrid TDDFT based on the Dirac-Coulomb Hamiltonian7 and the two-component hybrid TDDFT within the ZORA approach.14 Significant improvement over GGA or LDA functionals is observed mainly for the average values of excitation energies in case of charge-transfer transitions, which comes along with an increase of CPU times by a factor of about 25, very similar to one-component treatments. Excitation energies for other cases and sublevel splittings in general are similar for B3LYP and the other functionals. ACKNOWLEDGMENTS
M.K. is funded by the Carl Zeiss Foundation. The author thanks TURBOMOLE GmbH for financial support. 1Time-Dependent
Density Functional Theory, Lecture Notes in Physics, edited by M. A. L. Marques, C. A. Ullrich, F. Nogueira, A. Rubio, K. Burke, and E. K. U. Gross (Springer, Heidelberg, 2006). 2F. Furche, J. Chem. Phys. 114, 5982 (2001). 3M. A. Baldo, D. F. O’Brien, Y. You, A. Shoustikov, S. Sibley, M. E. Thompson, and S. R. Forrest, Nature 395, 151 (1998). 4R. C. Evans, P. Douglas, and C. J. Winscom, Coord. Chem. Rev. 250, 2093 (2006). 5J. Gao, W. Liu, B. Song, and C. Liu, J. Chem. Phys. 121, 6658 (2004). 6J. Gao, W. Zou, W. Liu, Y. Xiao, D. Peng, B. Song, and C. Liu, J. Chem. Phys. 123, 054102 (2005). 7R. Bast, H. J. A. Jensen, and T. Saue, Int. J. Quantum Chem. 109, 2091 (2009).
J. Chem. Phys. 142, 034116 (2015) 8F.
Wang, T. Ziegler, E. van Lenthe, S. van Gisbergen, and E. J. Baerends, J. Chem. Phys. 122, 204103 (2005). Peng, W. Zou, and W. Liu, J. Chem. Phys. 123, 144101 (2005). 10See http://www.scm.com for Amsterdam Density Functional (ADF) version 2013.01a. 11W. Liu, G. Hong, D. Dai, L. Li, and M. Dolg, Theor. Chem. Acc. 96, 75 (1997). 12W. Liu, F. Wang, and L. Li, J. Theor. Comput. Chem. 2, 257 (2003). 13W. Liu, F. Wang, and L. Li, Recent Advances in Relativistic Molecular Theory, Recent Advances in Computational Chemistry edited by K. Hirao and Y. Ishikawa (World Scientific, Singapore, 2004), Vol. 5, p. 257. 14A. Nakata, T. Tsuneda, and K. Hirao, J. Chem. Phys. 135, 224106 (2011). 15M. W. Schmidt, K. K. Baldridge, J. A. Boatz, S. T. Elbert, M. S. Gordon, J. H. Jensen, S. Koseki, N. Matsunaga, K. A. Nguyen, S. J. Su, T. L. Windus, M. Dupuis, and J. A. Montgomery, J. Comput. Chem. 14, 1347 (1993). 16S. Hirata and M. Head-Gordon, Chem. Phys. Lett. 314, 291 (1999). 17E. van Lenthe, E. J. Baerends, and J. G. Snijders, J. Chem. Phys. 99, 4597 (1993). 18E. van Lenthe, E. J. Baerends, and J. G. Snijders, J. Chem. Phys. 101, 9783 (1994). 19E. van Lenthe, A. W. Ehlers, and E. J. Baerends, J. Chem. Phys. 110, 8943 (1999). 20W. Kutzelnigg and W. Liu, J. Chem. Phys. 123, 241102 (2005). 21W. Kutzelnigg and W. Liu, Mol. Phys. 104, 2225 (2006). 22W. Liu and W. Kutzelnigg, J. Chem. Phys. 126, 114107 (2007). 23W. Liu and D. Peng, J. Chem. Phys. 125, 044102 (2006). 24D. Peng, D. Liu, Y. Xiao, and L. Cheng, J. Chem. Phys. 127, 104106 (2007). 25M. Iliaš and T. Saue, J. Chem. Phys. 126, 064102 (2007). 26D. Peng and M. Reiher, Theor. Chem. Acc. 131, 1081 (2012). 27D. Peng and M. Reiher, J. Chem. Phys. 136, 244108 (2012). 28D. Peng, N. Middendorf, F. Weigend, and M. Reiher, J. Chem. Phys. 138, 184105 (2013). 29W. Liu and D. Peng, J. Chem. Phys. 131, 031104 (2009). 30W. Xu, J. Ma, D. Peng, W. Zou, W. Liu, and V. Staemmler, Chem. Phys. 356, 219 (2009). 31W. Xu, Y. Zhang, and W. Liu, Sci. China, Ser. B: Chem. 52, 1945 (2009). 32Y. Zhang, W. Xu, Q. Sun, W. Zou, and W. Liu, J. Comput. Chem. 31, 532 (2010). 33H. J. Aa. Jensen, R. Bast, T. Saue, L. Visscher, V. Bakken, K. G. Dyall, S. Dubillard, U. Ekström, E. Eliav, T. Enevoldsen, T. Fleig, O. Fossgaard, A. S. P. Gomes, T. Helgaker, J. K. Lærdahl, Y. S. Lee, J. Henriksson, M. Iliaš, Ch. R. Jacob, S. Knecht, S. Komorovský, O. Kullie, C. V. Larsen, H. S. Nataraj, P. Norman, G. Olejniczak, J. Olsen, Y. C. Park, J. K. Pedersen, M. Pernpointner, K. Ruud, P. Sałek, B. Schimmelpfennig, J. Sikkema, A. J. Thorvaldsen, J. Thyssen, J. van Stralen, S. Villaume, O. Visser, T. Winther, and S. Yamamoto, DIRAC, a relativistic ab initio electronic structure program, Release DIRAC12 (2012), see http://www.diracprogram.org. 34M. Kühn and F. Weigend, J. Chem. Theory Comput. 9, 5341 (2013). 35Local version of TURBOMOLE V6.6 2014, a development of University of Karlsruhe and Forschungszentrum Karlsruhe GmbH, 1989-2007, TURBOMOLE GmbH since 2007. 36K. A. Peterson, D. Figgen, M. Dolg, and H. Stoll, J. Chem. Phys. 126, 124101 (2007). 37D. Figgen, K. A. Peterson, M. Dolg, and H. Stoll, J. Chem. Phys. 130, 164108 (2009). 38D. Figgen, G. Rauhut, M. Dolg, and H. Stoll, Chem. Phys. 311, 227 (2005). 39B. Metz, H. Stoll, and M. Dolg, J. Chem. Phys. 113, 2563 (2000). 40K. A. Peterson, D. Figgen, E. Goll, H. Stoll, and M. Dolg, J. Chem. Phys. 119, 11113 (2003). 41I. S. Lim, P. Schwerdtfeger, B. Metz, and H. Stoll, J. Chem. Phys. 122, 104103 (2005). 42I. S. Lim, H. Stoll, and P. Schwerdtfeger, J. Chem. Phys. 124, 034107 (2006). 43D. E. Jiang, M. Kühn, Q. Tang, and F. Weigend, J. Phys. Chem. Lett. 5, 3286 (2014). 44M. Kühn and F. Weigend, J. Chem. Phys. 141, 224302 (2014). 45P. Tecmer, R. Bast, K. Ruud, and L. Visscher, J. Phys. Chem. A 116, 7397 (2012). 46M. E. Casida, in Recent Advances in Density Functional Methods, Part I, edited by D. P. Chong (World Scientific, Singapore, 1995), pp. 155–192. 47F. Wang and W. Liu, J. Chin. Chem. Soc. 50, 597 (2003). 48F. Wang and T. Ziegler, J. Chem. Phys. 121, 12191 (2004). 49R. Bauernschmitt and R. Ahlrichs, Chem. Phys. Lett. 256, 454 (1996). 50F. Furche and D. Rappoport, Computational Photochemistry, Theoretical and Computational Chemistry edited by M. Olivucci (Elsevier, Amsterdam, 2005), Vol. 16, p. 93. 9D.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
034116-8 51M.
M. Kühn and F. Weigend
Kühn, J. Chem. Theory Comput. 10, 623 (2014). Oddershede and P. Jørgensen, J. Chem. Phys. 66, 1541 (1977). 53T. D. Bouman and A. E. Hansen, J. Chem. Phys. 66, 3460 (1977). 54M. E. Casida and M. Huix-Rotllant, Annu. Rev. Phys. Chem. 63, 287 (2012). 55M. Grüning, A. Marini, and X. Gonze, Nano Lett. 9, 2820 (2009). 56D. Peng, J. Ma, and W. Liu, Int. J. Quantum Chem. 109, 2149 (2009). 57See supplementary material at http://dx.doi.org/10.1063/1.4905829 for detailed description of the implemented 2c TDA algorithm and coordinates of all structures. 58C. van Wüllen, J. Comput. Chem. 32, 1195 (2010). 59A. D. Becke, Phys. Rev. A 38, 3098 (1988). 60J. P. Perdew, Phys. Rev. B 33, 8822 (1986). 61F. Weigend and A. Baldes, J. Chem. Phys. 133, 174102 (2010). 62K. A. Peterson, B. C. Shepler, D. Figgen, and H. Stoll, J. Phys. Chem. A 110, 13877 (2006). 63A. D. Becke, J. Chem. Phys. 98, 5648 (1993). 64F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys. 7, 3297 (2005). 65D. Andrae, U. Häußermann, M. Dolg, H. Stoll, and H. Preuß, Theor. Chim. Acta 77, 123 (1990). 66J. P. Perdew and Y. Wang, Phys. Rev. B 45, 13244 (1992). 52J.
J. Chem. Phys. 142, 034116 (2015) 67M.
Kühn and F. Weigend, ChemPhysChem 12, 3331 (2011). Eichkorn, F. Weigend, O. Treutler, and R. Ahlrichs, Theor. Chem. Acc. 97, 119 (1997). 69R. Bauernschmitt, M. Häser, O. Treutler, and R. Ahlrichs, Chem. Phys. Lett. 264, 573 (1997). 70F. Weigend, Phys. Chem. Chem. Phys. 8, 1057 (2006). 71M. K. Armbruster, W. Klopper, and F. Weigend, Phys. Chem. Chem. Phys. 8, 4862 (2006). 72O. Treutler and R. Ahlrichs, J. Chem. Phys. 102, 346 (1995). 73M. A. Baldo, S. Lamansky, P. E. Burrows, M. E. Thompson, and S. R. Forrest, Appl. Phys. Lett. 75, 4 (1999). 74W. J. Finkenzeller and H. Yersin, Chem. Phys. Lett. 377, 299 (2003). 75J. Fine, K. Diri, A. I. Krylov, C. Nemirow, Z. Lu, and C. Wittig, Mol. Phys. 110, 1849 (2012). 76C. E. Moore, Atomic Energy Levels Vol. 1-3, Natl. Bur. Stand. (U.S.) Circ. No. 467 (U.S. Government Printing Office, Washington, DC, 1949, 1952, 1958). 77J. Tellinghuisen, J. Chem. Phys. 58, 2821 (1973). 78J. Tellinghuisen, J. Chem. Phys. 76, 4736 (1982). 79K. P. Huber and G. Herzberg, Molecular Spectra and Molecular Structure (Van Nostrand and Reinhold, New York, 1979), Vol. IV. 68K.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 205.175.97.18 On: Fri, 13 Mar 2015 08:22:12
