The modeling of optical spectra of plasmonic nanoparticles via first-principles approaches is computationally expensive; thus, methods with high accuracy/computational cost ratio are required. Here, we show that the Time-Dependent Density Functional Theory (TDDFT) approach can be strongly simplified if only one s-type function per atom is employed in the auxiliary basis set, with a properly optimized exponent. This approach (named TDDFT-as, for auxiliary s-type) predicts excitation energies for silver nanoparticles with different sizes and shapes with an average error of only 12 meV compared to reference TDDFT calculations. The TDDFT-as approach resembles tight-binding approximation schemes for the linear-response treatment, but for the atomic transition charges, which are here computed exactly (i.e., without approximation from population analysis). We found that the exact computation of the atomic transition charges strongly improves the absorption spectra in a wide energy range.

The Time-Dependent Density Functional Theory (TDDFT) is the most used first-principles tool for the calculation of excitation energies in electronic systems.1–5 TDDFT is nowadays implemented in many free or commercial codes and can, in principle, predict the exact excitation energies. In practice, its accuracy depends on the exchange–correlation (XC) functional. For molecular systems, TDDFT has an accuracy of about 0.2 eV–0.5 eV, depending on the accuracy of the XC functional.6 More recently, TDDFT has been applied to plasmonic systems,7–14 e.g., silver and gold nanoparticles (NPs), which play a key role for applications in nano-optics, photovoltaics, surface-enhanced Raman scattering, and imaging.15–19 The absorption spectrum of metal NPs is mainly characterized by a localized surface plasmon resonance (LSPR), which can be accurately described using classical electromagnetism,20,21 but only for large systems. TDDFT, instead, includes all quantum effects, which are not considered in classical models or only approximated in other theoretical approaches for quantum plasmonics.12,22–27

Thanks to the computational efficiency of TDDFT, calculations for systems with more than 100 atoms have been reported for silver28–39 and gold29,31,34,35,39–45 nanosystems. In one of the largest application reported so far,42 the TDDFT absorption spectrum of Au1414 has been computed, but using more than 60 000 cores. Despite its efficiency, thus, TDDFT can be hardly applied to systems with thousands of atoms. To this end, different methods and algorithms have been developed46–50 in order to reduce the computational scaling. TDDFT implementations can be distinguished into two main classes: (i) real-time TDDFT and (ii) linear-response TDDFT. In the real-time approach,51 the evolution of the time-dependent Kohn–Sham (KS) states is computed following a perturbation. Absorption spectra are, then, obtained as the Fourier transform of the time-dependent dipole moment of the system. In the linear-response approach, a Dyson equation is solved in the frequency space,3,5 and excitation energies are obtained as poles of the interacting density–density response, or alternatively, the excitation energies can be obtained directly solving an eigenvalue problem as first pointed out by Casida.2 The spectra of both methods are equivalent.51,52

In this work, we will focus on the Casida formalism and, in particular, on the construction of the response matrix, which requires the calculation of Coulomb and XC kernel integrals. A standard approach to reduce the TDDFT computational cost is the resolution-of-the-identity (RI) method, in which an auxiliary basis set is used to expand the electronic density for a faster computation of the Coulomb integrals.53–61 More recently, following ideas of Time-Dependent Density Functional Tight Binding (TD-DFTB) methods,62–66 approximated TDDFT approaches have been proposed: in the simplified Tamm-Dancoff approximation (sTDA),67 sTDDFT,68 and TDDFT + TB69 methods, the computation of kernel integrals is avoided completely using semiempirical formulas. Although the applications of those methods for organic molecules are well established, investigations for more complex systems such as metal NPs have appeared only recently.70–72 

In this work, we introduce a modification of the first-principles RI-TDDFT approach, in which the auxiliary basis set is composed of a single s-type function per atom. In the conventional RI approach, the size of the auxiliary basis set is very large as it is usually optimized to reproduce the total Coulomb energy with chemical accuracy.60,61,73

Which is the impact of such minimal auxiliary basis set on the absorption spectra?

We show that our method, named TDDFT-as, for auxiliary s-type, can predict excitation energies for different silver nanoparticles with an average accuracy of only 12 meV, which is really negligible considering the much larger errors related to the orbital basis set and the choice of the XC functional.

The TDDFT-as method resembles the TDDFT + TB approach. In this work, we analyze and compare these two approaches in details. We found that the accuracy of TDDFT-as for silver nanoparticles is distinctively superior, and the origin on the TDDFT + TB lower accuracy is clarified.

This article is organized as follows: In Sec. II, we review the RI-TDDFT and the TDDFTB + TB theories and discuss the conditions under which these approaches can be considered equivalent; in Sec. III, we describe the systems under investigation, the computational details of the TDDFT-as approach, and the procedure we follow to compare absorption spectra computed by two different methods; in Sec. IV, we discuss the results of the exponent optimization and the comparison of the TDDFT-as results with reference TDDFT and the approximated TDDFT + TB ones; finally, in Sec. V, conclusions and future perspectives are drawn.

The TDDFT Casida eigenvalue equation2,4,5 is (considering for simplicity, closed-shell systems and singlet excitations)

(AB)1/2(A+B)(AB)1/2ZI=ωI2ZI,
(1)

where ωI are excitation energies and

(AB)as,bt=ωasδabδst,
(2)
(A+B)as,bt=ωasδabδst+4(as|bt)+2(as|fxc|bt).
(3)

In Eqs. (2) and (3), a, b (s, t) denote occupied (unoccupied) KS orbitals with eigenvalues ϵ, the single-particle gaps are ωas = ϵsϵa, and fxc=δ2Excδρδρ+δ2Excδρδρ is the XC kernel. The integrals of the Coulomb kernel, (as|bt)=ϕarϕsrϕbrϕtr/rrd3rd3r, appearing in Eq. (3) can be simplified using the RI approach, where an orbital basis function product can be approximated as

χμrχνrPΔμνPχPr,
(4)

where μ, ν, … indicate orbital basis functions, P, Q indicate the auxiliary basis functions, and the coefficients ΔμνP minimizing the error in the Coulomb norm are given by56,59

ΔμνP=Q(P|Q)1(Q|μν).
(5)

Here and hereafter, we use the notation (a|b)=χarχbr/|rr|d3rd3r. Orbital (and auxiliary) basis functions are atom-centered functions of the type χμr=χAlmr=RAl(|rRA|)Ylm((rRA)/|rRA|), where Ylm are real spherical harmonics and RA is the position of atom A.

In the RI approach, the Coulomb integrals appearing in Eq. (3) become

(as|bt)μνκλCμaCνsCκbCλtPQΔμνP(P|Q)ΔκλQ
(6)
=PQΓasP(P|Q)ΓbtQ,
(7)

where Cμa are the coefficients of the KS orbital ϕa and

ΓasP=μνCμaCνsΔμνP
(8)

is the contribution of the KS orbital product ϕaϕs to the auxiliary basis function χP.

The eigenvectors of Eq. (1) are related to the transition density, i.e., the first-order change of the electronic density: for excitation I, the transition density is

