Q-Chem 4.4 User’s Manual

6.3 Time-Dependent Density Functional Theory (TDDFT)

6.3.1 Brief Introduction to TDDFT

Excited states may be obtained from density functional theory by time-dependent density functional theory [367, 368], which calculates poles in the response of the ground state density to a time-varying applied electric field. These poles are Bohr frequencies or excitation energies, and are available in Q-Chem [369], together with the CIS-like Tamm-Dancoff approximation [188]. TDDFT is becoming very popular as a method for studying excited states because the computational cost is roughly similar to the simple CIS method (scaling as roughly the square of molecular size), but a description of differential electron correlation effects is implicit in the method. The excitation energies for low-lying valence excited states of molecules (below the ionization threshold, or more conservatively, below the first Rydberg threshold) are often remarkably improved relative to CIS, with an accuracy of roughly 0.1–0.3 eV being observed with either gradient corrected or local density functionals.

However, standard density functionals do not yield a potential with the correct long-range Coulomb tail, owing to the so-called self-interaction problem, and therefore excitation energies corresponding to states that sample this tail (e.g., diffuse Rydberg states and some charge transfer excited states) are not given accurately [187, 370, 156]. The extent to which a particular excited state is characterized by charge transfer can be assessed using an a spatial overlap metric proposed by Peach, Benfield, Helgaker, and Tozer (PBHT) [371]. (However, see Ref. Richard:2011 for a cautionary note regarding this metric.)

It is advisable to only employ TDDFT for low-lying valence excited states that are below the first ionization potential of the molecule. This makes radical cations a particularly favorable choice of system, as exploited in Ref. Hirata:2001. TDDFT for low-lying valence excited states of radicals is in general a remarkable improvement relative to CIS, including some states, that, when treated by wave function-based methods can involve a significant fraction of double excitation character [369]. The calculation of the nuclear gradients of full TDDFT and within the Tamm-Dancoff approximation is also implemented [373].

Standard TDDFT also does not yield a good description of static correlation effects (see Section 5.9), because it is based on a single reference configuration of Kohn-Sham orbitals. Recently, a new variation of TDDFT called spin-flip density functional theory (SFDFT) was developed by Yihan Shao, Martin Head-Gordon and Anna Krylov to address this issue [94]. SFDFT is different from standard TDDFT in two ways:

SF-DFT can describe the ground state as well as a few low-lying excited states, and has been applied to bond-breaking processes, and di- and tri-radicals with degenerate or near-degenerate frontier orbitals. Recently, we also implemented[374] a SFDFT method with a non-collinear exchange-correlation potential from Tom Ziegler et al. [375, 376], which is in many case an improvement over collinear SFDFT [94]. Recommended functionals for SF-DFT calculations are 5050 and PBE50 (see Ref. [374] for extensive benchmarks). See also Section 6.7.3 for details on wave function-based spin-flip models.

6.3.2 TDDFT within a Reduced Single-Excitation Space

Much of chemistry and biology occurs in solution or on surfaces. The molecular environment can have a large effect on electronic structure and may change chemical behavior. Q-Chem is able to compute excited states within a local region of a system through performing the TDDFT (or CIS) calculation with a reduced single excitation subspace [377]. This allows the excited states of a solute molecule to be studied with a large number of solvent molecules reducing the rapid rise in computational cost. The success of this approach relies on there being only weak mixing between the electronic excitations of interest and those omitted from the single excitation space. For systems in which there are strong hydrogen bonds between solute and solvent, it is advisable to include excitations associated with the neighboring solvent molecule(s) within the reduced excitation space.

The reduced single excitation space is constructed from excitations between a subset of occupied and virtual orbitals. These can be selected from an analysis based on Mulliken populations and molecular orbital coefficients. For this approach the atoms that constitute the solvent needs to be defined. Alternatively, the orbitals can be defined directly. The atoms or orbitals are specified within a $solute block. These approach is implemented within the TDA and has been used to study the excited states of formamide in solution [378], CO on the Pt(111) surface [379], and the tryptophan chromophore within proteins [380].

6.3.3 Job Control for TDDFT

Input for time-dependent density functional theory calculations follows very closely the input already described for the uncorrelated excited state methods described in the previous section (in particular, see Section 6.2.7). There are several points to be aware of:

TRNSS

Controls whether reduced single excitation space is used


TYPE:

LOGICAL


DEFAULT:

FALSE

Use full excitation space


OPTIONS:

TRUE

Use reduced excitation space


RECOMMENDATION:

None


TRTYPE

Controls how reduced subspace is specified


TYPE:

INTEGER


DEFAULT:

1


OPTIONS:

1

Select orbitals localized on a set of atoms

2

Specify a set of orbitals

3

Specify a set of occupied orbitals, include excitations to all virtual orbitals


RECOMMENDATION:

None


N_SOL

Specifies number of atoms or orbitals in $\$ solute$


TYPE:

INTEGER


DEFAULT:

No default


OPTIONS:

User defined


RECOMMENDATION:

None


CISTR_PRINT

Controls level of output


TYPE:

LOGICAL


DEFAULT:

FALSE

Minimal output


OPTIONS:

TRUE

Increase output level


RECOMMENDATION:

None


CUTOCC

Specifies occupied orbital cutoff


TYPE:

INTEGER: CUTOFF=CUTOCC/100


DEFAULT:

50


OPTIONS:

0-200


RECOMMENDATION:

None


CUTVIR

Specifies virtual orbital cutoff


TYPE:

INTEGER: CUTOFF=CUTVIR/100


DEFAULT:

0

No truncation


OPTIONS:

0-100


RECOMMENDATION:

None


PBHT_ANALYSIS

Controls whether overlap analysis of electronic excitations is performed.


TYPE:

LOGICAL


DEFAULT:

FALSE


OPTIONS:

FALSE

Do not perform overlap analysis

TRUE

Perform overlap analysis


RECOMMENDATION:

None


PBHT_FINE

Increases accuracy of overlap analysis


TYPE:

LOGICAL


DEFAULT:

FALSE

 

OPTIONS:

FALSE

 

TRUE

Increase accuracy of overlap analysis


RECOMMENDATION:

None


SRC_DFT

Selects form of the short-range corrected functional


TYPE:

INTEGER


DEFAULT:

No default


OPTIONS:

1

SRC1 functional

2

SRC2 functional


RECOMMENDATION:

None


OMEGA

Sets the Coulomb attenuation parameter for the short-range component.


TYPE:

INTEGER


DEFAULT:

No default


OPTIONS:

$n$

Corresponding to $\omega = n/1000$, in units of bohr$^{-1}$


RECOMMENDATION:

None


OMEGA2

Sets the Coulomb attenuation parameter for the long-range component.


TYPE:

INTEGER


DEFAULT:

No default


OPTIONS:

$n$

Corresponding to $\omega 2 = n/1000$, in units of bohr$^{-1}$


RECOMMENDATION:

None


HF_SR

Sets the fraction of Hartree-Fock exchange at r$_{12}$=0.


TYPE:

INTEGER


DEFAULT:

No default


OPTIONS:

$n$

Corresponding to HF_SR = $n/1000$


RECOMMENDATION:

None


HF_LR

Sets the fraction of Hartree-Fock exchange at r$_{12}$=$\infty $.


TYPE:

INTEGER


DEFAULT:

No default


OPTIONS:

$n$

Corresponding to HF_LR = $n/1000$


RECOMMENDATION:

None


WANG_ZIEGLER_KERNEL

Controls whether to use the Wang-Ziegler non-collinear exchange-correlation kernel in a SFDFT calculation.


TYPE:

LOGICAL


DEFAULT:

FALSE


OPTIONS:

FALSE

Do not use non-collinear kernel

TRUE

Use non-collinear kernel


RECOMMENDATION:

None


6.3.4 TDDFT Coupled with C-PCM for Excitation Energies and Properties Calculations

As described in Section 11.2 (and especially Section 11.2.2), continuum solvent models such as C-PCM allow one to include solvent effect in the calculations. TDDFT/C-PCM allows excited-state modeling in solution. Q-Chem also features TDDFT coupled with C-PCM which extends TDDFT to calculations of properties of electronically-excited molecules in solution. In particular, TDDFT/C-PCM allows one to perform geometry optimization and vibrational analysis [381].

