Development of a TDDFT-Based Protocol with Local Hybrid

Sep 1, 2017 - Chromophores suitable for singlet fission need to meet specific requirements regarding the relative energies of their S0, S1, and T1 (an...
2 downloads 15 Views 704KB Size
Article pubs.acs.org/JCTC

Development of a TDDFT-Based Protocol with Local Hybrid Functionals for the Screening of Potential Singlet Fission Chromophores Robin Grotjahn,† Toni M. Maier,† Josef Michl,‡,¶ and Martin Kaupp*,† †

Institut für Chemie, Theoretische Chemie/Quantenchemie, Sekr. C7, Technische Universität Berlin, Straße des 17. Juni 135, D-10623 Berlin, Germany ‡ Institute of Organic Chemistry and Biochemistry, Academy of Sciences of the Czech Republic, Flemingovo nám. 2, 16110 Prague 6, Czech Republic ¶ Department of Chemistry and Biochemistry, University of Colorado, Boulder, Colorado 80309-0215, United States S Supporting Information *

ABSTRACT: Chromophores suitable for singlet fission need to meet specific requirements regarding the relative energies of their S0, S1, and T1 (and T2) electronic states. Accurate quantumchemical computations of the corresponding energy differences are thus highly desirable for materials design. Methods based on density functional theory (DFT) have the advantage of being applicable to larger, often more relevant systems compared to more sophisticated post-Hartree−Fock methods. However, most exchange−correlation functionals do not provide the needed accuracy, in particular, due to an insufficient description of the T1 state. Here we use a recent singlet fission chromophore test set (Wen, J.; Havlas, Z.; Michl, J. J. Am. Chem. Soc. 2015, 137, 165−172) to evaluate a wide range of DFT-based methods, with an emphasis on local hybrid functionals with a position-dependent exact-exchange admixture. New reference vertical CC2/CBS benchmark excitation energies for the test set have been generated, which exhibit somewhat more uniform accuracy than the previous CASPT2-based data. These CC2 reference data have been used to evaluate a wide range of functionals, comparing full linear-response TDDFT, the Tamm−Dancoff approximation (TDA), and ΔSCF calculations. Two simple two-parameter local hybrid functionals and the more empirical M06-2X global meta-GGA hybrid provide the overall best accuracy. Due to its lower empiricism and wide applicability, the Lh12ct-SsifPW92 local hybrid is suggested as the main ingredient of an efficient computational protocol for prediction of the relevant excitation energies in singlet fission chromophores. Full TDDFT for the S1, S2, and T2 excitations is combined with ΔSCF for the T1 excitations. Making use also of some error compensation with suitable DFT-optimized structures, even the most critical T1 excitations can be brought close to the target accuracy of 0.20 eV, while the other excitation energies are obtained even more accurately. This fully DFT-based protocol should become a useful tool in the field of singlet fission.

1. INTRODUCTION In the singlet fission process,1−3 a part of the excitation energy of a singlet excited chromophore (S1 state) is transferred to a neighboring chromophore in its S0 ground state, yielding both chromophores in their first triplet excited state (T1 state). Thus, theoretically enabling triplet quantum yields close to 200%, exploitation of the singlet fission process allows, in principle, the design of single-junction solar cells (e.g., dye-sensitized solar cells, DSSCs4) with an efficiency beyond the Shockley− Queisser limit of 32%.5 For that purpose, it is crucial that the process of singlet fission be faster than competing processes depopulating the S1 state of the excited chromophore. Furthermore, annihilation of triplet excited chromophores needs to be sufficiently slow to enable their use for charge separation. © XXXX American Chemical Society

To identify chromophores that exhibit efficient singlet fission, certain criteria have been formulated.6 Besides the determination of an optimal structure for interchromophore coupling and the identification of conditions that minimize the quenching of chromophores in the T1 state, the energy ordering of the electronic states in the monomer involved in the singlet fission process is of major importance. In particular, the first singlet excitation is required to lie approximately within 0.1−0.2 eV of twice the excitation energy from the S0 to the T1 state. In addition, the T2 excitation energy should be larger than twice the T1 excitation energy to ensure endothermicity and thus slow rates of the T1 + T2 to T1 annihilation process.2 Received: July 6, 2017 Published: September 1, 2017 A

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

out to be the method of choice regarding computational efficiency.30 In particular, a recent seminumerical linearresponse TDDFT implementation of local hybrid functionals has been shown to be even more efficient than conventional analytical integration (as measured for global hybrids), already for relatively small molecules and basis sets.31 The objective of the present work is to investigate the performance of several state-of-the-art local hybrid functionals for the calculation of excitation energies relevant for singlet fission of potentially suitable chromophores. In particular, vertical excitations from the S0 ground state to the first two singlet and triplet excited states, that is, the S1, S2, T1, and T2 states, will be considered here for that purpose. Although not relevant for the singlet fission process, S2 excitation energies will be considered for comparison. Regarding the quantumchemical methodology, we aim to profit from the beneficial properties of local hybrid functionals, both the significantly improved description of triplet valence excitations compared to conventional XC functionals and the high efficiency of the recent seminumerical local hybrid TDDFT implementation. Therefore, we make use of the chromophore set by Wen et al.,6 which contains 11 captodatively stabilized biradicaloids modified by donor and acceptor groups to fine tune the energy levels of the excited states. Lewis structures and numbering of the chromophores are given in Figure 1.

There is a rather small number of chromophores known to satisfy these energy requirements.2 Common examples are polycyclic aromatic hydrocarbons such as pentacene, which has already been used in a DSSC prototype.7 Another class of molecules that are able to satisfy the energy requirements are certain captodatively stabilized biradicaloids.6 The development and investigation of further novel chromophores meeting the formulated energy requirements is of fundamental importance because there are additional physical, technological, and economical prerequisites for the successful industrial application of singlet fission in solar cells.2 An efficient way to screen for chromophores that satisfy the energy requirements is calculation of the respective excitation energies using quantum-chemical methods. Besides multireference methods, which explicitly consider excited-state determinants, such as multireference configuration interaction8 and complete active space (CAS) perturbation theory,9 linearresponse methods based on coupled-cluster (CC) theory10,11 and density functional theory (DFT)12,13 have earned appreciable popularity for the calculation of excited states.14,15 Among these methods, time-dependent density functional theory (TDDFT) in its linear-response formalism exhibits the best ratio between computational effort and achievable accuracy and thus became the method of choice when it comes to calculations on large molecular systems or when computational efficiency is crucial. Nevertheless, the accuracy of linear-response TDDFT calculations strongly depends on the choice of the employed exchange−correlation (XC) kernel and of the underlying potential and energy functional.15,16 Apart from the fact that LR-TDDFT describes only single excitations within the commonly used adiabatic approximation,17 that is, the XC kernel is defined as the second functional derivative of the XC functional with respect to the ground-state electron density, Rydberg,18,19 charge-transfer,20,21 core,22 and triplet valence excitations15 represent known cases that are difficult to describe with common XC functionals. A relatively recent class of XC functionals that have been shown to give promising results for some of these cases, in particular, triplet valence, core, and Rydberg excitations,23 are local hybrid functionals (local hybrids).24−26 Local hybrid functionals feature a local admixture of exact exchange managed by the real-space-dependent local mixing function (LMF) g(r), giving the XC functional ExcLH = Exex + Ecsl +

∑ ∫ [1 − gσ (r)]·[ex,slσ(r) − ex,exσ(r)] dr

Figure 1. Lewis structures of the chromophores studied, taken from ref 6. The inset shows a schematic structure of captodatively stabilized 1,2-biradicaloids.

σ

(1)

with eslx,σ(r) being a (semi)local exchange energy density, Eslc a (semi)local correlation functional, Exex the exact-exchange energy, and eex x,σ(r) the exact-exchange energy density in its conventional gauge27,28 e x,exσ(r) = −

1 2

occ.

∑∫ kl

φk*, σ (r)φl , σ (r)φl*, σ (r′)φk , σ (r′) |r − r′|

The paper is organized as follows: First, we present new CC2 reference values for the vertical excitation energies of the chromophore set and compare them to previous CASPT2 results by Wen et al.6 In this regard, special attention is paid to potential multireference character of the chromophores, which is essential for the reliability of the CC2 results. Subsequently, the new reference values are employed for benchmarking stateof-the-art local hybrid functionals against a large variety of common XC functionals. Besides full TDDFT, the Tamm− Dancoff approximation (TDA)32 as well as the ΔSCF approach33 are investigated. Finally, the influence of the underlying ground-state structures is studied in more detail. On the basis of these investigations, a DFT-based protocol for

dr′ (2)

The LMF contains information about the inhomogeneity of the system’s electron density and thus allows a more flexible exactexchange admixture than the mixing constant in conventional global hybrid functionals.29 For the determination of nonstandard integrals appearing due to the weighting of exact exchange with the LMF, seminumerical integration has turned B

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation the efficient screening of potential singlet fission chromophores is formulated.

CCSD(T), were carried out with the ricc2 module55−57 of the TURBOMOLE 6.6 program package.58 For the determination of vertical excitation energies within the DFT framework, three different methods were applied. Besides linear-response TDDFT within the adiabatic approximation,12 which has been applied with and without employing the TDA,32 the ΔSCF33 approach has been investigated. In the latter case, excited states are treated within the conventional ground-state Kohn−Sham formalism, yielding excitation energies as the energy differences to the ground state. Because the ΔSCF method thus requires convergence of the excited states during the self-consistent field (SCF) calculation, usually only a few low-lying excited states, which are often distinguished by molecular or spin symmetry, are accessible. In contrast to TDDFT, with which we consider the first two singlet and triplet excited states, only the S1 and T1 states of the chromophores are easily accessible with the ΔSCF method and are thus considered. However, due to the applied unrestricted Kohn−Sham formalism, the calculated S1 states exhibit unphysical spin contamination caused by mixing with the T1 state. For determination of the energy of the spin-pure S1 state, the energy projection technique described by Labanowski et al.59 was used. To evaluate the performance of local hybrid functionals, we have chosen five local hybrids known to give good results for thermochemistry and the calculation of vertical excitation energies.23,25,60 Two types of LMFs have been used: one of them is based on the ratio t = τW/τ of the von Weizsäcker kinetic energy density τW and the local kinetic energy density τ26 (“t-LMF”), and the other (“s-LMF”) is based on the dimensionless density gradient s.61 Specifically, we investigated three local hybrids based on Slater exchange62 and VWN correlation,63 that is, Lh07t-SVWN26 employing a spin channel t-LMF with a prefactor of 0.48, Lh07s-SVWN25,61 with a spinchannel s-LMF, that is, gσ(r) = erf(β·sσ(r)) with β = 0.277 (note that the parameter β from ref 25 has been adjusted to fit the original definition of sσ(r)64), and Lh12ct-SVWN60 using a common t-LMF with a prefactor of 0.534. In addition, two local hybrids based on self-correlation-corrected PW92 correlation and a common t-LMF were considered.65 They included a selfcorrelation-reduced variant, Lh12ct-SsirPW92,60 with a LMF prefactor of 0.646 and the self-correlation-free Lh12ctSsifPW9260 functional with a LMF prefactor of 0.709. For comparison, various established XC functionals have also been considered. Apart from SVWN,62,63 the GGA functionals BLYP,66,67 BP86,66,68 and PBE,64 as well as the meta-GGA functional TPSS,34,35 this includes the global hybrid functionals (listed in order of increasing exact-exchange admixture) TPSSh (10%),36 B3LYP (20%),37,38 PBE0 (25%),39 BLYP35 (35%),40 BMK (42%),69 BHLYP (50%),29 and M06-2X (54%).70 In addition, we also evaluated the range-separated hybrid functionals CAM-B3LYP (19−65%, ω = 0.33),71 LC-ωPBE (0−100%, ω = 0.40),72 and ωB97X-D (22.2036−100%, ω = 0.20),73 where the range separation parameter ω and the corresponding short-range and long-range exact-exchange admixture are given in parentheses. For all DFT and TDDFT calculations, the RI-J approximation47,74 and the def2-QZVPD basis set75 were used. The basis set has been chosen so as to keep the expected basis set error for all investigated excitation energies below 0.01 eV. As an example, this was checked for molecule 1 with TDDFT calculations using the B3LYP functional. Excitation energies at the CBS limit were determined using eq 3 and Jensen’s aug-pc-

2. COMPUTATIONAL DETAILS Unless stated otherwise, CASPT2/ANO-S-VDZP chromophore structures from ref 6, optimized with point group symmetry, were employed. In addition, DFT ground-state structures optimized with different XC functionals, in particular, TPSS,34,35 TPSSh,36 B3LYP,37,38 PBE0,39 and BLYP35,40 have been evaluated. These optimizations were carried out using the def2-TZVP basis set,41 a medium-sized integration grid (TURBOMOLE grid 3),42 and standard analytical integration for the determination of four-center integrals. The influence of Grimme’s DFT-D3 dispersion corrections43 on the structure optimization has been tested for the B3LYP functional but was found to result in negligible changes in structure and in the subsequent excitation energy calculations and will therefore not be further considered in this work. Potential multireference character of the chromophores in their ground states was evaluated using three different diagnostic tools. D1 diagnostics, proposed by Janssen and Nielsen,44 were calculated at the MP2/aug-cc-pVQZ, CC2/augcc-pVQZ, and CCSD/aug-cc-pVDZ levels of theory. Fractional orbital density (FOD) plots45 were generated with the program FODEN46 on the basis of B3LYP/def2-QZVPD calculations with fractional occupation at a Fermi temperature of 9000 K employing the resolution-of-the-identy (RI) approximation for Coulomb integrals (RI-J)47 and seminumerical integration of exact exchange48−50 on a small grid (Turbomole grid 1).42 Additionally, %TAEe[(T)] values,51 that is, the percentage of the perturbative triplet correction (T) to the total CCSD(T) atomization energy, were evaluated. Because %TAEe[(T)] values have been shown to exhibit a weak basis set dependence,51 which has been exemplarily confirmed for molecule 2 with the aug-cc-pVTZ and aug-cc-pVQZ basis sets, the corresponding CCSD(T) calculations used the aug-ccpVDZ basis set. For the calculation of vertical CC2 excitation energies,10 Dunning’s augmented, correlation-consistent aug-cc-pVXZ (X = D,T,Q) basis sets52 together with their corresponding auxiliary basis sets53 were used. Excitation energies at the complete basis set (CBS) limit were calculated by separately extrapolating reference and correlation energies of the ground state and the excited states, respectively. While the Hartree− Fock energy provides the reference energy for the ground state, the reference energy of the excited states was determined as a sum of the ground-state Hartree−Fock energy and the coupledcluster singles (CCS) excitation energy. For basis set extrapolation, the three-point exponential fitting formula54 E(l) = a ·exp( −b·l) + E CBS

(3)

was employed, where l is the highest angular momentum quantum number of the basis set (i.e., here l = 2, 3, 4) and ECBS is the extrapolated energy in the basis set limit. Additionally, CC2 excitation energies were evaluated with the def2-TZVP basis set.41 While four-center integrals were determined analytically for the Hartree−Fock reference in CC2 and CCSD(T) calculations, the RI approximation was applied in the respective coupled-cluster steps.55 Double excitation character of CC2 excitations is measured by the %;2 diagnostics.56 Coupled-cluster calculations, including CC2 and C

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation X (X = 2,3,4) basis sets,76−78 where ground- and excited-state energies were again extrapolated separately. For the evaluation of basis set dependencies, TDDFT excitation energies were calculated additionally with the def2-QZVP, def2-TZVPD, def2-TZVP, def2-SVPD, and def2-SVP basis sets75 using the Lh12ct-SsifPW92 functional. The current density response in XC functionals depending on the kinetic energy density79 was not considered, because the influence on the investigated excitations is expected to be negligible.23 Except for BMK, M06-2X, and the range-separated hybrid functionals, all DFT and TDDFT calculations were carried out with a developers’ version of TURBOMOLE 6.658 containing the recent seminumerical implementation of local hybrid functionals.30,31 Seminumerical evaluation of exact-exchange integrals has been used for local hybrid functionals and also for global hybrid functional calculations in TURBOMOLE (except for structure optimizations) using small molecular integration grids (TURBOMOLE grid 1).42 Calculations with BMK, M06-2X, and the range-separated hybrid functionals were done with Gaussian 09 (without RI-J).80 For statistical analysis of the various XC functionals, we give the mean absolute error (MAE) and the mean signed error (MSE) of DFT/TDDFT excitation energies with respect to CC2 excitation energies at the CBS limit.

exhibit larger energies than their corresponding CASPT2 counterparts. While absolute deviations are nonetheless mostly below 0.40 eV (illustrated by the dashed lines in Figure 2), the deviation in molecule 4 significantly exceeds this value with a difference of 0.52 eV. For S2 and T2 excitations, CC2 mostly gives systematically larger excitation energies than CASPT2, with absolute deviations commonly lying below 0.40 eV. However, the S2 excitation energy of molecule 10 as well as the T2 excitation energies of molecules 1 and 7 are considerably lower with CC2. Except for these outliers, CC2 and CASPT2 excitation energies agree well. To identify reasons for the outliers and to evaluate the reliability of the new CC2/CBS excitation energies, the methodological details of both the CASPT2 and CC2 calculations need to be scrutinized: first, the effect of the basis set extrapolation shall be examined. Therefore, we have compared CC2 excitation energies at the CBS limit with those using the def2-TZVP basis set (see Table S1 in the Supporting Information). In general, CC2/CBS and CC2/def2-TZVP excitation energies differ by less than 0.04 eV, which is significant for a quantitative evaluation but does not change the qualitative description of the respective excitations. Here, it appears that singlet excitation energies are usually lowered in energy upon going from def2-TZVP to the CBS limit, while triplet excitation energies are mostly increased. However, the deviations are significantly larger for the S2 and T2 excitations of molecules 1 and 2, with differences of 0.22 and 0.29 eV, and 0.07 and 0.05 eV, respectively. To a large extent, this behavior can be attributed to the structure of these two molecules, which differ from those of the other considered chromophores, specifically by the presence of two nitrile groups. In particular, it appears that the nitrile groups promote the mixing with low-lying Rydberg states, which results in an enhanced diffuseness of the respective T2 and S2 states. Hence, diffuse basis sets are inherently required for an accurate description. Accordingly, the basis set errors of CC2/aug-ccpVTZ for the excitations in question are below 0.04 eV, so that the larger basis set errors with the def2-TZVP basis set can be directly attributed to the absence of sufficiently diffuse basis functions. CC2 and CASPT2 differ significantly in how excitation energies are calculated. In CASPT2 calculations, excited states are considered by explicitly determining higher-lying roots of the configuration-interaction (CI) matrix, which is built in the underlying CASSCF calculation. Hence, excitations are restricted to those that are included in the CAS. Additionally, a larger CAS might be required for an accurate description compared to the ground state to account for potential multireference character of excited states or the partial mixing of higher-lying states. Hence, the calculation of excitation energies with CASPT2 is highly sensitive to the choice of CAS orbitals.81 In the case that the CAS is adequately chosen, CASPT2 is able to accurately treat degeneracy effects due to potential multireference character of the investigated molecules and thus also gives an accurate description of static correlation contributions. On the other hand, dynamic correlation is mainly described with second-order perturbation theory (PT2), which is similar to the MP2 method. An improved description of dynamical correlation within CASPT2 is only possible by increasing the CAS size, which is often not feasible due to the exponential scaling of the calculation cost. Excitation energies within the CC2 method are calculated within a linear-response approach. That is, excitation vectors describing the transition from the ground state to the excited

3. RESULTS 3.1. CC2 Reference Excitation Energies. Before turning to the discussion of the performance of the various XC functionals, the new CC2 vertical excitation energies, which will be used as reference for the subsequent DFT results, shall be evaluated and compared to the previous CASPT2 excitation energies by Wen et al.6 The corresponding correlation plot between CASPT2/ANO-L-VTZP and CC2/CBS vertical excitation energies is shown in Figure 2. For the CC2 excitations, S1, S2, T1, and T2 states were newly assigned in each calculation, while the assignment from ref 6 was taken for the CASPT2 results. For T1 excitations, differences between CC2 and CASPT2 excitation energies are relatively small and approximately statistically distributed, with a maximum absolute deviation of 0.13 eV. In contrast, CC2 S1 excitations usually

Figure 2. Correlation plot between vertical CC2/CBS reference excitation energies ECC2 of the present work and previous vertical CASPT2/ANO-L-VTZP excitation energies ECASPT2 by Wen et al.6 for the S1, S2, T1, and T2 excitations of the chromophores 1−11. The dashed lines mark deviations of ±0.4 eV. For excitations outside of this area, the chromophore number is given. D

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

While the factor between the highest and lowest D1 diagnostics, irrespective of the method on which it is based, amounts to around 1.4; the corresponding factor for the NFOD values is around 3.4. This suggests that the integrated FOD values are not directly correlated with the amount of static correlation for the present systems. Nevertheless, on the basis of the FOD plots, the π-systems could be identified as the main “hot spots” within chromophores 1−11. For the D1 diagnostic itself, thresholds of 0.04 and 0.05, based on MP2 and CC2 or CCSD calculations, respectively,44 are considered to indicate an inadequate method. In fact, this is the case for all investigated molecules. Like the FOD, this suggests distinct multireference cases. As a third multireference indicator, we used the %TAEe[(T)] diagnostic proposed by Martin and co-workers.51 Calculated values range from 1.960 to 2.853. According to the originally proposed and later confirmed thresholds,83 this indicates only low or mild static correlation and thus no pronounced multireference character. This is in distinct contrast to the other diagnostic tools and thus has to be evaluated carefully. As had been already pointed out by Karton et al.,83 D1 diagnostics may fail to identify multireference cases and are thus not suited as a sole indicator. Evaluation of several multireference diagnostics for transition-metal complexes by Wilson and coworkers84,85 further indicated that D1 thresholds need to be adjusted depending on the systems investigated. In particular, D1 thresholds (for CCSD) of 0.15 and 0.12 had been proposed for 3d and 4d transition-metal complexes, respectively. While these thresholds are not directly applicable to the present chromophore set, they suggest that the calculated D1 values might not indicate multireference character. FOD diagnostics are somewhat doubtful here, too, due to their low correlation with D1 and %TAEe[(T)] values. On the other hand, %TAEe[(T)] values are directly related to the reliability of CC2 calculations,84 which is in fact at the heart of the present methodological evaluation. Because %TAEe[(T)] values have turned out to be a robust and reliable diagnostic in various studies,51,83−85 their present relatively low values suggest that molecules 1−11 do not exhibit significant multireference character and that CC2 ground-state calculations are reliable. To make sure we have no excessive double-excitation character in the relevant excitations, which would render the reliability of CC2 excitation energies insufficient, we calculated %;2 diagnostics for the S1 excitations (Table 1). Except for molecule 10, which exhibits a %;2 value of 15.49 for the S2 excitation, S1 excitations appear to exhibit the largest values. In view of a proposed threshold of 15%,56 all S1 excitations appear to be free from excessive double-excitation character. Even the S2 excitation of molecule 10 exceeds this threshold only slightly. Due to a different assignment of the S2 state of molecule 10 within CC2 and CASPT2, there is also no contradiction with the large deviation between both methods for the corresponding S2 excitation energy (see Figure 2; we note that inclusion of the S3 state into the state averaging of the CASPT2 calculation lowers the “S2” excitation energy by more than 1 eV). Hence, we conclude that CC2 is reliable for calculation of the vertical excitation energies investigated. Furthermore, the predominant single-excitation character also suggests well-behaved TDDFT results within the adiabatic approximation used here. In summary, CC2 generally appears to be reliable for calculation of the excitation energies in question. Despite the biradicaloid character of the singlet fission chromophores studied here, neither significant multireference nor extensive

states are optimized, instead of explicitly considering the excited states as done in CASPT2. Hence, there are no pitfalls comparable to the choice of an appropriate CAS. On the other hand, CC2 is clearly a single-reference method, so that degeneracy effects due to potential multireference character in the ground state of the investigated molecule are not covered. Hence, CC2 is not applicable to such multireference molecules. Nevertheless, mixing of various determinants in the excited states represents no problem as long as the singleexcitation character of the respective excitation lies above 85%.56 Dynamic as well as parts of the static correlations are described within a conventional CC singles and doubles (CCSD) approach, with the double-excitation operator being approximated to first order. While static correlation is thus usually treated less well compared to CASPT2, dynamical correlation is generally better described than that with a conventional PT2 approach. Single and double excitations are correct to second and zeroth order, respectively. For valence excitations with predominant single-excitation character, CC2 and CCSD are closely correlated,14 with expected CC2 errors compared to experiments of around 0.2−0.3 eV.82 Accordingly, CC2 calculations are only reliable if multireference character in the ground state can be excluded and the vertical excitations considered do not exhibit extensive doubleexcitation character. Therefore, we used several diagnostics (cf. the Computational Details section), that is, the integrated FOD value NFOD, D1 diagnostics, and the %TAEe[(T)] for evaluation of the multireference character (i.e., the importance of static correlation) as well as the %;2 values of the S1 excitations for the assessment of the double-excitation character. The results are shown in Table 1. FOD plots can be found in Figure S1 in the Supporting Information. Table 1. Comparison of Different Diagnostics for the Evaluation of Potential Multireference Character of Molecules 1−11a D1 diagnostic mol.

NFOD

MP2

CC2

CCSD

%TAEe[(T)]

%;2 (S1)

1 2 3 4 5 6 7 8 9 10 11

0.416 0.311 0.717 0.618 0.562 0.995 0.861 0.903 0.625 1.057 0.786

0.045 0.044 0.056 0.056 0.058 0.061 0.058 0.058 0.057 0.060 0.056

0.084 0.077 0.093 0.104 0.101 0.110 0.100 0.096 0.103 0.102 0.091

0.059 0.056 0.070 0.077 0.075 0.079 0.073 0.072 0.075 0.075 0.070

2.218 2.853 2.422 1.966 1.960 2.056 2.201 2.299 2.220 2.248 2.257

11.18 10.83 12.23 12.00 11.94 13.03 12.76 12.40 12.13 13.67 11.98

a

Fractional orbital densities integrated over space NFOD, D1 diagnostics for MP2, CC2, and CCSD calculations, as well as % TAEe[(T)] values are presented. The double-excitation character of S1 excitations is evaluated using the %;2 diagnostic.

Calculated NFOD values for the chromophores lie within a range of 0.311 for molecule 2 to 1.057 for molecule 10, suggesting a relatively wide distribution of multireference character. For all molecules, the FOD plots are distinctly delocalized, which, according to ref 45, indicates true multireference cases and thus strong static correlation. Turning to the more established D1 diagnostic, the relatively weak correlation (R2 = 0.33) with the NFOD values is most noticeable. E

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

Table 2. MAEs (MSEs) in eV for Excitation Energies Relative to CC2/CBS Reference Values Calculated with TDDFT, TDATDDFT, and ΔSCF, with Different XC Functionalsa TDDFT functional SVWN BLYP BP86 PBE TPSS TPSSh B3LYP PBE0 BLYP35 BMK BHLYP M06-2X CAM-B3LYP LC-ωPBE ωB97X-D Lh07s-SVWN Lh07t-SVWN Lh12ct-SVWN Lh12ct-SsirPW92 Lh12ct-SsifPW92

ΔSCF

TDA-TDDFT

S1

S2

T1

T2

S1

S2

T1

T2

S1

T1

0.83 (−0.83) 0.75 (−0.75) 0.76 (−0.76) 0.77 (−0.77) 0.67 (−0.67) 0.48 (−0.48) 0.37 (−0.37) 0.30 (−0.29) 0.15 (−0.09) 0.12 (−0.06) 0.17 (0.17) 0.08 (0.05) 0.09 (−0.02) 0.25 (0.24) 0.10 (−0.05) 0.38 (−0.38) 0.30 (−0.30) 0.25 (−0.24) 0.20 (−0.17) 0.18 (−0.13)

0.94 (−0.94) 0.85 (−0.85) 0.86 (−0.86) 0.88 (−0.88) 0.72 (−0.72) 0.47 (−0.47) 0.34 (−0.34) 0.26 (−0.26) 0.09 (−0.01) 0.07 (−0.07) 0.30 (0.27) 0.06 (−0.06) 0.12 (0.05) 0.26 (0.26) 0.09 (0.03) 0.38 (−0.38) 0.26 (−0.26) 0.19 (−0.19) 0.12 (−0.12) 0.08 (−0.08)

0.77 (−0.77) 0.78 (−0.78) 0.83 (−0.83) 0.81 (−0.81) 0.82 (−0.82) 0.81 (−0.81) 0.74 (−0.74) 0.82 (−0.82) 0.83 (−0.83) 0.59 (−0.59) 1.12 (−1.12) 0.47 (−0.47) 0.79 (−0.79) 1.06 (−1.06) 0.71 (−0.71) 0.71 (−0.71) 0.67 (−0.67) 0.56 (−0.56) 0.52 (−0.52) 0.50 (−0.50)

1.00 (−1.00) 0.95 (−0.95) 0.99 (−0.99) 1.01 (−1.01) 0.89 (−0.89) 0.68 (−0.68) 0.53 (−0.53) 0.50 (−0.50) 0.34 (−0.34) 0.26 (−0.26) 0.29 (−0.28) 0.17 (−0.17) 0.25 (−0.25) 0.20 (−0.20) 0.20 (−0.20) 0.52 (−0.52) 0.43 (−0.43) 0.30 (−0.30) 0.22 (−0.22) 0.17 (−0.17)

0.57 (−0.57) 0.50 (−0.50) 0.51 (−0.51) 0.52 (−0.52) 0.42 (−0.42) 0.30 (−0.24) 0.25 (−0.13) 0.22 (−0.06) 0.18 (0.13) 0.18 (0.17) 0.39 (0.39) 0.31 (0.31) 0.20 (0.20) 0.48 (0.48) 0.19 (0.18) 0.25 (−0.14) 0.25 (−0.06) 0.24 (0.00) 0.22 (0.07) 0.21 (0.10)

0.92 (−0.92) 0.83 (−0.83) 0.85 (−0.85) 0.86 (−0.86) 0.71 (−0.71) 0.45 (−0.45) 0.31 (−0.31) 0.23 (−0.23) 0.11 (0.02) 0.05 (−0.02) 0.34 (0.32) 0.07 (0.02) 0.14 (0.09) 0.31 (0.31) 0.12 (0.06) 0.36 (−0.36) 0.24 (−0.24) 0.17 (−0.17) 0.09 (−0.09) 0.06 (−0.05)

0.74 (−0.74) 0.72 (−0.72) 0.75 (−0.75) 0.75 (−0.75) 0.70 (−0.70) 0.62 (−0.62) 0.54 (−0.54) 0.53 (−0.53) 0.41 (−0.41) 0.32 (−0.32) 0.31 (−0.31) 0.22 (−0.22) 0.38 (−0.38) 0.30 (−0.30) 0.36 (−0.36) 0.53 (−0.53) 0.48 (−0.48) 0.40 (−0.40) 0.35 (−0.35) 0.32 (−0.32)

0.98 (−0.98) 0.92 (−0.92) 0.96 (−0.96) 0.97 (−0.97) 0.84 (−0.84) 0.62 (−0.62) 0.48 (−0.48) 0.44 (−0.44) 0.19 (−0.19) 0.20 (−0.20) 0.11 (0.05) 0.11 (−0.11) 0.13 (−0.13) 0.06 (0.02) 0.12 (−0.12) 0.49 (−0.49) 0.39 (−0.39) 0.27 (−0.27) 0.19 (−0.19) 0.14 (−0.14)

0.83b (−0.83)b 0.75 (−0.75) 0.75 (−0.75) 0.76 (−0.76) 0.66 (−0.66) 0.45 (−0.45) 0.33 (−0.33) 0.25 (−0.25) 0.11 (−0.03) 0.11a (−0.03)a 0.29 (0.28) 0.11b (0.08)b 0.12b (0.04)b 0.36b (0.36)b 0.12b (0.02)b 0.34 (−0.34) 0.30 (−0.30) 0.23 (−0.23) 0.16 (−0.15) 0.13 (−0.11)

0.40 (−0.40) 0.42 (−0.42) 0.47 (−0.47) 0.46 (−0.46) 0.47 (−0.47) 0.48 (−0.48) 0.42 (−0.42) 0.47 (−0.47) 0.44 (−0.44) 0.33 (−0.33) 0.45 (−0.45) 0.27 (−0.27) 0.42 (−0.42) 0.46 (−0.46) 0.39 (−0.39) 0.39 (−0.39) 0.37 (−0.37) 0.30 (−0.30) 0.28 (−0.28) 0.27 (−0.27)

a Functionals are grouped as (semi)local, global hybrid, range-separated hybrid, and local hybrid (from top to bottom). bStructure 7 is not included due to convergence problems in some of the ΔSCF calculations for the S1 state (mainly with Gaussian09).

double-excitation character emerged. In view of the modest static correlation effects, benefits from the use of CASPT2 instead of CC2 are expected to be minor. In fact, the more sophisticated consideration of dynamical correlation and the ability to use large, diffuse basis sets to calculate excitation energies at the CBS limit (see section 2 for details) make CC2 seem even preferable. The few notable deviations between the previous CASPT2 excitation energies by Wen et al.6 and the present CC2 results likely arise from insufficient active space, state-averaging, dynamical correlation errors, or misassignment of S2 or T2 states. However, the previous CASPT2 results for the S1 and T1 excitations, which are most important for evaluation of potential singlet fission chromophores, mostly agree well with the new CC2 excitation energies. The overall findings of ref 6 thus remain untouched. In the following

evaluation of local hybrid functionals, we will use the new CC2/CBS values as reference. 3.2. Evaluation of DFT Methods. Here we assess the performance of several state-of-the-art local hybrid functionals, in comparison with some (semi)local functionals, as well as global and range-separated hybrids, based on the new CC2/ CBS reference data of section 3.1. We make use of three different DFT-based methods for the evaluation of vertical excitation energies: TDDFT, TDA-TDDFT, and ΔSCF (see section 2). S1, S2, T1, and T2 excitations of the chromophores are analyzed separately. We restrict ourselves to a discussion of MAEs and MSEs (Table 2) and note that maximum errors are in most cases highly correlated with the MAEs (individual excitation energies can be found in Tables S2−S4 in the Supporting Information). F

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

agreement with the TDDFT results, but the agreement depends somewhat on the EXX admixture. While the ΔSCF and TDDFT MAEs differ only by 0.03 eV for TPSSh (10%), the difference is 0.12 eV for BHLYP (50%; differences are between 0.01 and 0.05 eV for the other global hybrids). The reason is that the energy projection scheme for the S1 state becomes less accurate for larger EXX admixtures due to mixing in of higher spin multiplets (M06-2X and BMK are exceptions; see also below). Dependence on the EXX admixture is not as straightforward for T1 excitations (Table 2). Notably, TDDFT errors for “simpler” global hybrid functionals, for example, TPSSh, B3LYP, PBE0, and BLYP35, appear to be comparable to those of (semi)local functionals, with a similar systematic underestimation. In fact, BHLYP performs even worse here than (semi)local functionals. In contrast, BMK and, particularly, M06-2X exhibit significantly lower errors. Most likely, this can be attributed to their high parametrization and to the inclusion of τ-dependent meta-GGA contributions rather than to the relatively large amount of exact exchange (note that the EXX admixtures of M06-2X and BHLYP, 54% vs 50%, are rather similar). Nevertheless, even M06-2X, the clearly best-performing global hybrid, shows a significant underestimation of 0.47 eV for T1 excitations at the TDDFT level. Improvement of the agreement for T1 excitations with the reference by using the TDA depends notably on the EXX admixture; the effects range from a very small reduction of the MAE for (semi)local functionals to much larger reductions with the largest EXX admixtures (e.g., by more than 0.80 eV for BHLYP). This is known to be an effect of overcoming triplet instabilities86,87 by neglecting the B matrix17 within the TDA. The instabilities and thus the effects of the TDA have been shown to increase with increasing amounts of EXX admixture.86,88 Interestingly, the effects of the TDA for the T1 excitation energies with BMK and, particularly, with M06-2X are smaller than expected for their large EXX admixtures. This suggests that these functionals buffer some of the adverse effects of the large EXX admixtures (triplet instability of the ground state) already at the TDDFT level. At the TDA level, M06-2X exhibits the lowest MAE of 0.22 eV, even though excitation energies are still systematically too low. Turning to the ΔSCF results for T1 excitations, differences between the global hybrids are less pronounced than those at the TDA level. The performance is overall much better than that with TDDFT and for low EXX admixtures also superior to TDA: only with BMK and M06-2X, TDA performs slightly better than ΔSCF. Trends of the ΔSCF T1 excitation energies concerning the EXX admixture are somewhat similar to those at the TDDFT level. This suggests that the different behavior of TDA-TDDFT might be an artifact of neglecting the B matrix. Again, M06-2X is the best-performing global hybrid functional with ΔSCF, giving a MAE of 0.27 eV for T1 excitations. Turning to the three range-separated hybrid functionals studied, we see good TDDFT performance for S1, S2, and T2 states, in particular, for ωB97X-D and CAM-B3LYP (Table 2). These two functionals feature an average EXX admixture comparable to BMK. LC-ωPBE, which interpolates between zero and full EXX admixture and consequently requires a larger range separation parameter (see the Computational Details section), may be viewed as having more average exact exchange, comparable to that of BHLYP. Indeed, the errors of the rangeseparated hybrids for these three excitation types can be compared to those global hybrids. However, at the TDDFT

Since we aim at establishing a predictive protocol for the screening of potential singlet fission chromophores, the accuracy of the calculated excitation energies is of primary importance. In particular, we require the S1 and T1 excitations to be described well because both states are directly involved in the singlet fission process. Therefore, we aim for an accuracy of 0.20 eV, which represents the expected precision of our CC2/ CBS reference values. We note in passing that DFT as well as CC2 results were ensured to be within 0.01 eV of the basis set limit, so that no significant basis set errors are expected to influence the outcome. In contrast, T2 states are usually just required to lie distinctly above the respective T1 states to exclude competing processes (see section 1). Hence, a lower accuracy is generally acceptable, although in borderline cases, for example, when the T2 excitation energy is close to twice the T1 excitation energy, higher accuracies might be necessary. Because the biradicaloid singlet fission chromophores investigated tend to place very specific demands on the quantumchemical methodology, and these have not been extensively studied in a DFT framework so far, results of the entire set of XC functionals will be discussed in detail. However, special emphasis is put on the comparison with the local hybrid functionals. The (semi)local functionals exhibit a systematic underestimation of all considered excitation energies, with MAEs up to or exceeding 1.00 eV (Table 2). In comparison to previous studies,15,16,23 which found MAEs of around 0.50 eV for common valence excitations for (semi)local functionals, errors are thus significantly larger here. This illustrates the larger methodological demands of the present chromophores. Systematically lower MSEs compared to previous studies have also been found for all other investigated XC functionals, which has a significant influence on the individual performances. Regarding differences between the methods employed, TDATDDFT shows the expected systematical blue shift of excitation energies compared to TDDFT.23,86 While S2 and T2 excitations appear to be generally less affected by this effect, especially S1 excitations show blue shifts of up to 0.26 eV. In the case of (semi)local XC functionals, the influence of the TDA on the T1 excitations is less pronounced. S1 excitation energies calculated with ΔSCF are in excellent agreement with the TDDFT results. This shows that the employed energy projection scheme (section 2) used for the S1 states, as well as the ΔSCF method in general, is able to give reliable results. On the other hand, the ΔSCF MAEs for T1 excitation energies are significantly lower than the respective TDDFT values. Because the calculated T1 states do not exhibit significant spin contamination, and S1 excitations, which use the respective T1 state for energy projection, appear to be well described, we tend to trust the ΔSCF values more for T1 excitations and attribute errors of the TDDFT results to the use of the adiabatic approximation. For the global hybrid functionals considered, a systematic increase of the MSEs for the S1, S2, and T2 excitations with an increasing amount of exact exchange is observed (Table 2). This agrees with trends in previous studies.15,16,23 Due to the generally lower MSEs for the chromophores investigated (see above), M06-2X with its 54% EXX admixture provides the lowest TDDFT errors, significantly below the target accuracy of 0.20 eV for the S1, S2, and T2 excitations. Due to the known blue shift, in this case, TDA in fact produces inferior results (except for the T1 excitations). For most global hybrid functionals, ΔSCF S1 excitation energies are again in good G

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

does not degrade these results, which contrasts with our findings for M06-2X. The MAE for T1 excitations is 0.27 even at the ΔSCF level, that is, still slightly above the desired threshold of 0.20 eV. This is nevertheless better than those with most of the other functionals studied, which remain much further above this threshold. The only exception is M06-2X, which provides a closely comparable performance, and BMK also comes close. Note also that the Lh12ct-SsirPW92 local hybrid, having a slightly smaller prefactor of 0.646 for the tLMF, performs almost as well as Lh12ct-SsifPW92. We see these two local hybrids thus as clear contenders for the best functional to be used for the purpose of studying the relevant singlet and triplet excitations of singlet fission chromophores. A clear advantage of Lh12ct-SsifPW92 over functionals like M062X is its much lower degree of empiricism. Only the t-LMF prefactor should be viewed as an adjustable parameter here, while the performance of the functional depends only marginally on the value of the range separation parameter in the modified PW92 correlation part.60 We also note that both Lh12ct-SsifPW92 and Lh12ct-SsirPW92 perform better than M06-2X for other vertical excitation test sets, including valence, Rydberg, and core excitations.23 Nevertheless, some excitations in the present test set, in particular, S0 → S1, are described slightly better by M06-2X. To complete our setup of a computational protocol to be used for larger-scale studies, we need not only select a functional but also decide on the use of TDDFT, TDA, or ΔSCF. TDDFT is the method of choice for S1, S2, and T2 excitations. It does not rely on potentially unfavorable approximations of the TDA and lacks the uncertainties associated with the energy projection scheme of ΔSCF for S1 states (which also sometimes may cause convergence problems). For T1 excitations, ΔSCF represents an invaluable tool because it eliminates errors arising from the adiabatic approximation in TDDFT. Due to spin symmetry, a SCF calculation for the T1 state is usually feasible for all kinds of chromophores with a singlet ground state, thus rendering ΔSCF a good choice for the calculation of vertical T1 excitation energies. While the TDA shows comparable accuracy when using Lh12ct-SsifPW92, it falls slightly behind ΔSCF. Besides, in a previous study of the Thiel test set, we found that the systematic blue shift of the TDA can also deteriorate the agreement for triplet excitations, in particular, for XC functionals such as Lh12ct-SsirPW92 and Lh12ct-SsifPW92, which already exhibit small errors within full TDDFT.23 That is, the blue shift provided by the TDA may be either favorable or unfavorable, depending on system and excitation. Here we prefer ΔSCF over the TDA for T1 excitations to ensure optimal applicability of our protocol. The proposed protocol is based on the use of the Lh12ctSsifPW92 local hybrid functional and uses TDDFT for S1 and T2 excitations and ΔSCF for T1 excitations. We note that Truhlar and co-workers recently used a similar strategy with full linear-response TDDFT for singlet and ΔSCF for triplet excitations when computing semiconductor band gaps and molecular excitation energies.89 While we have used the very large def2-QZVPD basis sets here to work essentially at the basis set limit, it is of interest to know what errors will be introduced by the use of smaller basis sets, which may be required in applications to very large chromophores. Further computations over the entire test set with Lh12ct-SsifPW92 at the TDDFT level afforded the following MAEs relative to the def2-QZVPD results (Table S5 provides individual excitation

level, they perform rather poorly for T1 excitations, LC-ωPBE in particular. The TDA provides similarly large improvements here as those for the corresponding global hybrids (much less so for the other excitations). In fact, large errors of LRC hybrids for some triplet excitations using full TDDFT and significant improvements with TDA-TDDFT have been observed previously.88 Similar trends are found for ΔSCF (Table 2). Overall, the range-separated hybrids perform well but lag somewhat behind the best global hybrids (in particular, M062X; see above) for T1 excitations. Finally turning to local hybrid functionals, we see trends (Table 2) that are similar to those that we found in a recent evaluation for valence singlet and triplet excitations of a typical test set of organic chromophores.23 The MSEs here are lower, as we found also above, resulting in systematically underestimated excitation energies for the present diradicaloid chromophores at the TDDFT level (and thus larger MAEs). The generally better performance of the t-LMF compared to the s-LMF is also in line with our previous results. It is related to differences in the real-space behavior of these two types of LMFs, especially in the valence region. Therefore, among the local hybrids, Lh07s-SVWN exhibits the largest MAEs for all excitations, with errors nonetheless already comparable to those of B3LYP. We also confirm that common LMFs generally perform better for triplet valence excitations than their spinchannel counterparts due to an improved implicit description of some static correlation contributions.23 That is, results for T1 and T2 excitations are inferior with Lh07s-SVWN and Lh07tSVWN (both using a spin-channel LMF) compared to those with the other three local hybrids (all using a common-t-LMF). This becomes most clearly apparent when comparing Lh07tSVWN and Lh12ct-SVWN: the t-LMF prefactors of these two functionals differ by only 0.054, but the MAEs for triplet excitations are noticeably lower with the common LMF. For all excitations, we observe increasing MSEs (and thus lower MAEs) with increasing prefactor of the t-LMF. While the latter is directly related to the amount of EXX admixture, which explains the dependence of S1, S2, and T2 excitations, also T1 excitations appear to be affected, which is not seen with global hybrid functionals (except when using the TDA; see above). It appears that a position-dependent EXX admixture with a common LMF, as opposed to a constant mixing parameter, introduces some additional static correlation effects needed to accurately describe the T1 excitations. An increased t-LMF prefactor extends the influence of the LMF and thus improves the description of the T1 excitations. Effects of the TDA are similar to those described above for global hybrids, that is, relatively small effects on S1, S2, and T2 excitations and a somewhat more pronounced increase of the MSE (and thus decrease of MAE) for T1 excitations. Similarly, ΔSCF provides only a little improvement for S1 excitations but a more noticeable one for T1 excitations, in fact even slightly superior to the TDA results (Table 2). Overall, the Lh12ct-SsifPW92 functional, which uses a common t-LMF with a large prefactor of 0.709, exhibits the lowest MAEs among all local hybrid functionals examined. While the MAEs for singlet excitations are comparable to those of BLYP35, triplet excitations are described similarly well as those with M06-2X, the best-performing global hybrid functional. Most notably, the MAEs for S1, S2, and T2 excitations determined with TDDFT are all below the threshold of 0.20 eV (Table 2). Due to the systematic underestimation of the respective excitation energies, TDA H

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

for the S1 and T1 excitations, which are much more important for the purpose of chromophore screening (see above), clearly improve along the series of DFT-optimized structures, suggesting PBE0 or BLYP35 as good choices. The slightly increased deviations for S2 and T2 at these structures are comparably unimportant. We note in passing that BLYP35 was set up originally as a pragmatic approach to quantify localization/delocalization in mixed-valence systems.40,92 Turning to a fully (TD)DFT-based protocol, we evaluate the same DFT-based structures also for vertical excitation energies obtained with Lh12ct-SsifPW92, our “best-choice” functional (section 3.2), against CC2/CBS reference values obtained at the CASPT2/ANO-S-VDZP structures. As pointed out above (section 3.2), S1, S2, and T2 excitations are calculated with TDDFT and T1 excitations with ΔSCF. MAEs are illustrated in Figure 4.

energies): 0.008 eV (def2-QZVP), 0.004 eV (def2-TZVPD), 0.021 eV (def2-TZVP), 0.056 eV (def2-SVPD), and 0.139 eV (def2-SVP). Clearly, inclusion of diffuse functions is decisive, and def2-TZVPD or def2-SVPD basis sets appear to be good choices for larger systems. 3.3. Effects of Input Structure. Our evaluations so far have relied on the CASPT2/ANO-S-VDZP-optimized groundstate structures of ref 6. To derive a purely (TD)DFT-based protocol for the screening of potential singlet fission candidates, we need to also investigate the use of DFT-optimized structures (so far, only for the ground state; studies of adiabatic excitation energies will also require excited-state structure optimizations). We have thus reoptimized the ground-state structures of all systems of this study using XC functionals with different amounts of EXX admixture, that is, TPSS, TPSSh, B3LYP, PBE0, and BLYP35 (Tables S6−S10 give Cartesian coordinates). We then investigated the effect of structural differences on the CC2/CBS vertical excitation energies (see Table S11 for explicit values). Local hybrid structure optimizations, which have been reported only recently,90 have not been performed here because no significant benefits were expected (these closed-shell ground-state structures are not particularly demanding methodologically, in contrast to other cases studied elsewhere91). MSEs of the CC2/CBS excitation energies at the DFT-optimized structures are compared with CC2/CBS results at CASPT2/ANO-S-VDZP structures in Figure 3.

Figure 4. MAEs in eV for TDDFT (S1, S2, and T2) and ΔSCF (T1) results at the Lh12ct-SsifPW92/def2-QZVPD level at structures optimized with various XC functionals and with CASPT2, relative to CC2/CBS//CASPT2/ANO-S-VDZP reference values.

Obviously, the use of structures obtained with low EXX admixtures adds to the deviations between the DFT- and CC2based calculations. While the effects are essentially negligible for S2 excitations, the use of a TPSS structure, for example, may increase the MAE for the S1 above our target threshold of 0.20 eV and also increase further the errors for T1 excitations. In contrast, PBE0 structures perform similary to the CASPT2 structures, leaving the comparison essentially unaffected. Interestingly, the use of BLYP35 structures appears to introduce some favorable error cancellation, which provides the overall lowest MAEs for S1, T1, and T2 excitations and also low MAEs for S2 excitations. Notably, the MAE for the T1 excitations has now moved close to the target threshold, while those for the other excitations remain clearly below. Obviously, the systematic underestimation of excitation energies (particularly to the T1 state) by our Lh12ct-SsifPW92-based protocol is partly compensated by slightly shorter bond lengths in the BLYP35 structures, which increases the excitation energies. We may thus make use of the non-negligible influence of the ground-state structure on the computed excitation energies and benefit from some favorable error compensation by using BLYP35 for structure optimization as part of a (TD)DFTbased computational protocol. We note in passing that we have also compared the MAEs for Lh12ct-SsifPW92 vs CC2/CBS when identical structures are

Figure 3. MSEs in eV of the CC2/CBS vertical S1, S2, T1, and T2 excitation energies calculated at various DFT-optimized structures, relative to CC2/CBS results computed at CASPT2/ANO-S-VDZP structures. The functionals have been ordered according to increasing EXX admixture.

With increasing EXX admixture in the functional used for structure optimization, the MSEs relative to the CASPT2 structures increase systematically for all four excitation types. Simultaneously, the bond lengths in the π-systems of the biradicaloid decrease. Interestingly, the S1 and T1 excitations exhibit very similar MSEs, as do the S2 and T2 excitations. This suggests a similar dependence on changes in the ground-state structure within each of the two groups of excitations. Differences between the structures may reflect errors of both CASPT2 optimizations (possibly due to state averaging or too small basis sets6) and the DFT optimizations. The differences between the two groups of excitations (about 0.10 eV in the MSEs) may in turn indicate a different structural dependence of the higher excited states and the lower-lying ones. The MSEs I

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

energies, the use of the BLYP35 functional gives even slightly shorter average bond lengths within the π-systems of the biradicaloid chromophores. This provides some favorable error compensation for the critical S1 and T1 excitation energies obtained with the Lh12ct-SsifPW92 functional, bringing even the most difficult T1 excitation energies almost within the target accuracy of 0.20 eV, while MAEs for the other excitation energies are clearly below this threshold. This leads us to propose the following DFT-based protocol for the calculation of vertical excitation energies in potential singlet fission chromophores: 1. Structure optimization at the BLYP35/def2-TZVP level. 2. TDDFT calculation of S1 and T2 excitation energies at the Lh12ct-SsifPW92/def2-QZVPD level. 3. ΔSCF calculation of the T1 excitation energy at the Lh12ct-SsifPW92/def2-QZVPD level. If smaller basis sets are required due to system size, the def2TZVPD basis set can be used with an essentially negligible loss of accuracy. We note in passing that seminumerical integration techniques render both DFT and TDDFT calculations with local hybrids extremely efficient, also for large basis sets. We expect that this protocol will provide useful accuracy for applications to a screening or design of singlet fission chromophores or of their aggregates, for example, in the construction of DSSCs. Further improvement can be expected from the ongoing development of additional sophisticated orbital-dependent XC functionals (e.g., by combining ideas of local hybrids and range separation). We furthermore plan to implement excited-state gradients for such functionals, providing access to more accurate adiabatic excitation energies.

used for both. In this case, the changes with different input structures are negligible and thus not discussed further (explicit values can be found in Table S12 in the Supporting Information).

4. CONCLUSIONS AND OUTLOOK A computationally expedient and reliable protocol for the calculation of vertical excitation energies of potential singlet fission chromophores has been developed and validated, based on (TD)DFT methodology. It promises to be applicable to larger chromophores and aggregates compared to more sophisticated approaches based on post-Hartree−Fock methodologies. Various methods were compared for the CASPT2based test set of Wen et al.6 We first established new CC2/CBS references for these systems as we came to the conclusion that CC2 provides more uniform high-quality data for the relevant excitations. Various diagnostics were used to ensure the absence of both significant multireference character of the ground state and of appreciable double-excitation character. The latter point is also important for the applicability of TDDFT with adiabatic kernels. CC2 is thus a good choice as a benchmark method. As a single-reference approach, it does not feature uncertainties regarding active spaces or starting orbitals, and its treatment of dynamical correlation is somewhat superior to that of CASPT2. Moreover, we have been able to obtain the CC2 data essentially at the complete basis set limit. We note in passing that the conclusions regarding the suitability of individual chromophores 1−11 for experimental work, as drawn by Wen et al.,6 ultimately remain unaffected, making an equivalent discussion of the new results redundant. On the basis of the generated CC2/CBS reference data, we have provided the first extensive evaluation of different XC functionals for the calculation of vertical excitation energies in potential singlet fission chromophores, comparing full TDDFT, TDDFT within the TDA, and ΔSCF calculations (for the S1 and T1 excitations). Special focus has been on evaluation of local hybrid functionals, which were recently found to exhibit promising results for a wide variety of excitation classes. In fact, the simple two-parameter local hybrid functional Lh12ctSsifPW92 (and the related Lh12ct-SsirPW92 functional) turned out to exhibit the best performance among the investigated XC functionals, except for the similarly accurate but much more highly parametrized (33 parameters) M06-2X global hybrid. In view of its lower empiricism and its excellent performance in previous wider validation, we chose Lh12ct-SsifPW92 over M06-2X as the central ingredient of the new protocol. The approach is based on full TDDFT for the S1, S2, and T2 excitations and on ΔSCF for the T1 excitations. While the TDA also improves the T1 excitations relative to full TDDFT, we prefer the ΔSCF method because it provides slightly more accurate results for the chosen local hybrids and does not rely on a systematic blue shift as in the TDA by neglecting the B matrix. In fact, the latter has been shown to systematically improve T1 excitations for the molecules considered in the current study, which however cannot be generalized to other cases. To arrive at a fully (TD)DFT-based approach suitable for larger systems, we have also evaluated different XC functionals for the optimization of the ground-state structures. The EXX admixture to the functional used in these optimizations has a non-negligible effect on the computed vertical excitation energies. While PBE0 structures provide closest agreement with the CASPT2 ones and thus very similar excitation



ASSOCIATED CONTENT

S Supporting Information *

This material is available free of charge via the Internet at http://pubs.acs.org/. The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.7b00699. Tables and figures with additional computational data (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Martin Kaupp: 0000-0003-1582-2819 Funding

This work has been funded by DFG project KA1187/14-1. J.M. thanks the Alexander von Humboldt Foundation for supporting a stay at TU Berlin and U.S. DOE BES DE-SC0007004 for funding. R.G. thanks the Studienstiftung des Deutschen Volkes and Fonds der chemischen Industrie for financial support. Notes

The authors declare no competing financial interest.



REFERENCES

(1) Singh, S.; Jones, W. J.; Siebrand, W.; Stoicheff, B. P.; Schneider, W. G. Laser Generation of Excitons and Fluorescence in Anthracene Crystals. J. Chem. Phys. 1965, 42, 330−342. (2) Smith, M. B.; Michl, J. Singlet Fission. Chem. Rev. 2010, 110, 6891−6936.

J

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation (3) Smith, M. B.; Michl, J. Recent Advances in Singlet Fission. Annu. Rev. Phys. Chem. 2013, 64, 361−386. (4) O’Regan, B.; Grätzel, M. A low-cost, high-efficiency solar cell based on dye-sensitized colloidal TiO2 films. Nature 1991, 353, 737− 740. (5) Shockley, W.; Queisser, H. J. Detailed Balance Limit of Efficiency of p-n Junction Solar Cells. J. Appl. Phys. 1961, 32, 510−519. (6) Wen, J.; Havlas, Z.; Michl, J. Captodatively Stabilized Biradicaloids as Chromophores for Singlet Fission. J. Am. Chem. Soc. 2015, 137, 165−172. (7) Congreve, D. N.; Lee, J.; Thompson, N. J.; Hontz, E.; Yost, S. R.; Reusswig, P. D.; Bahlke, M. E.; Reineke, S.; Van Voorhis, T.; Baldo, M. A. External Quantum Efficiency Above 100% in a Singlet-ExcitonFission-Based Organic Photovoltaic Cell. Science 2013, 340, 334−337. (8) Siegbahn, P. E. M. Generalizations of the direct CI method based on the graphical unitary group approach. II. Single and double replacements from any set of reference configurations. J. Chem. Phys. 1980, 72, 1647−1656. (9) Roos, B. O.; Andersson, K.; Fülscher, M. P. Towards an accurate molecular orbital theory for excited states: the benzene molecule. Chem. Phys. Lett. 1992, 192, 5−13. (10) Christiansen, O.; Koch, H.; Jørgensen, P. The second-order approximate coupled cluster singles and doubles model CC2. Chem. Phys. Lett. 1995, 243, 409−418. (11) Koch, H.; Jørgensen, P. Coupled cluster response functions. J. Chem. Phys. 1990, 93, 3333−3344. (12) Casida, M. E. In Recent Advances in Density Functional Methods, Part I; Chong, D. P., Ed.; World Scientific: Singapore, 1995; pp 155− 192. (13) Gross, E. K. U.; Maitra, N. T. Fundamentals of Time-Dependent Density Functional Theory; Springer, 2012; Chapter 4. (14) Schreiber, M.; Silva-Junior, M. R.; Sauer, S. P.; Thiel, W. Benchmarks for electronically excited states: CASPT2, CC2, CCSD, and CC3. J. Chem. Phys. 2008, 128, 134110. (15) 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. (16) Laurent, A. D.; Jacquemin, D. TD-DFT Benchmarks: A Review. Int. J. Quantum Chem. 2013, 113, 2019−2039. (17) Dreuw, A.; Head-Gordon, M. Single-Reference ab Initio Methods for the Calculation of Excited States of Large Molecules. Chem. Rev. 2005, 105, 4009−4037. (18) Tozer, D. J.; Handy, N. C. Improving virtual Kohn-Sham orbitals and eigenvalues: Application to excitation energies and static polarizabilities. J. Chem. Phys. 1998, 109, 10180−10189. (19) Tozer, D. J.; Handy, N. C. On the determination of excitation energies using density functional theory. Phys. Chem. Chem. Phys. 2000, 2, 2117−2121. (20) 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. (21) Dreuw, A.; Head-Gordon, M. Failure of Time-Dependent Density Functional Theory for Long-Range Charge-Transfer Excited States: The Zincbacteriochlorin-Bacteriochlorin and Bacteriochlorophyll-Spheroidene Complexes. J. Am. Chem. Soc. 2004, 126, 4007− 4016. (22) Besley, N.; Peach, M. J. G.; Tozer, D. J. Time-dependent density functional theory calculations of near-edge X-ray absorption fine structure with short-range corrected functionals. Phys. Chem. Chem. Phys. 2009, 11, 10350−10358. (23) Maier, T. M.; Bahmann, H.; Arbuznikov, A. V.; Kaupp, M. Validation of local hybrid functionals for TDDFT calculations of electronic excitation energies. J. Chem. Phys. 2016, 144, 074106. (24) Jaramillo, J.; Scuseria, G. E.; Ernzerhof, M. Local hybrid functionals. J. Chem. Phys. 2003, 118, 1068−1073.

(25) Kaupp, M.; Bahmann, H.; Arbuznikov, A. V. Local hybrid functionals: An assessment for thermochemical kinetics. J. Chem. Phys. 2007, 127, 194102. (26) Bahmann, H.; Rodenberg, A.; Arbuznikov, A. V.; Kaupp, M. A thermochemically competitive local hybrid functional without gradient corrections. J. Chem. Phys. 2007, 126, 011103. (27) Tao, J.; Springborg, M.; Perdew, J. P. Properties of the exchange hole under an appropriate coordinate transformation. J. Chem. Phys. 2003, 119, 6457−6464. (28) Tao, J.; Staroverov, V. N.; Scuseria, G. E.; Perdew, J. P. Exactexchange energy density in the gauge of a semilocal density-functional approximation. Phys. Rev. A: At., Mol., Opt. Phys. 2008, 77, 012509. (29) Becke, A. D. A new mixing of Hartree-Fock and local densityfunctional theories. J. Chem. Phys. 1993, 98, 1372−1377. (30) Bahmann, H.; Kaupp, M. Efficient Self-Consistent Implementation of Local Hybrid Functionals. J. Chem. Theory Comput. 2015, 11, 1540−1548. (31) Maier, T. M.; Bahmann, H.; Kaupp, M. Efficient Semi-numerical Implementation of Global and Local Hybrid Functionals for TimeDependent Density Functional Theory. J. Chem. Theory Comput. 2015, 11, 4226−4237. (32) Hirata, S.; Head-Gordon, M. Time-dependent density functional theory within the Tamm-Dancoff approximation. Chem. Phys. Lett. 1999, 314, 291−299. (33) Jones, R. O.; Gunnarsson, O. The density functional formalism, its applications and prospects. Rev. Mod. Phys. 1989, 61, 689−746. (34) Perdew, J. P.; Tao, J.; Staroverov, V. N.; Scuseria, G. E. Metageneralized gradient approximation: Explanation of a realistic nonempirical density functional. J. Chem. Phys. 2004, 120, 6898−6911. (35) Tao, J.; Perdew, J. P.; Staroverov, V. N.; Scuseria, G. E. Climbing the Density Functional Ladder: Nonempirical Meta-Generalized Gradient Approximation Designed for Molecules and Solids. Phys. Rev. Lett. 2003, 91, 146401. (36) Staroverov, V. N.; Scuseria, G. E.; Tao, J.; Perdew, J. Comparative assessment of a new nonempirical density functional: Molecules and hydrogen-bonded complexes. J. Chem. Phys. 2003, 119, 12129−12137. (37) Becke, A. D. Density functional thermochemistry. III. The role of exact exchange. J. Chem. Phys. 1993, 98, 5648−5652. (38) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. Ab Initio Calculation of Vibrational Absorption and Circular Dichroism Spectra Using Density Functional Force Fields. J. Phys. Chem. 1994, 98, 11623−11627. (39) Adamo, C.; Barone, V. Toward reliable density functional methods without adjustable parameters: The PBE0 model. J. Chem. Phys. 1999, 110, 6158−6170. (40) Renz, M.; Theilacker, K.; Lambert, C.; Kaupp, M. A Reliable Quantum-Chemical Protocol for the Characterization of Organic Mixed-Valence Compounds. J. Am. Chem. Soc. 2009, 131, 16292− 16302. (41) Weigend, F.; Ahlrichs, R. Balanced basis sets of split valence, triple zeta valence and quadruple zeta valence quality for H to Rn: Design and assessment of accuracy. Phys. Chem. Chem. Phys. 2005, 7, 3297−3305. (42) Treutler, O.; Ahlrichs, R. Efficient molecular numerical integration schemes. J. Chem. Phys. 1995, 102, 346−353. (43) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu. J. Chem. Phys. 2010, 132, 154104. (44) Janssen, C. L.; Nielsen, I. M. B. New diagnostics for coupledcluster and Møller-Plesset perturbation theory. Chem. Phys. Lett. 1998, 290, 423−430. (45) Grimme, S.; Hansen, A. A Practicable Real-Space Measure and Visualization of Static Electron-Correlation Effects. Angew. Chem., Int. Ed. 2015, 54, 12308−12313. (46) FODplot Tools, a development of Mulliken Center for Theoretical Chemistry, University of Bonn. http://www.thch.unibonn.de/tc/ (2015). K

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation (47) Eichkorn, K.; Treutler, O.; Ö hm, H.; Häser, M.; Ahlrichs, R. Auxiliary basis sets to approximate Coulomb potentials. Chem. Phys. Lett. 1995, 240, 283−290. (48) Neese, F.; Wennmohs, F.; Hansen, A.; Becker, U. Efficient, approximate and parallel Hartree-Fock and hybrid DFT calculations. A ’chain-of-spheres’ algorithm for the Hartree-Fock exchange. Chem. Phys. 2009, 356, 98−109. (49) Izsak, R.; Neese, F.; Klopper, W. Robust fitting techniques in the chain of spheres approximation to the Fock exchange: The role of the complementary space. J. Chem. Phys. 2013, 139, 094111. (50) Plessow, P.; Weigend, F. Seminumerical Calculation of the Hartree-Fock Exchange Matrix: Application to Two-Component Procedures and Efficient Evaluation of Local Hybrid Density Functionals. J. Comput. Chem. 2012, 33, 810−816. (51) Karton, A.; Daon, S.; Martin, J. M. L. W4−11: A highconfidence benchmark dataset for computational thermochemistry derived from first-principles W4 data. Chem. Phys. Lett. 2011, 510, 165−178. (52) Dunning, T. H., Jr. Gaussian basis sets for use in correlated molecular calculations. I. The atoms boron through neon and hydrogen. J. Chem. Phys. 1989, 90, 1007−1023. (53) Weigend, F.; Köhn, A.; Hättig, C. Efficient use of the correlation consistent basis sets in resolution of the identity MP2 calculations. J. Chem. Phys. 2002, 116, 3175−3183. (54) Feller, D. Benchmarks of improved complete basis set extrapolation schemes designed for standard CCSD(T) atomization energies. J. Chem. Phys. 2013, 138, 074103. (55) Hättig, C.; Weigend, F. CC2 excitation energy calculations on large molecules using the resolution of the identity approximation. J. Chem. Phys. 2000, 113, 5154−5161. (56) Hättig, C.; Köhn, A.; Hald, K. First-order properties for triplet excited states in the approximated coupled cluster model CC2 using an explicitly spin coupled basis. J. Chem. Phys. 2002, 116, 5401−5410. (57) Hättig, C.; Hald, K. Implementation of RI-CC2 triplet excitation energies with an application to trans-azobenzene. Phys. Chem. Chem. Phys. 2002, 4, 2111−2118. (58) TURBOMOLE V6.6 2014, a development of University of Karlsruhe and Forschungszentrum Karlsruhe GmbH, 1989−2007. TURBOMOLE GmbH http://www.turbomole.com (2007). (59) Ovchinnikov, A. A.; Labanowski, J. K. Simple spin correction of unrestricted density-functional calculation. Phys. Rev. A: At., Mol., Opt. Phys. 1996, 53, 3946−3952. (60) Arbuznikov, A. V.; Kaupp, M. Importance of the correlation contribution for local hybrid functionals: Range separation and selfinteraction corrections. J. Chem. Phys. 2012, 136, 014111. (61) Arbuznikov, A. V.; Kaupp, M. Local hybrid exchange-correlation functionals based on the dimensionless density gradient. Chem. Phys. Lett. 2007, 440, 160−168. (62) Slater, J. C. A Simplification of the Hartree-Fock Method. Phys. Rev. 1951, 81, 385−390. (63) 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. (64) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865−3868. (65) Perdew, J. P.; Wang, Y. Accurate and simple analytic representation of the electron-gas correlation energy. Phys. Rev. B: Condens. Matter Mater. Phys. 1992, 45, 13244−13249. (66) Becke, A. D. Density-functional exchange-energy approximation with correct asymptotic behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (67) 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: Condens. Matter Mater. Phys. 1988, 37, 785−789. (68) Perdew, J. P. Density-functional approximation for the correlation energy of the inhomogeneous electron gas. Phys. Rev. B: Condens. Matter Mater. Phys. 1986, 33, 8822−8824.

(69) Boese, A. D.; Martin, J. M. L. Development of density functionals for thermochemical kinetics. J. Chem. Phys. 2004, 121, 3405−3416. (70) Zhao, Y.; Truhlar, D. G. The M06 suite of density functionals for main group thermochemistry, thermochemical kinetics, noncovalent interactions, excited states, and transition elements: two new functionals and systematic testing of four M06-class functionals and 12 other functionals. Theor. Chem. Acc. 2008, 120, 215−241. (71) Yanai, T.; Tew, D. P.; Handy, N. C. A new hybrid exchangecorrelation functional using the Coulomb-attenuating method (CAMB3LYP). Chem. Phys. Lett. 2004, 393, 51−57. (72) Vydrov, O. A.; Scuseria, G. E. Assessment of a long-range corrected hybrid functional. J. Chem. Phys. 2006, 125, 234109. (73) Chai, J.-D.; Head-Gordon, M. Long-range corrected hybrid density functionals with damped atom-atom dispersion corrections. Phys. Chem. Chem. Phys. 2008, 10, 6615−6620. (74) Bauernschmitt, R.; Häser, M.; Treutler, O.; Ahlrichs, R. Calculation of excitation energies within time-dependent density functional theory using auxiliary basis set expansions. Chem. Phys. Lett. 1997, 264, 573−578. (75) Rappoport, D.; Furche, F. Property-optimized Gaussian basis sets for molecular response calculations. J. Chem. Phys. 2010, 133, 134105. (76) Jensen, F. Polarization consistent basis sets: Principles. J. Chem. Phys. 2001, 115, 9113−9125. (77) Jensen, F. Polarization consistent basis sets. III. The importance of diffuse functions. J. Chem. Phys. 2002, 117, 9234−9240. (78) Jensen, F. Polarization consistent basis sets. II. Estimating the Kohn-Sham basis set limit. J. Chem. Phys. 2002, 116, 7372−7379. (79) Bates, J. E.; Furche, F. Harnessing the meta-generalized gradient approximation for time-dependent density functional theory. J. Chem. Phys. 2012, 137, 164105. (80) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J.; J, A.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, O.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian09, revision E.01; Gaussian Inc.: Wallingford, CT, 2009. (81) Veryazov, V.; Malmqvist, P. Å.; Roos, B. O. How to Select Active Space for Multiconfigurational Quantum Chemistry? Int. J. Quantum Chem. 2011, 111, 3329−3338. (82) Christiansen, O.; Koch, H.; Jørgensen, P.; Helgaker, T. Integral direct calculation of CC2 excitation energies: singlet excited states of benzene. Chem. Phys. Lett. 1996, 263, 530−539. (83) Karton, A.; Rabinovich, E.; Martin, J. M. L.; Ruscic, B. W4 theory for computational thermochemistry: In pursuit of confident sub-kJ/mol predictions. J. Chem. Phys. 2006, 125, 144108. (84) Jiang, W.; DeYonker, N. J.; Wilson, A. K. Multireference Character for 3d Transition-Metal-Containing Molecules. J. Chem. Theory Comput. 2012, 8, 460−468. (85) Wang, J.; Manivasagam, S.; Wilson, A. K. Multireference Character for 4d Transition Metal-Containing Molecules. J. Chem. Theory Comput. 2015, 11, 5865−5872. (86) Peach, M. J. G.; Tozer, D. J. Overcoming Low Orbital Overlap and Triplet Instability Problems in TDDFT. J. Phys. Chem. A 2012, 116, 9783−9789. (87) Casida, M. E.; Gutierrez, F.; Guan, J.; Gadea, F.-X.; Salahub, D.; Daudey, J.-P. Charge-transfer correction for improved time-dependent local density approximation excited-state potential energy curves: L

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation Analysis within the two-level model with illustration for H2 and LiH. J. Chem. Phys. 2000, 113, 7062−7071. (88) Sears, J. S.; Koerzdoerfer, T.; Zhang, C.-R.; Brédas, J.-L. Communication: Orbital instabilities and triplet states from timedependent density functional theory and long-range corrected functionals. J. Chem. Phys. 2011, 135, 151103. (89) Verma, P.; Truhlar, D. G. HLE17: An Improved Local Exchange-Correlation Functional for Computing Semiconductor Band Gaps and Molecular Excitation Energies. J. Phys. Chem. C 2017, 121, 7144−7154. (90) Klawohn, S.; Bahmann, H.; Kaupp, M. Implementation of Molecular Gradients for Local Hybrid Density Functionals Using Seminumerical Integration Techniques. J. Chem. Theory Comput. 2016, 12, 4254−4262. (91) Kaupp, M.; Karton, A.; Bischoff, F. A. [Al2O4]−, a Benchmark Gas-Phase Class II Mixed-Valence Radical Anion for the Evaluation of Quantum-Chemical Methods. J. Chem. Theory Comput. 2016, 12, 3796−3806. (92) Parthey, M.; Kaupp, M. Quantum-chemical insights into mixedvalence systems: within and beyond the Robin-Day scheme. Chem. Soc. Rev. 2014, 43, 5067−5088.

M

DOI: 10.1021/acs.jctc.7b00699 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX