Article pubs.acs.org/JCTC
Comparison of Real-Time and Linear-Response Time-Dependent Density Functional Theories for Molecular Chromophores Ranging from Sparse to High Densities of States Samat Tussupbayev,† Niranjan Govind,*,‡ Kenneth Lopata,§ and Christopher J. Cramer*,† †
Department of Chemistry, Supercomputing Institute, and Chemical Theory Center, University of Minnesota, Minneapolis, Minnesota 55455, United States ‡ Environmental Molecular Sciences Laboratory, Pacific Northwest National Laboratory, Richland, Washington 99338, United States § Department of Chemistry and Center for Computation & Technology, Louisiana State University, Baton Rouge, Louisiana 70803, United States S Supporting Information *
ABSTRACT: We assess the performance of real-time timedependent density functional theory (RT-TDDFT) for the calculation of absorption spectra of 12 organic dye molecules relevant to photovoltaics and dye-sensitized solar cells with 8 exchange-correlation functionals (3 traditional, 3 global hybrids, and 2 range-separated hybrids). We compare the calculations with traditional linear-response (LR) TDDFT and experimental spectra. In addition, we demonstrate the efficacy of the RT-TDDFT approach to calculate wide absorption spectra of two large chromophores relevant to photovoltaics and molecular switches. RT-TDDFT generally requires longer simulation times, compared to LR-TDDFT, for absorption spectra of small systems. However, it becomes more effective for the calculation of wide absorption spectra of large molecular complexes and systems with very high densities of states.
1. INTRODUCTION The time-dependent response of molecules and materials under the influence of external fields is at the center of several fundamental processes, such as photodissociation, organic and inorganic photovoltaic systems, electron transport, charge injection, higher-order nonlinear response, etc. In a nutshell, light-matter interactions can be classified as either weak, when the interaction with the field is much smaller than the magnitude of the intermolecular interactions, and strong otherwise. In the former, the excitation is weak and results in a small perturbation away from the ground state. In this regime, perturbation theories such as linear-response or frequency domain density functional theory (LR-TDDFT) provide reliable descriptions of the linear-response properties of molecules and materials (e.g., electronic absorption spectra).1,2 Strong interactions, on the other hand, require methodologies that go beyond linear response. Over the past few years, there has been growing interest in real-time formulations and real-time time-dependent density functional theory (RT-TDDFT) in particular. It is an alternative approach to the traditional response formulation which allows one to capture the effects of external fields of molecules and materials in both the weak and strong interaction regimes in a seamless manner. Unlike the weak limit, where frequency domain (or linear response) perturbative © XXXX American Chemical Society
approaches are sufficient, the strong limit is best captured via a real-time, real-space, and time-correlation approach. In addition to being able to simulate transient spectroscopies, where full dynamical information about the excitation process is needed, the nonlinear spectral information of molecules and materials beyond perturbation theory can also be calculated. Due in part to its growing popularity, many implementations of RT-TDDFT have been reported in the literature3−14 and have been used in several applications such as molecular conductance,12,15 excited states at metal surfaces,16 absorption properties of silicon clusters,9 double excitations,17 singlet− triplet transition,18 electric circular dichroism (ECD) spectroscopy,19 and magnetic circular dichroism (MCD) spectroscopy.20 There has also been extensive work in developing schemes that go beyond the Born−Oppenheimer (BO) approximation to explicitly treat nuclear motion, such as Ehrenfest dynamics,21,22 Liouville−von Neumann (LvN) molecular dynamics with real-time tight binding,6 the use of surface hopping between adiabatic states,23−25 and correlated electron−ion dynamics.26 Recent implementations have also been reported in combination with solvation models,27−29 applications to core-level spectra,30 optical spectra of metalReceived: August 23, 2014
A
DOI: 10.1021/ct500763y J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation oxides31 and with non-Hermitian von Neumann density matrix propagation with absorbing potentials to model near and above ionization electronic excitations and core spectra.32,63 An important feature of RT-TDDFT is its ability to capture excitations across broad spectral regions (5−20 eV) in large systems and systems with high densities of states (DOSs). For example, the optical absorption spectrum of transition-metal oxide materials such as α-Fe2O3 typically involves closely spaced excitations between the O 2p and unoccupied transition-metal (Fe) 3d states.31 Traditional methods such as LR-TDDFT may be used to address the lowest excited states. However, these “bottom-up” iterative approaches become computationally prohibitive if large numbers of excitations (∼103−104) are involved, because of the tetradic nature of the RPA matrix.33 Formally, the numerical cost to diagonalize the full LR-TDDFT equations scales as O(N6). Similar issues are encountered in the calculation of core-level spectra (XANES).30 RT-TDDFT, in combination with a delta-function field, addresses this by simultaneously exciting all the electronic modes with the length of the simulation determining the resolution of the spectral features. In this context, we should also mention the complex polarization,34 damped response approaches,35,36 and multishift linear solvers,37 which have also been shown to be viable approaches to address wide absorption spectra. Other approaches include the self-consistent constricted variational DFT approach by Ziegler and co-workers38 and the simplified Tamm-Dancoff39 and simplified TDDFT reported by Grimme and co-workers.40 In this work, we present a systematic benchmark study, using our previously reported RT-TDDFT implementation in NWChem,13 of the absorption spectra (or weak field limit) of 12 organic dye molecules relevant for photovoltaics and dyesensitized solar cells using 8 exchange-correlation functionals (3 traditional, 3 global hybrids, and 2 range-separated hybrids). These calculations are compared with traditional LR-TDDFT. We also compare the RT-TDDFT calculated spectrum with the gas phase experimental spectrum for the Zn-phtalocyanine (ZnPc) system using the B3LYP and M06 functionals. Finally, to the best of our knowledge, for the first time, we calculate the wide absorption spectra of two large organic molecular complexes41 relevant to photovoltaics and molecular switches. The primary goal of this paper is to demonstrate the consistency of the two approaches over a broad range of systems, exchange-correlation functionals and the utility of the RT-TDDFT approach for the calculation of spectra of large systems. We do not address the shortcomings of the different exchange-correlation functionals. The rest of the paper is organized as follows: we first present a brief overview of the method, followed by the computational details and a discussion of our results.
i
αii(ω) =
i
μi (ω) Ei(ω)
i = x, y, z
(4)
where μi(ω) and Ei(ω) are the Fourier transforms of the dipole moment and applied electric fields in the i direction (the offdiagonals are assumed to be zero). The absorption cross section tensor is then given by 4πω σ (ω) = Im[α(ω)] (5) c and the resulting dipole strength function is then S(ω) =
1 Tr[σ (ω)] 3
(6)
It is worth emphasizing that the RT-TDDFT approach is dependent only on the exchange-correlation potential via the time-dependent Fock matrix, unlike the LR-TDDFT approach, which requires the full second derivative of the exchangecorrelation functional as part of the response equations. This not only simplifies the implementation, but also removes the dependence on higher derivatives and potential numerical instabilities that may stem from them. As such, this approach can also be used to test the accuracy of response properties and excited states from potential-driven density functional theory, i.e., for exchange-correlation functionals which have no energy functional or for which higher derivatives have not yet been implemented. For further details of the implementation and information about time propagators, we refer the interested reader to our original paper.13 Computational Details. All calculations were performed with the NWChem 6.3 software suite,42 except for linearresponse (LR-TDDFT) calculations employing meta-GGA functionals, which were done using the Gaussian 09 software suite.43 In total, eight exchange-correlation functionals were considered and were divided into three groups, depending on the amount of HF exchange. Three local functionals were considered (BLYP,44,45 M06-L,46 and M11-L47 (which has dual-range DFT exchange)), along with three global hybrids (B3LYP (VWN5 version),45,48−50 M06, and M06-2X,51 with 20%, 27%, and 54% Hartree−Fock (HF) exchange, respectively), and two range-separated hybrids (CAM-B3LYP52 and M11,53 with 19% and 42.8% HF exchange, respectively, in the short range and 65% and 100% HF exchange, respectively, in the long range). The dyes considered in this study vary in size from 14 atoms for the smallest (2,3,5-trifluorobenzaldehyde) to 162 atoms for the largest ( f-coronene) and cover both organic and organometallic compounds. They include 2,3,5-trifluorobenzaldehyde (1), dibenzonaphthyridindione (2), coronene (3), coumarin
NMO
∑ C*μi(t )Cνi(t )
(3)
where F′(t) is the time-dependent Fock matrix, which also includes the applied field, in the MO basis. For every time step, P′ is propagated forward in time using the second-order Magnus scheme.13 For convenience, relative intensities are typically computed directly from the Fourier transform of the dipole moment; however, in principle, RT-TDDFT spectra yield meaningful absolute values of the oscillator strengths identical to those of LR-TDDFT. First, the on-diagonals of the complex polarizability tensor are computed from x,y,z kick-type excitation simulations:
2. METHODOLOGY In our Gaussian basis set implementation,13 we use the singleparticle reduced density matrix representation, which is given by P′μν(t ) =
∂P′(t ) = [F′(t ), P′(t )] ∂t
(2)
where Cμi(t) is the time-dependent molecular orbital coefficient matrix. In the molecular orbital (MO) basis, the time evolution of the density matrix is governed by the von Neumann equation: B
DOI: 10.1021/ct500763y J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
calculated. The ground-state single-point DFT calculations to generate the initial density for RT-TDDFT were run with a tighter density convergence (RMS difference of