When TDDFT/C-PCM is applied to calculate vertical excitation energies, the solvent around vertically excited solute is out of equilibrium. While the solvent electron density equilibrates fast to the density of the solute (electronic response), the relaxation of nuclear degrees of freedom (e.g., orientational polarization) takes place on a slower timescale. To describe this situation, an optical dielectric constant is employed. To distinguish between equilibrium and non-equilibrium calculations, two dielectric constants are used in these calculations: a static constant ($\varepsilon _0$), equal to the equilibrium bulk value, and a fast constant ($\varepsilon _{fast}$) related to the response of the medium to high frequency perturbations. For vertical excitation energy calculations (corresponding to the unrelaxed solvent nuclear degrees of freedom), it is recommended to use the optical dielectric constant for $\varepsilon _{fast}$), whereas for the geometry optimization and vibrational frequency calculations, the static dielectric constant should be used [381].

The example below illustrates TDDFT/C-PCM calculations of vertical excitation energies. More information concerning the C-PCM and the various PCM job control options can be found in Section 11.2.

Example 6.101  TDDFT/C-PCM low-lying vertical excitation energy

$molecule
0 1
C           0  0  0.0
O           0  0  1.21
$end

$rem
EXCHANGE            B3lyp
CIS_N_ROOTS         10
cis_singlets        true
cis_triplets        true
RPA                 TRUE
BASIS               6-31+G*
XC_GRID             1
solvent_method      pcm
$end

$pcm
Theory      CPCM
Method      SWIG
Solver      Inversion
Radii       Bondi
$end

$solvent
Dielectric         78.39
OpticalDielectric  1.777849
$end

6.3.5 Analytical Excited-State Hessian in TDDFT

To carry out vibrational frequency analysis of an excited state with TDDFT [382, 383], an optimization of the excited-state geometry is always necessary. Like the vibrational frequency analysis of the ground state, the frequency analysis of the excited state should be also performed at a stationary point on the excited state potential surface. The $rem variable CIS_STATE_DERIV should be set to the excited state for which an optimization and frequency analysis is needed, in addition to the $rem keywords used for an excitation energy calculation.

Compared to the numerical differentiation method, the analytical calculation of geometrical second derivatives of the excitation energy needs much less time but much more memory. The computational cost is mainly consumed by the steps to solve both the CPSCF equations for the derivatives of molecular orbital coefficients C$^{\mathrm{x}}$ and the CP-TDDFT equations for the derivatives of the transition vectors, as well as to build the Hessian matrix. The memory usages for these steps scale as ${\cal {O}}({3mN^2})$, where $N$ is the number of basis functions and m is the number of atoms. For large systems, it is thus essential to solve all the coupled-perturbed equations in segments. In this case, the $rem variable CPSCF_NSEG is always needed.

In the calculation of the analytical TDDFT excited-state Hessian, one has to evaluate a large number of energy-functional derivatives: the first-order to fourth-order functional derivatives with respect to the density variables as well as their derivatives with respect to the nuclear coordinates. Therefore, a very fine integration grid for DFT calculation should be adapted to guarantee the accuracy of the results.

Analytical TDDFT/C-PCM Hessian has been implemented in Q-Chem. Normal mode analysis for a system in solution can be performed with the frequency calculation by TDDFT/C-PCM method. the $rem and $pcm variables for the excited state calculation with TDDFT/C-PCM included in the vertical excitation energy example above are needed. When the properties of large systems are calculated, you must pay attention to the memory limit. At present, only a few exchange correlation functionals, including Slater+VWN, BLYP, B3LYP, are available for the analytical Hessian calculation.

Example 6.102  A B3LYP/6-31G* optimization in gas phase, followed by a frequency analysis for the first excited state of the peroxy radical

$molecule
0 2
C  1.004123 -0.180454  0.000000
O -0.246002  0.596152  0.000000
O -1.312366 -0.230256  0.000000
H  1.810765  0.567203  0.000000
H  1.036648 -0.805445 -0.904798
H  1.036648 -0.805445  0.904798
$end

$rem
jobtype opt
exchange b3lyp
cis_state_deriv 1
basis 6-31G*
cis_n_roots 10
cis_singlets true
cis_triplets false
xc_grid 000075000302
RPA 0
$end

@@@

$molecule 
Read
$end

$rem
jobtype freq
exchange b3lyp
cis_state_deriv 1
basis 6-31G*
cis_n_roots 10
cis_singlets true
cis_triplets false
RPA 0
xc_grid 000075000302
$end

Example 6.103  The optimization and Hessian calculation for low-lying excited state with TDDFT/C-PCM