δρIr=as2ωasωIZasIϕarϕsr
(9)
as2ωasωIZasIPΓasPχPr.
(10)

It is interesting to show that Eq. (7) can be recast to the approximation, which is used in TB schemes.62–69 In fact, if in Eq. (5) we do a Mulliken approximation [i.e., χμrχμ|χνχνr=Sμνχνr, and symmetrizing], we obtain

ΔμνP12Q(P|Q)1Sμν(Q|μμ)+(Q|νν)12SμνQ(P|Q)1(Q|Aμ)+(Q|Aν)12Sμν(δAμ,P+δAν,P),
(11)

where in the second line, we assume that each |χμr|2=|χAlmr|2 corresponds to the auxiliary basis function χAμ [i.e., by keeping only the monopolar term of |χμr|2 so that every orbital basis function on the atom A corresponds to the same auxiliary basis function χAμ]. Using Eq. (11) in Eq. (8), we obtain (using the atom index A to name the auxiliary basis function as there is just one auxiliary basis function per atom)

ΓasPΓ̃asA=12μA,νCμaCνsSμν+CνaCμsSνμ,
(12)

which equals the TD-DFTB Mulliken atomic transition charges.62 

It is important to observe that atomic transition charges in Eq. (12) are “population analysis,” i.e., they depend only on the overlap between atomic orbitals and do not depend at all on the auxiliary basis function. This is very different from Eq. (8), where ΔμνP clearly depends on the auxiliary function P, i.e., on the exponent of the Gaussian function, see Eq. (5). Recently, another population analysis (Löwdin) has been used in STDA and TDDFT + TB67,69 instead of the Mulliken population.

Finally, combining Eq. (12) with Eq. (7), we obtain

(as|bt)ABΓ̃asA(A|B)Γ̃btB,
(13)

which is the TD-DFTB approximation for the Coulomb kernel integrals.62–69 

Actually, Eq. (13) corresponds to the full TDDFT kernel, i.e., with Coulomb and XC effects, if the latter are treated with semilocal XC functionals. In fact, in this case, the last term in Eq. (3) can be approximated as

(as|fxc|bt)ABΓasA(A|fxc|B)ΓbtBAΓasAWAΓbtA,
(14)

where WA is an onsite parameter, which can be, thus, incorporated into Eq. (13). In fact, in the conventional TD-DFTB approaches and related methods, the two-center integrals (A|B) are usually not computed from the basis set but with an analytic function of the internuclear distance (RAB) and Hubbard parameters (U) of the atom types. For example, the sTDA/sTDDFT methods67,68 (for the same atom type, as in the present work) employ

(A|B)=γ(RAB)=1RABa+Ua1/a,
(15)

where a is a parameter. Thus, we have that U = γ(0) can be computed from isolated atoms and includes both the Coulomb and XC effects. Note also that for the metallic nanoparticles considered in this work, the role of the XC kernel is very small.74 

In summary, we have that the TB integrals in Eq. (13) employ three approximations:

  • a single and spherical function for each atom is used to expand the transition densities;

  • an atomic Hubbard parameter U is used to define the single and spherical function;

  • a population analysis (Mulliken or Löwdin) approximation is used to compute the atomic transition charges Γ̃asA.

On the other hand, no approximations (beside the size of the auxiliary basis set) are used for the general expression in Eqs. (5), (7), and (8).

Despite the quite good accuracy of the conventional TD-DFTB approaches and related methods,62–69 no detailed study of the role of approximations (i)–(ii)–(iii) is present in the literature.

In this work, we propose to keep only one Gaussian s-type orbital function per atom in the auxiliary basis set, i.e., approximation (i), but the Gaussian exponent is fitted to reproduce the absorption spectra of test systems, i.e., we do not use an atomic Hubbard parameter as in approximation (ii). Moreover, the atomic transition charges are computed exactly from Eqs. (5) and (8), without using the population approximation (iii). In this way, we can test directly the accuracy of using only approximation (i). We named such a method with one s-type function per atom in the auxiliary basis set as “TDDFT-as.”

In TDDFT-as, the two-center integral between two normalized s-type functions, i.e., gsr=(α/π)3/2exp(αr2) with gsrd3r=1, centered on atoms A and B is75,76

γ(RAB)=gs(rRA)gs(rRB)|rr|d3rd3r
(16)
=erf(RABα/2)RAB,
(17)

where RAB is the distance between atoms A and B. When RAB = 0, we have

γ(0)=2α/π=Uα,
(18)

which relates the exponent of the Gaussian s-type function to the Hubbard parameter. Note, however, that, in general, the auxiliary basis set is Coulomb normalized, i.e., gCsr=Uα1/2gsr=(23/4α5/4)/(2π5/4)exp(αr2). If α = 0, all the kernel contributions vanish and a single-particle absorption spectrum is obtained. In Sec. IV A, we will verify if an optimal α exists for silver NPs.

In this work, we considered eight silver NPs of different sizes and shapes: Ag20, Ag35+, and Ag120 (tetrahedron, Td symmetry); Ag32 (cubic rod, D4d symmetry); Ag551 and Ag555+ (cuboctahedron, Oh symmetry); and Ag553 and Ag555+ (icosahedron, Ih symmetry).

All the simulations including geometry relaxations and ground state and excited state calculations have been performed using the TURBOMOLE program package.77 For all these tasks, we used the Perdew–Burke–Ernzerhof (PBE)78 functional for the XC energy, the default effective core potential (ECP) with 28 core electrons,79 and the split valence plus polarization (def-SVP) basis set.80 This choice represents a good compromise between accuracy and computational efficiency.

For each system, we first performed a reference TDDFT calculation with the conventional auxiliary81 basis set, i.e., 7s4p3d3f2g, which corresponds to 97 Cartesian basis functions per atom. Then, we ran TDDFT-as calculations with just one s-type auxiliary basis function and neglecting the XC kernel.

All TDDFT calculations have also been done freezing the 4s and 4p occupied orbitals (i.e., with eigenvalues below 2 a.u.) and a number of virtual orbitals corresponding to the number of atoms (i.e., with eigenvalues above 50 a.u.). The number of considered occupied and virtual orbitals is reported in Table I. Freezing these orbitals has no effect on the spectra in the considered range (i.e., up to about 5 eV).

TABLE I.

Details of all systems considered. From left to right: symmetry, the number of occupied orbitals considered, the number of virtual orbitals considered, optical active space for TDDFT, highest energy considered in eV (ωmax), and the corresponding number of states computed (Nstates).

Syst.Symm.Occ.Virt.SpaceωmaxNstates
Ag20 Td 110 250 3 503 t2 5.30 100 t2 
Ag32 D4h 176 400 4 988 a2u4.59 66 a2u
    9 517 eu  121 eu 