$comment
9-Fluorenone + 2 methanol in methanol solution
$end

$molecule
0 1
          6   -1.987249    0.699711    0.080583
          6   -1.987187   -0.699537   -0.080519
          6   -0.598049   -1.148932   -0.131299
          6    0.282546    0.000160    0.000137
          6   -0.598139    1.149219    0.131479
          6   -0.319285   -2.505397   -0.285378
          6   -1.386049   -3.395376   -0.388447
          6   -2.743097   -2.962480   -0.339290
          6   -3.049918   -1.628487   -0.186285
          6   -3.050098    1.628566    0.186246
          6   -2.743409    2.962563    0.339341
          6   -1.386397    3.395575    0.388596
          6   -0.319531    2.505713    0.285633
          8    1.560568    0.000159    0.000209
          1    0.703016   -2.862338   -0.324093
          1   -1.184909   -4.453877   -0.510447
          1   -3.533126   -3.698795   -0.423022
          1   -4.079363   -1.292006   -0.147755
          1    0.702729    2.862769    0.324437
          1   -1.185378    4.454097    0.510608
          1   -3.533492    3.698831    0.422983
          1   -4.079503    1.291985    0.147594
          8    3.323150    2.119222    0.125454
          1    2.669309    1.389642    0.084386
          6    3.666902    2.489396   -1.208239
          1    4.397551    3.298444   -1.151310
          1    4.116282    1.654650   -1.759486
          1    2.795088    2.849337   -1.768206
          1    2.669205   -1.389382   -0.084343
          8    3.322989   -2.119006   -0.125620
          6    3.666412   -2.489898    1.207974
          1    4.396966   -3.299023    1.150789
          1    4.115800   -1.655485    1.759730
          1    2.794432   -2.850001    1.767593
$end

$rem
jobtype             OPT
EXCHANGE            B3lyp
CIS_N_ROOTS         10
cis_singlets        true
cis_triplets        true
cis_state_deriv     1   Lowest TDDFT state
RPA                 TRUE
BASIS               6-311G**
XC_GRID             000075000302
solvent_method      pcm
$end

$pcm
Theory      CPCM
Method      SWIG
Solver      Inversion
Radii       Bondi
$end

$solvent
Dielectric 32.613
$end

@@@

$molecule
read
$end

$rem
jobtype             freq
EXCHANGE            B3lyp
CIS_N_ROOTS         10
cis_singlets        true
cis_triplets        true
RPA                 TRUE
cis_state_deriv     1   Lowest TDDFT state
BASIS               6-311G**
XC_GRID             000075000302
solvent_method      pcm
mem_static          4000
mem_total           24000
cpscf_nseg    3
$end

$pcm
Theory      CPCM
Method      SWIG
Solver      Inversion
Radii       Bondi
$end

$solvent
Dielectric 32.613
$end

6.3.6 Various TDDFT-Based Examples

Example 6.104  This example shows two jobs which request variants of time-dependent density functional theory calculations. The first job, using the default value of RPA = FALSE, performs TDDFT in the Tamm-Dancoff approximation (TDA). The second job, with RPA = TRUE performs a both TDA and full TDDFT calculations.

$comment
   methyl peroxy radical
   TDDFT/TDA and full TDDFT with 6-31+G*
$end

$molecule
   0  2
   C  1.00412  -0.18045    0.00000
   O -0.24600   0.59615    0.00000
   O -1.31237  -0.23026    0.00000
   H  1.81077   0.56720    0.00000
   H  1.03665  -0.80545   -0.90480
   H  1.03665  -0.80545    0.90480
$end

$rem
   EXCHANGE          b
   CORRELATION       lyp
   CIS_N_ROOTS       5
   BASIS             6-31+G*
   SCF_CONVERGENCE   7
$end

@@@

$molecule
   read
$end

$rem
   EXCHANGE          b
   CORRELATION       lyp
   CIS_N_ROOTS       5
   RPA               true
   BASIS             6-31+G*
   SCF_CONVERGENCE   7
$end

Example 6.105  This example shows a calculation of the excited states of a formamide-water complex within a reduced excitation space of the orbitals located on formamide

$comment
   formamide-water
   TDDFT/TDA in reduced excitation space
$end