Ag35+1 Td 192 438 10 651 t2 4.74 165 t2 
Ag551 Oh 303 687 13 180 t1u 4.72 186 t1u 
Ag55+5 Oh 300 690 13 112 t1u 5.06 221 t1u 
Ag553 Ih 304 686 5 368 t1u 4.92 101 t1u 
Ag555+ Ih 300 690 5 332 t1u 5.50 115 t1u 
Ag120 Td 660 1500 124 481 t2 4.18 1000 t1u 
Syst.Symm.Occ.Virt.SpaceωmaxNstates
Ag20 Td 110 250 3 503 t2 5.30 100 t2 
Ag32 D4h 176 400 4 988 a2u4.59 66 a2u
    9 517 eu  121 eu 
Ag35+1 Td 192 438 10 651 t2 4.74 165 t2 
Ag551 Oh 303 687 13 180 t1u 4.72 186 t1u 
Ag55+5 Oh 300 690 13 112 t1u 5.06 221 t1u 
Ag553 Ih 304 686 5 368 t1u 4.92 101 t1u 
Ag555+ Ih 300 690 5 332 t1u 5.50 115 t1u 
Ag120 Td 660 1500 124 481 t2 4.18 1000 t1u 

For the exponent optimization (see Sec. IV A), we considered all excited states up to the energy ωmax reported in Table I, with the corresponding number of (optically allowed) states considered. In this way, we considered not only the main plasmon peak but also structures in the high-energy region.

We also performed TDDFT + TB calculations using a locally modified version of the STDA code,82 in which the kernel integrals have been replaced by Eq. (17). The STDA code, in contrast to TURBOMOLE, does not use any symmetry: in order to reduce the computational cost, we have removed from the transition space all single-particle transitions, which are optically dark by symmetry.

To estimate the accuracy of TDDFT-as results with respect to the reference TDDFT ones, we considered the following root mean square (rms) error indicators for excitation energies and oscillator strengths (fi):

Eene[α]=i=1Nfiref(ωirefωj(i)[α])2i=1Nfiref1/2,
(19)
Eosc[α]=i=1N(fireffj(i)[α])2N1/2,
(20)

where N is Nstates in Table I and the function j(i) associates the excited state i of the reference TDDFT calculation with the excited state j in the TDDFT-as calculation. Assuming j(i) = i means that the order of excited states is the same in TDDFT-as and TDDFT, which is often not the case. The association is done considering the quantity

Oij=asZasi,refZasj,αexp|fireffj[α]|
(21)

and, for a given i, selecting the j-state with the largest Oij. Note that 0 < Oij < 1. In Eq. (21), the first part is the overlap between the eigenvectors of Eq. (1), whereas the second exponential factor is included to select states with closer oscillator strengths. In fact, in some cases, the overlap between different eigenvectors can be similar, and the second exponential factor helps to select the right state.

The indicators Eene[α] and Eosc[α] are very accurate, but also quite complicated and require that the two sets of excitation energies and oscillator strengths are close each other otherwise the association will fail (i.e., Oij will be always too small). A simpler indicator with a more general applicability can be defined as

Espe[α]=0ωmaxσref(ω)σ(ω)[α]2dω0ωmax|σref(ω)|2dω1/2,
(22)

where σref(ω) and σ(ω)[α] are the photoabsorption cross sections for TDDFT and TDDFT-as, respectively, and ωmax is defined in Table I.

The photoabsorption cross section is obtained from the computed excitation energies and oscillator strengths using a Lorentzian broadening, namely,

σ(ω)=2π2ciNfiπg(ωωi)2+g2.
(23)

The function Espe[α] does not include any direct information about the error on excitation energies; moreover, it depends on the broadening g. However, it is very simple to compute and, as we will show in this article, can give quite accurate information.

To optimize the exponent of the auxiliary basis set, we considered seven silver NPs: all the ones in Table I but Ag120, which has been considered only for timings, see Sec. IV D.

In Fig. 1(a), we report the error Eene[α] for all systems considered. The first interesting result is that all the curves have a similar shape with a single, quite flat minimum at around α = 0.03, more precisely in the range from α = 0.026 for Ag553(Ih) to α = 0.040 for Ag32(D4h). Although the values of optimal α change a little among different systems, most of the curves are quite flat, meaning that accurate results can be obtained within a significant range of α values.

FIG. 1.

(a) Error Eene for the TDDFT-as method as a function of exponent α for all systems. (b) Averaged Eene and averaged Eosc as a function of exponent α.

FIG. 1.

(a) Error Eene for the TDDFT-as method as a function of exponent α for all systems. (b) Averaged Eene and averaged Eosc as a function of exponent α.

Close modal

The second important result is that the value of Eene[α] at the minima is very small, in the range 7 meV–17 meV, meaning that TDDFT-as predicts all excitation energies with an average rms error less than 20 meV, which means a very high accuracy (less than the orbital basis set error).

Then we considered, see Fig. 1(b), the quantities ⟨Eene⟩ and ⟨Eosc⟩, which are the average among all systems. The ⟨Eene⟩ minimum is for α = 0.036 with an average rms error among all systems of only 12 meV.

For the oscillator strength, the ⟨Eosc⟩ curve does not show any well-defined minimum: it is very flat, with an average rms error of 0.25. Thus, we can consider the optimized value for ⟨Eene⟩,

α=0.036,
(24)

as a fixed exponent for all TDDFT-as calculations of silver NPs. An error distribution analysis (see Fig. S1 of the supplementary material) shows that all energies are underestimated by TDDFT-as. The results for the Espe indicator (see Fig. S2 of the supplementary material) give similar conclusion but with a smaller exponent (α = 0.034). Although the absolute value of Espe depends on the broadening parameter, it turned out that the optimized α is not.

The accuracy of TDDFT-as in reproducing the TDDFT reference spectra can also be directly recognized in Fig. 2, where three structures (the smallest one and the two extreme cases for the optimal α) have been reported. Figure 2 clearly shows a very high accuracy over the whole energy range. More in detail, we note that the intensity of the main plasmon peak is a bit reduced in TDDFT-as. We found that the intensity on the main peak in TDDFT-as can be increased reducing the α value; however, in this case, also the intensity of all other peaks at high-energy is increased, thus decreasing the overall accuracy.

FIG. 2.

Absorption cross section for the three systems studied: Ag20 (Td), Ag32 (D4h), and Ag553 (Ih). The curves end at ωmax of Table I. The black solid curves correspond to the standard TDDFT calculation, while the red-dashed ones to TDDFT-as with α = 0.036. The spectra are obtained by applying a Lorentzian broadening of 0.05 eV.

FIG. 2.

Absorption cross section for the three systems studied: Ag20 (Td), Ag32 (D4h), and Ag553 (Ih). The curves end at ωmax of Table I. The black solid curves correspond to the standard TDDFT calculation, while the red-dashed ones to TDDFT-as with α = 0.036. The spectra are obtained by applying a Lorentzian broadening of 0.05 eV.

Close modal

In Fig. 3, we report the reference transition densities, i.e., computed using KS orbitals, see Eq. (9), and the ones using the auxiliary basis set, see Eq. (10), considering just one optimized s-type function.

FIG. 3.

Transition densities at the main plasmon peak for (a) Ag20(Td), (b) Ag32(D4h), and (c) Ag553(Ih), as computed (left) from KS orbitals, i.e., using Eq. (9), or (right) using Eq. (10) with only one s-type function in the auxiliary basis set.

FIG. 3.

Transition densities at the main plasmon peak for (a) Ag20(Td), (b) Ag32(D4h), and (c) Ag553(Ih), as computed (left) from KS orbitals, i.e., using Eq. (9), or (right) using Eq. (10) with only one s-type function in the auxiliary basis set.

Close modal

The reference transition densities show a very rich structure due to atomistic details, and it is not very simple to distinguish the direction of the dipole. On the other hand, the TDDFT-as approximated transition densities are much smoother and without any atomistic structure, but both having similar (deviations about 10%, see Table S1 of the supplementary material) oscillator strength of the reference one. Thus, a much simpler interpretation of dipole moment and the charge redistribution in the excited state is obtained; moreover, the calculation of transition density plots is much cheaper.

However, we note that the transition dipole moments [μI=rδρIrd3r], which are used to compute the TDDFT-as absorption spectrum, are computed directly from orbitals, i.e., using Eq. (9), and not from Eq. (10). In fact, the auxiliary basis set is used in TDDFT-as only to approximate the Coulomb integrals in the TDDFT kernel matrix. An alternative approach is to use directly the expansion in Eq. (10) also for the computation of the transition dipole moments. This is, indeed, very efficient when only one s-type basis set is involved; in fact, in this case, the transition dipole moment is

μITDDFTas=as2ωasωIZasIAΓasARANA,
(25)

where the sum in A is over all atoms and

NA=gsrd3r=Uα1/2=1.119515α1/4.
(26)

Thus, the calculation of transition dipole moment between each pairs of atomic orbitals is not required anymore. A comparison of the oscillator strength for the main plasmon peaks for all systems is reported in the supplementary material (Table S1). We have also recomputed Eosc for all systems with this alternative definition of the transition dipole moment. Results show (see Fig. S3 of the supplementary material) that a well-defined minimum is present, in contrast to Fig. 1(b). In this case, in fact, there are stronger constraints (both Coulomb integrals and dipole moments) on the optimization of the α exponent. We found that the optimized value is α = 0.039, thus very close to the one obtained in Sec. IV A. However, the rms error on the oscillator strength increases by about 40% (⟨Eosc⟩ = 0.034), and the resulting spectra are not very accurate at high energies. Thus, while transition density plots at the plasmon peak can be safely used, Eq. (9) is more accurate for the calculation of the absorption spectrum, and it is used by default in TDDFT-as calculations.

The optimized α value for the TDDFT-as approach (α = 0.036) for the silver atom corresponds to an onsite Hubbard parameter U = 0.151 a.u., see Eq. (18). This value is much smaller than the one used in DFTB parameterization (U = 0.2191 a.u.) or in the sTDA and sTDDFT methods (U = 0.250 a.u.).

In Fig. 4, we report the absorption spectra computed with TDDFT + TB with different values of Hubbard parameters, for three selected systems Ag20(Td), Ag32(D4h), and Ag553(Ih).

FIG. 4.

Absorption spectra for the three systems studied, Ag20(Td), Ag32(D4h), and Ag551(Ih), from TDDFT + TB with different Hubbard parameters. The black curve is the TDDFT reference.

FIG. 4.

Absorption spectra for the three systems studied, Ag20(Td), Ag32(D4h), and Ag551(Ih), from TDDFT + TB with different Hubbard parameters. The black curve is the TDDFT reference.

Close modal

Recall that while in TDDFT-as atomic transition charges are computed from Eqs. (5) and (8), in TDDFT + TB atomic transition charges are computed with the Löwdin population approximation.69 For two-center integrals, we also use Eq. (17) in TDDFT + TB so that our implementation of TDDFT + TB is a little bit different from the original scheme.69 A modification of the two-center formula has some effects on the spectrum, but much less than the modification of the Hubbard parameter: recall, in fact, that all two-center formulas have γ(0) = U and γ(R) = 1/R for R large; thus, modifications can occur only for first neighboring atoms. In both TDDFT-as and TDDFT + TB, transitions dipole moments are computed directly from KS orbitals, i.e., computing the dipole moment of Eq. (9). Thus, the only difference between TDDFT-as and our TDDFT + TB is in the atomic transition charge computational scheme, besides the actual value of the Hubbard parameter (or α exponent). The top panel of Fig. 4 clearly shows that if the Hubbard parameter from our TDDFT-as approach is combined with a population analysis approximation for the atomic transition charges, very bad absorption spectra are obtained for Ag20 (Td). The TDDFT + TB spectrum using U = 0.151 a.u. is also very bad for the other two systems, and it is not reported for graphical reasons.

Using a larger Hubbard parameter, improved spectra are obtained for Ag20(Td). Yet, the default TDDFT + TB still overestimates the intensity of the main plasmon peak. A better accuracy is obtained using the optimized value (U = 0.337 a.u.), obtained minimizing Espe, which improved the description of the plasmon peak. Nevertheless, the intensities of high-energy part (>4 eV) of the spectrum are overestimated.

For Ag553(Ih), bottom panel of Fig. 4, the situation is similar. In this case, the optimized Hubbard parameter is even larger (U = 0.4 a.u.), and the agreement with the reference TDDFT is better. Instead, the behavior of Ag32(D4h) is very different. With default TDDFT + TB, very bad absorption spectra are obtained, with a big peak at around 3.8 eV, not present in the reference TDDFT. In this case, increasing the Hubbard parameter does not help because in this way, the main plasmon peak at 2.9 eV disappears.

In Ref. 72, TDDFT and TDDFT + TB absorption spectra for silver and gold clusters have been reported. Silver clusters are similar to the ones of this work, but the TDDFT reference employs a different XC functional83 and TDDFT + TB a different formula for two-center integrals. In any case, results of Ref. 72 indicate a similar accuracy and overestimation of the oscillator strengths of TDDFT + TB with respect to TDDFT at high energy, which is consistent with our findings.

Comparing Fig. 4 with Fig. 2, the much higher accuracy of TDDFT-as with respect to TDDFT + TB is, thus, evident. A more quantitative analysis is reported in Fig. 5 for all systems. The default TDDFT + TB yields an error for Espe in the range 0.70–1.08. Optimizing the Hubbard parameter for each system, the Espe can be reduced down to 0.25–0.69. Clearly, these values of Espe are much larger than the ones that can be obtained using the TDDFT-as with α = 0.036, with an error of Espe in the range 0.09–0.20.

FIG. 5.

Error on the spectrum Espe for all systems as obtained from TDDFT + TB with different Hubbard parameters and TDDFT-as. In TDDFT + TB (U opt.), the Hubbard parameter is optimized for each system.

FIG. 5.

Error on the spectrum Espe for all systems as obtained from TDDFT + TB with different Hubbard parameters and TDDFT-as. In TDDFT + TB (U opt.), the Hubbard parameter is optimized for each system.

Close modal

This direct comparison between TDDFT + TB and TDDFT-as clearly indicates that the origin of the reduced accuracy of TDDFT + TB is related to the population analysis approximation (as all other quantities are the same in the two methods).