$molecule
   0 1
   H  1.13  0.49 -0.75
   C  0.31  0.50 -0.03
   N -0.28 -0.71  0.08
   H -1.09 -0.75  0.67
   H  0.23 -1.62 -0.22
   O -0.21  1.51  0.47
   O -2.69  1.94 -0.59
   H -2.59  2.08 -1.53
   H -1.83  1.63 -0.30
$end

$rem
   EXCHANGE          b3lyp
   CIS_N_ROOTS       10
   BASIS             6-31++G**
   TRNSS             TRUE
   TRTYPE            1
   CUTOCC            60
   CUTVIR            40
   CISTR_PRINT       TRUE
$end

$solute
1
2
3
4
5
6
$end

Example 6.106  This example shows a calculation of the core-excited states at the oxygen $K$-edge of CO with a short-range corrected functional.

$comment
   TDDFT with short-range corrected (SRC1) functional for the 
   oxygen K-edge of CO
$end

$molecule
0 1
C           0.000000    0.000000   -0.648906
O           0.000000    0.000000    0.486357
$end

$rem
exchange gen
basis 6-311(2+,2+)G**
cis_n_roots 6
cis_triplets false
trnss true
trtype 3
n_sol 1
src_dft 1
omega  560
omega2 2450
HF_SR  500
HF_LR  170
$end

$solute
1
$end

$XC_Functional
X HF    1.00
X B     1.00
C LYP   0.81
C VWN   0.19
$end

Example 6.107  This example shows a calculation of the core-excited states at the phosphorus $K$-edge with a short-range corrected functional.

$comment
   TDDFT with short-range corrected (SRC2) functional for the 
   phosphorus K-edge of PH3 
$end

$molecule
0 1
H        1.196206      0.000000     -0.469131
P        0.000000      0.000000      0.303157
H       -0.598103     -1.035945     -0.469131
H       -0.598103      1.035945     -0.469131
$end

$rem
exchange gen
basis 6-311(2+,2+)G**
cis_n_roots 6
cis_triplets false
trnss true
trtype 3
n_sol 1
src_dft 2
omega  2200
omega2 1800
HF_SR  910
HF_LR  280
$end

$solute
1
$end

$XC_Functional
X HF    1.00
X B     1.00
C LYP   0.81
C VWN   0.19
$end

Example 6.108  SF-TDDFT SP calculation of the 6 lowest states of the TMM diradical using recommended 50-50 functional

$molecule
0 3
C
C 1 CC1
C 1 CC2   2 A2
C 1 CC2   2 A2    3 180.0
H 2 C2H   1 C2CH  3 0.0
H 2 C2H   1 C2CH  4 0.0
H 3 C3Hu  1 C3CHu 2 0.0
H 3 C3Hd  1 C3CHd 4 0.0
H 4 C3Hu  1 C3CHu 2 0.0
H 4 C3Hd  1 C3CHd 3 0.0

CC1    = 1.35
CC2    = 1.47
C2H    = 1.083
C3Hu   = 1.08
C3Hd   = 1.08
C2CH   = 121.2
C3CHu  = 120.3
C3CHd  = 121.3
A2    = 121.0
$end

$rem
jobtype            SP
EXCHANGE           GENERAL            Exact exchange
BASIS              6-31G*
SCF_GUESS          CORE
SCF_CONVERGENCE    10
MAX_SCF_CYCLES     100
SPIN_FLIP       1
CIS_N_ROOTS        6
CIS_CONVERGENCE   10
MAX_CIS_CYCLES  = 100
$end

$xc_functional
 X HF 0.5
 X S  0.08
 X B  0.42
 C VWN 0.19
 C LYP 0.81
$end

Example 6.109  SFDFT with non-collinear exchange-correlation functional for low-lying states of $\rm CH_2$

$comment
non-collinear SFDFT calculation for CH2
at 3B1 state geometry from EOM-CCSD(fT) calculation
$end

$molecule
0 3
C
H  1 rCH
H  1 rCH  2 HCH

rCH    = 1.0775
HCH    = 133.29 
$end

$rem
JOBTYPE             SP
UNRESTRICTED        TRUE
EXCHANGE            PBE0
BASIS               cc-pVTZ
SPIN_FLIP           1
WANG_ZIEGLER_KERNEL TRUE
SCF_CONVERGENCE     10
CIS_N_ROOTS         6
CIS_CONVERGENCE     10
$end