In Fig. 6, we compared the atomic transition charges computed with our TDDFT-as method, i.e., using Eq. (8), with the ones from the TDDFT + TB approach, for two different single-particle transitions of Ag20. The plots clearly show that the TDDFT + TB atomic transition charges are much smaller than the TDDFT-as ones; thus, a larger Hubbard parameter is required in the former to obtain similar Coulomb integrals. This finding explains the results in Fig. 4 where an incorrect spectrum is obtained with U = 0.151 a.u.

FIG. 6.

Atomic transition charges for all atoms (ordered according to their position in the z axis) for the Ag20 system, for two single-particle transitions. The black curves are results from TDDFT-as [i.e., Eq. (8)] and the red curves from TDDFT + TB (i.e., Löwdin population analysis).

FIG. 6.

Atomic transition charges for all atoms (ordered according to their position in the z axis) for the Ag20 system, for two single-particle transitions. The black curves are results from TDDFT-as [i.e., Eq. (8)] and the red curves from TDDFT + TB (i.e., Löwdin population analysis).

Close modal

Finally, we considered a larger nanoparticle, Ag120 (Td), not included in the fitting procedure. The computed photoabsorption cross sections from TDDFT and TDDFT-as are reported in Fig. 7. We used a log-scale for the y axis and a small broadening (g = 20 meV) so that the accuracy for all peaks can be well distinguished. Clearly, the agreement with the reference TDDFT calculation is excellent as the difference between the two curves can be hardly distinguished, for all the energy range considered.

FIG. 7.

Photoabsorption cross section for Ag120(Td) using TDDFT and TDDFT-as. Note the log-scale on the y axis. The broadening is g = 20 meV.

FIG. 7.

Photoabsorption cross section for Ag120(Td) using TDDFT and TDDFT-as. Note the log-scale on the y axis. The broadening is g = 20 meV.

Close modal

We also considered the timings for the Ag120(Td) nanoparticle. We report in Table II the time (in seconds) for the Coulomb, XC, and RI part, i.e., the calculation in memory of the three-center integrals (Q|μν), see Eq. (5). The TDDFT-as computational cost for all these steps is less than 0.5% of the full TDDFT results, which is consistent with the reduction in the auxiliary basis size of two orders of magnitude and the reduced number of Davidson iteration steps required.

TABLE II.

Computational cost for the Ag120(Td) nanoparticle with different options for the number of states and symmetries considered. From left to right: the number of steps of the Davidson diagonalization routine, seconds for the calculation of the Coulomb integrals for the TDDFT kernel, seconds for the calculation of the XC integrals for the TDDFT kernel, seconds for the calculation of all three-center integrals. Timings have been performed on a single core of Xeon Gold 6132 CPU 2.6 GHz.

States/Symm.MethodStepsCoul.XCRI
17/C1 TDDFT 10 11 463 13 836 11 576 
17/C1 TDDFT-as 53 87 
100/Td TDDFT 5 945 7 986 2 043 
100/Td TDDFT-as 19 21 
400/Td TDDFT 23 320 25 989 23 323 
400/Td TDDFT-as 75 77 
States/Symm.MethodStepsCoul.XCRI
17/C1 TDDFT 10 11 463 13 836 11 576 
17/C1 TDDFT-as 53 87 
100/Td TDDFT 5 945 7 986 2 043 
100/Td TDDFT-as 19 21 
400/Td TDDFT 23 320 25 989 23 323 
400/Td TDDFT-as 75 77 

However, we point out that we did not optimize the code for efficiency. In fact, the main focus of the present work was to verify the accuracy of a single s-type function in the auxiliary basis set.

In this article, we have shown that it is possible to reproduce TDDFT results for silver nanoparticles with an accuracy of 12 meV using only one s-type function per atom in the auxiliary basis set. In addition, the exponent α of this s-type function is optimized to model XC effects too. In this approach, named TDDFT-as, the computational cost for the calculation of the TDDFT kernel matrix elements is reduced by two orders of magnitude, with respect to TDDFT with a default auxiliary basis set. The TDDFT-as method is based on two observations: (i) to reproduce the absorption spectrum of silver NPs, there is no need to reproduce exactly all the atomistic details of the transition density: they can be averaged out for the calculation of Coulomb integrals and also for the dipole moment and (ii) the effects of a (semi-)local XC kernel are much smaller than the Coulomb one;74 thus, the former can be included as an onsite correction of the latter.

In this work, we have optimized the exponent among seven silver NPs with different sizes and symmetries, showing the transferability of the α exponent. The TDDFT-as approach accurately reproduces the main plasmon peak as well as the high-energy part of the absorption spectrum.

The TDDFT-as approach resembles the TDDFT + TB scheme: both are based on the KS eigenvalues and eigenvectors and use a monopolar approximation to compute the Coulomb-XC kernel matrix elements. TDDFT + TB employs two other additional approximations, i.e., the population analysis to compute the atomic transition charges and a Hubbard atomic parameter. In this work, we have shown that TDDFT-as is by far more accurate than any TDDFT + TB approach with a population analysis approximation and optimized Hubbard parameter. For some systems, e.g., a cubic silver rod, TDDFT + TB even predicts a spurious peak in the high-energy part of the spectrum.

The TDDFT-as approach is very simple, and it is easily available in all RI-TDDFT codes (changing the auxiliary basis set and switching-off XC contributions).

There are many directions to be followed in the future:

In this work, we have, in fact, only considered weakly charged systems and the PBE XC functional. It will be interesting to study if the optimal α will change with the KS potential, which can be largely modified using a different XC functional (i.e., the LB94 functional83) or by a large system charge. In this case, an alternative approach could be used to extract a system-dependent optimal α, using information from the ground-state KS orbitals.

Another path to follow is to further optimize the auxiliary basis set in order to verify if it is possible to obtain even better accuracy, without increasing the computational cost. It will be interesting, for example, to use a contraction of more s-type functions or different s-type functions for different atoms (e.g., on the surface or inside the nanoparticle).

Finally, TDDFT-as can be extended to other systems, e.g., other metals, metal alloys, or even molecular systems. We expect that TDDFT-as can be quite accurate for all these systems as it is closer to TDDFT than TD-DFTB and TDDFT + TB. Preliminary calculations show that TDDFT-as is quite accurate in reproducing the absorption spectra of gold NPs. However, the overall accuracy is lower for silver NPs; in fact, gold NPs have a more complex absorption spectrum with stronger contributions from the d-band,84 which cannot be fully described with just one s-type function. When systems with different atom types are considered, it must be kept in mind that the optimized α for a given atom might also depend on the local environment as it aims to describe transition densities, which can be quite delocalized. Extensions of the TDDFT-as method (e.g., including two basis functions per atom) are also possible.

See the supplementary material for TDDFT-as error distribution for excitation energies and oscillator strengths, excitation energies and oscillator strengths at the plasmon peak for all systems, error Espe for all systems and averaged value, and error Eosc for all systems and averaged value using transition dipole moments from the auxiliary basis set.

The data that support the findings of this study are available from the corresponding author upon reasonable request.

1.
E.
Runge
and
E. K. U.
Gross
, “
Density-functional theory for time-dependent systems
,”
Phys. Rev. Lett.
52
,
997
(
1984
).
2.
M. E.
Casida
, “
Time-dependent density functional response theory for molecules
,” in
Recent Advances in Density Functional Methods
, edited by
D. P.
Chong
(
World Scientific
,
Singapore
,
1995
), Vol. 1, pp.
155
192
.
3.
M. A. L.
Marques
and
E. K. U.
Gross
, “
Time-dependent density functional theory
,”
Ann. Rev. Phys. Chem.
55
,
427
(
2004
).
4.
M. E.
Casida
and
M.
Huix-Rotllant
, “
Progress in time-dependent density-functional theory
,”
Annu. Rev. Phys. Chem.
63
,
287
(
2012
).
5.
C. A.
Ullrich
,
Time-Dependent Density-Functional Theory: Concepts and Applications
, Oxford Graduate Texts (
Oxford University Press
,
Oxford, New York
,
2012
).
6.
A. D.
Laurent
and
D.
Jacquemin
, “
TD-DFT benchmarks: A review
,”
Int. J. Quantum Chem.
113
,
2019
(
2013
).
7.
S. M.
Morton
,
D. W.
Silverstein
, and
L.
Jensen
, “
Theoretical studies of plasmonics using electronic structure methods
,”
Chem. Rev.
111
,
3962
(
2011
).
8.
P.
Zhang
,
J.
Feist
,
A.
Rubio
,
P.
García-González
, and
F. J.
García-Vidal
, “
Ab initio nanoplasmonics: The impact of atomic structure
,”
Phys. Rev. B
90
,
161407
(
2014
).
9.
A.
Varas
,
P.
García-González
,
F. J.
García-Vidal
, and
A.
Rubio
, “
Anisotropy effects on the plasmonic response of nanoparticle dimers
,”
J. Phys. Chem. Lett.
6
,
1891
(
2015
).
10.
M.
Barbry
,
P.
Koval
,
F.
Marchesin
,
R.
Esteban
,
A. G.
Borisov
,
J.
Aizpurua
, and
D.
Sánchez-Portal
, “
Atomistic near-field nanoplasmonics: Reaching atomic-scale resolution in nanooptics
,”
Nano Lett.
15
,
3410
(
2015
).
11.
A.
Varas
,
P.
García-González
,
J.
Feist
,
F. J.
García-Vidal
, and
A.
Rubio
, “
Quantum plasmonics: From jellium models to ab initio calculations
,”
Nanophotonics
5
,
409
(
2016
).
12.
W.
Zhu
,
R.
Esteban
,
A. G.
Borisov
,
J. J.
Baumberg
,
P.
Nordlander
,
H. J.
Lezec
,
J.
Aizpurua
, and
K. B.
Crozier
, “
Quantum mechanical effects in plasmonic structures with subnanometre gaps
,”
Nat. Commun.
7
,
11495
(
2016
).
13.
P.
Zhang
,
W.
Jin
, and
W.
Liang
, “
Size-dependent optical properties of aluminum nanoparticles: From classical to quantum description
,”
J. Phys. Chem. C
122
,
10545
(
2018
).
14.
F.
Della Sala
,
M.
Pezzolla
,
S.
D’Agostino
, and
E.
Fabiano
, “
Ab initio plasmonics of externally doped silicon nanocrystals
,”
ACS Photonics
6
,
1474
(
2019
).
15.
S. A.
Maier
,
Plasmonics, Fundamentals and Applications
(
Springer
,
2007
).
16.
Modern Plasmonics
, Handbook of Surface Science, edited by
A. A.
Maradudin
,
J. R.
Sambles
, and
W. L.
Barnes
(
Elsevier
,
Amsterdam
,
2014
), Vol. 4.
17.
H. A.
Atwater
and
A.
Polman
, “
Plasmonics for improved photovoltaic devices
,”
Nat. Mater.
9
,
205
213
(
2010
).
18.
S.
Nie
and
S. R.
Emory
, “
Probing single molecules and single nanoparticles by surface-enhanced Raman scattering
,”
Science
275
,
1102
(
1997
).
19.
S.
Kawata
,
Y.
Inouye
, and
P.
Verma
, “
Plasmonics for near-field nano-imaging and superlensing
,”
Nat. Photonics
3
,
388
(
2009
).
20.
M.
Quinten
,
Optical Properties of Nanoparticle Systems: Mie and beyond
(
Wiley-VCH Verlag GmbH & Co. KGaA
,
Weinheim, Germany
,
2011
).
21.
C. F.
Bohren
and
D. R.
Huffman
,
Absorption and Scattering of Light by Small Particles
(
Wiley VCH
,
Weinheim
,
2004
).
22.
Y.
Luo
,
A. I.
Fernandez-Dominguez
,
A.
Wiener
,
S. A.
Maier
, and
J. B.
Pendry
, “
Surface plasmons and nonlocality: A simple model
,”
Phys. Rev. Lett.
111
,
093901
(
2013
).
23.
T.
Christensen
,
W.
Yan
,
A.-P.
Jauho
,
M.
Soljačić
, and
N. A.
Mortensen
, “
Quantum corrections in nanoplasmonics: Shape, scale, and material
,”
Phys. Rev. Lett.
118
,
157402
(
2017
).
24.
W.
Yan
,
M.
Wubs
, and
N.
Asger Mortensen
, “
Projected dipole model for quantum plasmonics
,”
Phys. Rev. Lett.
115
,
137403
(
2015
).
25.
R.
Esteban
,
A. G.
Borisov
,
P.
Nordlander
, and
J.
Aizpurua
, “
Bridging quantum and classical plasmonics with a quantum-corrected model
,”
Nat. Commun.
3
,
825
(
2012
).
26.
G.
Toscano
,
J.
Straubel
,
A.
Kwiatkowski
,
C.
Rockstuhl
,
F.
Evers
,
H.
Xu
,
N.
Asger Mortensen
, and
M.
Wubs
, “
Resonance shifts and spill-out effects in self-consistent hydrodynamic nanoplasmonics
,”
Nat. Commun.
6
,
7132
(
2015
).
27.
C.
Ciracì
and
F.
Della Sala
, “
Quantum hydrodynamic theory for plasmonics: Impact of the electron density tail
,”
Phys. Rev. B
93
,
205405
(
2016
).
28.
C. M.
Aikens
,
S.
Li
, and
G. C.
Schatz
, “
From discrete electronic states to plasmons: TDDFT optical absorption properties of agn (n = 10, 20, 35, 56, 84, 120) tetrahedral clusters
,”
J. Phys. Chem. C
112
,
11272
(
2008
).
29.
H.-C.
Weissker
and
C.
Mottet
, “
Optical properties of pure and core-shell noble-metal nanoclusters from TDDFT: The influence of the atomic structure
,”
Phys. Rev. B
84
,
165443
(
2011
).
30.
G.-T.
Bae
and
C. M.
Aikens
, “
Time-dependent density functional theory studies of optical properties of Ag nanoparticles: Octahedra, truncated octahedra, and icosahedra
,”
J. Phys. Chem. C
116
,
10356
(
2012
).
31.
X.
López-Lozano
,
H.
Barron
,
C.
Mottet
, and
H.-C.
Weissker
, “
Aspect-ratio- and size-dependent emergence of the surface-plasmon resonance in gold nanorods—An ab initio TDDFT study
,”
Phys. Chem. Chem. Phys.
16
,
1820
(
2014
).
32.
G.
Barcaro
,
L.
Sementa
,
A.
Fortunelli
, and
M.
Stener
, “
Optical properties of silver nanoshells from time-dependent density functional theory calculations
,”
J. Phys. Chem. C
118
,
12450
(
2014
).
33.
M.
Kuisma
,
A.
Sakko
,
T. P.
Rossi
,
A. H.
Larsen
,
J.
Enkovaara
,
L.
Lehtovaara
, and
T. T.
Rantala
, “
Localized surface plasmon resonance in silver nanoparticles: Atomistic first-principles time-dependent density-functional theory calculations
,”
Phys. Rev. B
91
,
115431
(
2015
).
34.
H.-C.
Weissker
and
X.
López-Lozano
, “
Surface plasmons in quantum-sized noble-metal clusters: TDDFT quantum calculations and the classical picture of charge oscillations
,”
Phys. Chem. Chem. Phys.
17
,
28379
(
2015
).
35.
O.
Baseggio
,
M.
De Vetta
,
G.
Fronzoni
,
M.
Stener
,
L.
Sementa
,
A.
Fortunelli
, and
A.
Calzolari
, “
Photoabsorption of icosahedral noble metal clusters: An efficient TDDFT approach to large-scale systems
,”
J. Phys. Chem. C
120
,
12773
(
2016
).
36.
P.
Koval
,
F.
Marchesin
,
D.
Foerster
, and
D.
Sánchez-Portal
, “
Optical response of silver clusters and their hollow shells from linear-response TDDFT
,”
J. Phys.: Condens. Matter
28
,
214001
(
2016
).
37.
T. P.
Rossi
,
M.
Kuisma
,
M. J.
Puska
,
R. M.
Nieminen
, and
P.
Erhart
, “
Kohn–Sham decomposition in real-time time-dependent density-functional theory: An efficient tool for analyzing plasmonic excitations
,”
J. Chem. Theory Comput.
13
,
4779
(
2017
).
38.
R.
Schira
and
F.
Rabilloud
, “
Localized surface plasmon resonance in free silver nanoclusters Agn, n = 20–147
,”
J. Phys. Chem. C
123
,
6205
(
2019
).
39.
R.
Sinha-Roy
,
P.
García-González
, and
H.-C.
Weissker
, “
How metallic are noble-metal clusters? Static screening and polarizability in quantum-sized silver and gold nanoparticles
,”
Nanoscale
12
,
4452
(
2020
).
40.
G.
Piccini
,
R. W. A.
Havenith
,
R.
Broer
, and
M.
Stener
, “
Gold nanowires: A time-dependent density functional assessment of plasmonic behavior
,”
J. Phys. Chem. C
117
,
17196
(
2013
).
41.
S.
Malola
,
L.
Lehtovaara
,
J.
Enkovaara
, and
H.
Häkkinen
, “
Birth of the localized surface plasmon resonance in monolayer-protected gold nanoclusters
,”
ACS Nano
7
,
10263
(
2013
).
42.
K.
Iida
,
M.
Noda
,
K.
Ishimura
, and
K.
Nobusada
, “
First-principles computational visualization of localized surface plasmon resonance in gold nanoclusters
,”
J. Phys. Chem. A
118
,
11317
(
2014
).
43.
R. W.
Burgess
and
V. J.
Keast
, “
TDDFT study of the optical absorption spectra of bare gold clusters
,”
J. Phys. Chem. C
118
,
3194
(
2014
).
44.
K.
Iida
,
M.
Noda
, and
K.
Nobusada
, “
Interface electronic properties between a gold core and thiolate ligands: Effects on an optical absorption spectrum in Au133(SPh-TBu)52
,”
J. Phys. Chem. C
120
,
2753
(
2016
).
45.
S.
Hernandez
,
Y.
Xia
,
V.
Vlček
,
R.
Boutelle
,
R.
Baer
,
E.
Rabani
, and
D.
Neuhauser
, “
First-principles spectra of au nanoparticles: From quantum to classical absorption
,”
Mol. Phys.
116
,
2506
(
2018
).
46.
D.
Rocca
,
R.
Gebauer
,
Y.
Saad
, and
S.
Baroni
, “
Turbo charging time-dependent density-functional theory with Lanczos chains
,”
J. Chem. Phys.
128
,
154105
(
2008
).
47.
T. J.
Zuehlsdorff
,
N. D. M.
Hine
,
J. S.
Spencer
,
N. M.
Harrison
,
D. J.
Riley
, and
P. D.
Haynes
, “
Linear-scaling time-dependent density-functional theory in the linear response formalism
,”
J. Chem. Phys.
139
,
064104
(
2013
).
48.
O.
Baseggio
,
G.
Fronzoni
, and
M.
Stener
, “
A new time dependent density functional algorithm for large systems and plasmons in metal clusters
,”
J. Chem. Phys.
143
,
024106
(
2015
).
49.
J.
Brabec
,
L.
Lin
,
M.
Shao
,
N.
Govind
,
C.
Yang
,
Y.
Saad
, and
E. G.
Ng
, “
Efficient algorithms for estimating the absorption spectrum within linear response TDDFT
,”
J. Chem. Theory Comput.
11
,
5197
(
2015
).
50.
Y.
Gao
,
D.
Neuhauser
,
R.
Baer
, and
E.
Rabani
, “
Sublinear scaling for time-dependent stochastic density functional theory
,”
J. Chem. Phys.
142
,
034106
(
2015
).
51.
K.
Yabana
and
G. F.
Bertsch
, “
Time-dependent local-density approximation in real time
,”
Phys. Rev. B
54
,
4484
(
1996
).
52.
K. L. D. M.
Weerawardene
and
C. M.
Aikens
, “
Comparison and convergence of optical absorption spectra of noble metal nanoparticles computed using linear-response and real-time time-dependent density functional theories
,”
Comput. Theor. Chem.
1146
,
27
(
2018
).
53.
E. J.
Baerends
,
D. E.
Ellis
, and
P.
Ros
, “
Self-consistent molecular Hartree–Fock–Slater calculations. I. The computational procedure
,”
Chem. Phys.
2
,
41
(
1973
).
54.
B. I.
Dunlap
,
J. W. D.
Connolly
, and
J. R.
Sabin
, “
On some approximations in applications of Xα theory
,”
J. Chem. Phys.
71
,
3396
(
1979
).
55.
K.
Eichkorn
,
O.
Treutler
,
H.
Öhm
,
M.
Häser
, and
R.
Ahlrichs
, “
Auxiliary basis sets to approximate Coulomb potentials
,”
Chem. Phys. Lett.
240
,
283
(
1995
).
56.
R.
Bauernschmitt
,
M.
Häser
,
O.
Treutler
, and
R.
Ahlrichs
, “
Calculation of excitation energies within time-dependent density functional theory using auxiliary basis set expansions
,”
Chem. Phys. Lett.
264
,
573
(
1997
).
57.
H. H.
Heinze
,
A.
Görling
, and
N.
Rösch
, “
An efficient method for calculating molecular excitation energies by time-dependent density-functional theory
,”
J. Chem. Phys.
113
,
2088
(
2000
).
58.
F.
Neese
and
G.
Olbrich
, “
Efficient use of the resolution of the identity approximation in time-dependent density functional calculations with hybrid density functionals
,”
Chem. Phys. Lett.
362
,
170
(
2002
).
59.
T. B.
Pedersen
,
F.
Aquilante
, and
R.
Lindh
, “
Density fitting with auxiliary basis sets from Cholesky decompositions
,”
Theor. Chem. Acc.
124
,
1
(
2009
).
60.
F.
Weigend
,
M.
Kattannek
, and
R.
Ahlrichs
, “
Approximated electron repulsion integrals: Cholesky decomposition versus resolution of the identity methods
,”
J. Chem. Phys.
130
,
164106
(
2009
).
61.
G. L.
Stoychev
,
A. A.
Auer
, and
F.
Neese
, “
Automatic generation of auxiliary basis sets
,”
J. Chem. Theory Comput.
13
,
554
(
2017
).
62.
T. A.
Niehaus
,
S.
Suhai
,
F.
Della Sala
,
P.
Lugli
,
M.
Elstner
,
G.
Seifert
, and
T.
Frauenheim
, “
Tight-binding approach to time-dependent density-functional response theory
,”
Phys. Rev. B
63
,
085108
(
2001
).
63.
T. A.
Niehaus
, “
Approximate time-dependent density functional theory
,”
J. Mol. Struct.
914
,
38
(
2009
).
64.
F.
Trani
,
G.
Scalmani
,
G.
Zheng
,
I.
Carnimeo
,
M. J.
Frisch
, and
V.
Barone
, “
Time-dependent density functional tight binding: New formulation and benchmark of excited states
,”
J. Chem. Theory Comput.
7
,
3304
(
2011
).
65.
A.
Humeniuk
and
R.
Mitrić
, “
Long-range correction for tight-binding TD-DFT
,”
J. Chem. Phys.
143
,
134120
(
2015
).
66.
A.
Domínguez
,
B.
Aradi
,
T.
Frauenheim
,
V.
Lutsker
, and
T. A.
Niehaus
, “
Extensions of the time-dependent density functional based tight-binding approach
,”
J. Chem. Theory Comput.
9
,
4901
(
2013
).
67.
S.
Grimme
, “
A simplified Tamm-Dancoff density functional approach for the electronic excitation spectra of very large molecules
,”
J. Chem. Phys.
138
,
244104
(
2013
).
68.
C.
Bannwarth
and
S.
Grimme
, “
A simplified time-dependent density functional theory approach for electronic ultraviolet and circular dichroism spectra of very large molecules
,”
Comput. Theor. Chem.
1040-1041
,
45
(
2014
).
69.
R.
Rüger
,
E.
van Lenthe
,
T.
Heine
, and
L.
Visscher
, “
Tight-binding approximations to time-dependent density functional theory—A fast approach for the calculation of electronically excited states
,”
J. Chem. Phys.
144
,
184103
(
2016
).
70.
S.
D’Agostino
,
R.
Rinaldi
,
G.
Cuniberti
, and
F.
Della Sala
, “
Density functional tight binding for quantum plasmonics
,”
J. Phys. Chem. C
122
,
19756
(
2018
).
71.
F.
Alkan
and
C. M.
Aikens
, “
TD-DFT and TD-DFTB investigation of the optical properties and electronic structure of silver nanorods and nanorod dimers
,”
J. Phys. Chem. C
122
,
23639
(
2018
).
72.
N.
Asadi-Aghbolaghi
,
R.
Rüger
,
Z.
Jamshidi
, and
L.
Visscher
, “
TD-DFT + TB: An efficient and fast approach for quantum plasmonic excitations
,”
J. Phys. Chem. C
124
,
7946
(
2020
).
73.
F.
Weigend
, “
Accurate coulomb-fitting basis sets for H to Rn
,”
Phys. Chem. Chem. Phys.
8
,
1057
(
2006
).
74.
J. C.
Idrobo
and
S. T.
Pantelides
, “
Origin of bulklike optical response in noble-metal Ag and Au nanoparticles
,”
Phys. Rev. B
82
,
085420
(
2010
).
75.
P.
Koskinen
and
V.
Mäkinen
, “
Density-functional tight-binding for beginners
,”
Comput. Mater. Sci.
47
,
237
(
2009
).
76.
I. V.
Bodrenko
,
M.
Sierka
,
E.
Fabiano
, and
F. D.
Sala
, “
A periodic charge-dipole electrostatic model: Parametrization for silver slabs
,”
J. Chem. Phys.
137
,
134702
(
2012
).
77.
TURBOMOLE V7.2, 2009, a development of University of Karlsruhe and Forschungszentrum Karlsruhe GmbH, 1989–2007, TURBOMOLE GmbH, since 2007; available from http://www.turbomole.com.
78.
J. P.
Perdew
,
K.
Burke
, and
M.
Ernzerhof
, “
Generalized gradient approximation made simple
,”
Phys. Rev. Lett.
77
,
3865
(
1996
).
79.
D.
Andrae
,
U.
Hääußermann
,
M.
Dolg
,
H.
Stoll
, and
H.
Preuß
, “
Energy-adjusted ab initio pseudopotentials for the second and third row transition elements
,”
Theor. Chim. Acta
77
,
123
(
1990
).
80.
A.
Schäafer
,
H.
Horn
, and
R.
Ahlrichs
, “
Fully optimized contracted Gaussian basis sets for atoms Li to Kr
,”
J. Chem. Phys.
97
,
2571
(
1992
).
81.
K.
Eichkorn
,
F.
Weigend
,
O.
Treutler
, and
R.
Ahlrichs
, “
Auxiliary basis sets for main row atoms and transition metals and their use to approximate Coulomb potentials
,”
Theor. Chem. Acc.
97
,
119
(
1997
).
82.
The sTDA code can be downloaded from https://github.com/grimme-lab/stda.
83.
R.
van Leeuwen
and
E. J.
Baerends
, “
Exchange-correlation potential with correct asymptotic behavior
,”
Phys. Rev. A
49
,
2421
(
1994
).
84.
G.-T.
Bae
and
C. M.
Aikens
, “
Time-dependent density functional theory studies of optical properties of Au nanoparticles: Octahedra, truncated octahedra, and icosahedra
,”
J. Phys. Chem. C
119
,
23127
(
2015
).

Supplementary Material