CP2K Open Source Molecular Dynamics

User Tools

Site Tools

Trace:
howto:tddft

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision		Previous revision
Next revisionBoth sides next revision
howto:tddft [2022/07/18 21:32] – ahehnhowto:tddft [2022/07/19 13:17] – ahehn
Line 1: Line 1:
 ====== How to run a LR-TDDFT calculation for absorption and emission spectroscopy ====== ====== How to run a LR-TDDFT calculation for absorption and emission spectroscopy ======
  
-This is a short tutorial on how to run linear-response TDDFT computations for absorption and emission spectroscopy. The TDDFT module enables a description of excitation energies and excited-state computations within the Tamm-Dancoff approximation (TDA) featuring GGA and hybrid functionals as well as semi-empirical simplified TDA kernels. The details of the implementation can be found in [[https://aip.scitation.org/doi/full/10.1063/1.5078682]] (Ref. [1]) and in [[https://pubs.acs.org/doi/10.1021/acs.jctc.2c00144]] for corresponding excited-state gradients.  +This is a short tutorial on how to run linear-response time-dependent density functional theory (LR-TDDFT) computations for absorption and emission spectroscopy. The TDDFT module enables a description of excitation energies and excited-state computations within the Tamm-Dancoff approximation (TDA) featuring GGA and hybrid functionals as well as semi-empirical simplified TDA kernels. The details of the implementation can be found in [[https://aip.scitation.org/doi/full/10.1063/1.5078682]] [1] and in [[https://pubs.acs.org/doi/10.1021/acs.jctc.2c00144]] [2] for corresponding excited-state gradients.  
-Note that the current module is based on an earlier TDDFT implementation [[https://chimia.ch/chimia/article/view/2005_499]]. +Note that the current module is based on an earlier TDDFT implementation [[https://chimia.ch/chimia/article/view/2005_499]] [3]. 
-Please cite these papers if you were to use the TDDFT module for the computation of excitation energies or excited-state gradients.+Please cite these papers if you were to use the TDDFT module for the computation of excitation energies [1,3] or excited-state gradients [2].
  
 ===== Brief theory recap ===== ===== Brief theory recap =====
Line 11: Line 11:
 \begin{equation} \label{tda_equation} \begin{equation} \label{tda_equation}
 \begin{aligned} \begin{aligned}
-      \mathbf{A} \mathbf{X}^p &= \Omega^p \mathbf{X}^p \, , \\ +      \mathbf{A} \mathbf{X}^p &= \Omega^p \mathbf{S} \mathbf{X}^p \, , \\ 
-      \sum_{\kappa k} [ F_{\mu \kappa \sigma} \delta_{ik} - F_{ik \sigma} S_{\mu \kappa} ] X^p_{\kappa k \sigma} + \sum_{\lambda} K_{\mu \lambda \sigma} [\mathbf{X}^p] C_{\lambda i \sigma} &=  \sum_{\kappa} \Omega^p S_{\mu \kappa} X^p_{\kappa i \sigma} \, . +      \sum_{\kappa k} [ F_{\mu \kappa \sigma} \delta_{ik} - F_{ik \sigma} S_{\mu \kappa} ] X^p_{\kappa k \sigma} + \sum_{\lambda} K_{\mu \lambda \sigma} [\mathbf{D}^{{\rm{\tiny{X}}}p}] C_{\lambda i \sigma} &=  \sum_{\kappa} \Omega^p S_{\mu \kappa} X^p_{\kappa i \sigma} \, . 
     \end{aligned}     \end{aligned}
 \end{equation} \end{equation}
Line 21: Line 21:
 \begin{aligned} \begin{aligned}
 F_{\mu \nu \sigma} [\mathbf{D}] &= h_{\mu \nu} + J_{\mu \nu \sigma} [\mathbf{D}] - a_{\rm{\tiny{EX}}}K^{\rm{\tiny{EX}}}_{\mu \nu \sigma} [\mathbf{D}] + V_{\mu \nu \sigma}^{\rm{\tiny{XC}}} \, , \\ F_{\mu \nu \sigma} [\mathbf{D}] &= h_{\mu \nu} + J_{\mu \nu \sigma} [\mathbf{D}] - a_{\rm{\tiny{EX}}}K^{\rm{\tiny{EX}}}_{\mu \nu \sigma} [\mathbf{D}] + V_{\mu \nu \sigma}^{\rm{\tiny{XC}}} \, , \\
-K_{\mu \nu \sigma} [\mathbf{D}^{\rm{\tiny{X}}}] &=  J_{\mu \nu \sigma} [\mathbf{D}^{\rm{\tiny{X}}}] - a_{\rm{\tiny{EX}}} K^{\rm{\tiny{EX}}}_{\mu \nu \sigma}[\mathbf{D}^{\rm{\tiny{X}}}] + \sum_{\kappa \lambda \sigma'} f^{\rm{\tiny{XC}}}_{\mu \nu \sigma,\kappa \lambda \sigma'} D_{\kappa \lambda \sigma'}^{\rm{\tiny{X}}} \, .+K_{\mu \nu \sigma} [\mathbf{D}^{{\rm{\tiny{X}}}p}] &=  J_{\mu \nu \sigma} [\mathbf{D}^{{\rm{\tiny{X}}}p}] - a_{\rm{\tiny{EX}}} K^{\rm{\tiny{EX}}}_{\mu \nu \sigma}[\mathbf{D}^{{\rm{\tiny{X}}}p}] + \sum_{\kappa \lambda \sigma'} f^{\rm{\tiny{XC}}}_{\mu \nu \sigma,\kappa \lambda \sigma'} D_{\kappa \lambda \sigma'}^{{\rm{\tiny{X}}}p} \, .
 \end{aligned} \end{aligned}
 \end{equation} \end{equation}
Line 30: Line 30:
 \begin{aligned} \begin{aligned}
 D_{\mu \nu \sigma} &= \sum_k C_{\mu k \sigma} C_{\nu k \sigma}^{\rm{T}} \, , \\ D_{\mu \nu \sigma} &= \sum_k C_{\mu k \sigma} C_{\nu k \sigma}^{\rm{T}} \, , \\
-D_{\mu \nu \sigma}^{\rm{\tiny{X}}} &= \frac{1}{2} \sum_{k} ( X^p_{\mu k \sigma} C_{\nu k \sigma}^{\rm{T}} + C_{\mu k \sigma} (X^p_{\nu k \sigma})^{\rm{T}} ) \, .+D_{\mu \nu \sigma}^{{\rm{\tiny{X}}}p} &= \frac{1}{2} \sum_{k} ( X^p_{\mu k \sigma} C_{\nu k \sigma}^{\rm{T}} + C_{\mu k \sigma} (X^p_{\nu k \sigma})^{\rm{T}} ) \, .
 \end{aligned} \end{aligned}
 \end{equation} \end{equation}
Line 49: Line 49:
 Within the current implementation, oscillator strengths can be calculated for molecular systems in the length form and for periodic systems using the velocity form (see Ref.[1]). Within the current implementation, oscillator strengths can be calculated for molecular systems in the length form and for periodic systems using the velocity form (see Ref.[1]).
  
-Based on Eq.\ (\ref{tda_equation}), excited-state gradients can be formulated based on a variational Lagrangian for each excited state $p$,+Based on the TDA eigenvalue problem, excited-state gradients can be formulated based on a variational Lagrangian for each excited state $p$,
  
 \begin{equation} \begin{equation}
Line 80: Line 80:
  
 ==== Excitation energies for acetone ==== ==== Excitation energies for acetone ====
 +
 +The following input is a standard input for calculating excitation energies with the hybrid functional PBE0. It should be noted that it is possible to compute the exact exchange integrals for hybrid functionals analytically at however high computational costs. It is therefore recommended to use the Auxiliary Density Matrix Method (ADMM) [[https://pubs.acs.org/doi/10.1021/ct1002225]] to approximate the exact exchange contribution. For GGAs, the corresponding ADMM sections are not required and it is sufficient to specify the chosen functional in the subsection ''&XC_FUNCTIONAL'' of the ''&XC'' section. 
  
 <code - acetone_tddft.inp> <code - acetone_tddft.inp>
Line 136: Line 138:
       REL_CUTOFF 80       REL_CUTOFF 80
     &END MGRID     &END MGRID
-    &AUXILIARY_DENSITY_MATRIX_METHOD       ! When choosing ADMM to approximate exact exchange, the +    &AUXILIARY_DENSITY_MATRIX_METHOD       ! For hybrid functionals, it is recommended to choose ADMM  
-      METHOD BASIS_PROJECTION              ! environment for ground and excited state has to be +      METHOD BASIS_PROJECTION              ! the ADMM environment for ground and excited state has to be 
       EXCH_SCALING_MODEL NONE              ! identical       EXCH_SCALING_MODEL NONE              ! identical
-      EXCH_CORRECTION_FUNC NONE +      EXCH_CORRECTION_FUNC NONE            ! Triple-zeta auxiliary basis sets are recommended (see below) 
-      ADMM_PURIFICATION_METHOD NONE +      ADMM_PURIFICATION_METHOD NONE        ! For periodic systems (see below), only specific ADMM options 
-    &END AUXILIARY_DENSITY_MATRIX_METHOD+    &END AUXILIARY_DENSITY_MATRIX_METHOD   ! are available
     &POISSON     &POISSON
        PERIODIC NONE        PERIODIC NONE
Line 274: Line 276:
 ==== Excited-state gradient for acetone ==== ==== Excited-state gradient for acetone ====
  
 +To perform an excited-state optimization for emission spectroscopy, the run type has to be set to ''GEO_OPT'' and the state to be optimized has to be specified in the ''&EXCITED_STATES'' section. Note that the number of excited states chosen in the TDDFPT section should be larger or at least equal to the number of the chosen excited state. 
 +
 +<code cp2k>
 +
 + &GLOBAL
 +  PROJECT S20Acetone
 +  RUN_TYPE ENERGY_FORCE                ! The run type has to be changed to ENERGY_FORCE of GEO_OPT
 +  PREFERRED_DIAG_LIBRARY SL
 +  PRINT_LEVEL medium
 + &END GLOBAL
 + &PROPERTIES
 +    &TDDFPT
 +       NSTATES 10                     
 +       MAX_ITER   100                        
 +       CONVERGENCE [eV] 1.0e-7        
 +       RKS_TRIPLETS F 
 +       ADMM_KERNEL_CORRECTION_SYMMETRIC T   ! required keyword when using hybrid functionals and ADMM for 
 +    &END TDDFPT                             ! the exact exchange contribution for ES gradients
 + &END PROPERTIES
 + ...
 + &DFT
 +    &QS
 +      METHOD GPW
 +      EPS_DEFAULT 1.0E-17
 +      EPS_PGF_ORB 1.0E-20
 +    &END QS
 +    &EXCITED_STATES T
 +     STATE 1                     ! the excited state to be optimized has to be specified in this section
 +    &END EXCITED_STATES
 +    &SCF
 +      SCF_GUESS restart
 +      &OT
 +         PRECONDITIONER FULL_ALL
 +         MINIMIZER DIIS
 +      &END OT
 +      &OUTER_SCF
 +         MAX_SCF 900
 +         EPS_SCF 1.0E-7
 +      &END OUTER_SCF
 +      MAX_SCF 10
 +      EPS_SCF 1.0E-7
 +    &END SCF
 + ...   
 + &END DFT
 +
 +</code>
 +
 +The resulting output file contains for each step of the geometry optimization the already explained output of the TDDFPT computation and additionally information on the excited state that is optimized as well as the iterative solution of the Z vector equations.  
 +
 +<code cp2k>
 +
 + !--------------------------- Excited State Energy ----------------------------!
 + Excitation Energy [Hartree]                                        0.1719188002
 + Total Energy [Hartree]                                           -36.3332484083
 + !-----------------------------------------------------------------------------!
 + !--------------------------- Excited State Forces ----------------------------!
 +
 +  Iteration    Method   Restart      Stepsize      Convergence         Time
 +  ------------------------------------------------------------------------------
 +        1        PCG       F         0.00E+00      0.0231470476        0.02
 +        2        PCG       F         0.30E+00      0.0003972862        1.56
 +        3        PCG       F         0.60E+00      0.0000361012        3.11
 +        4        PCG       F         0.62E+00      0.0000031865        4.65
 +        5        PCG       F         0.63E+00      0.0000002183        6.19
 +        6        PCG       F         0.73E+00      0.0000000147        7.73
 +        7        PCG       F         0.70E+00      0.0000000014        9.27
 +        8        PCG       F         0.72E+00      0.0000000001       10.82
 +        9        PCG       F         0.71E+00      0.0000000000       12.35
 + !-----------------------------------------------------------------------------!
 +
 + ENERGY| Total FORCE_EVAL ( QS ) energy [a.u.]:              -36.333248408276823
 +
 + -------------------------------------------------------------------------------
 +
 +</code> 
 +
 +==== Choosing a semi-empirical kernel ====
 +
 +To speed up computation times for broad-band absorption spectra, a semi-empirical simplified Tamm-Dancoff (sTDA) kernel can be chosen by setting ''KERNEL sTDA'' [[https://www.chemie.uni-bonn.de/pctc/mulliken-center/software/stda/stda]]. The semi-empirical electron repulsion operators depend on several empirical parameters, which need to be adjusted depending on the system under investigation. Most importantly, the amount of exact exchange, scaled by adjusting the parameter $a_{\rm{\tiny{EX}}}$ needs to be chosen carefully and is really crucial to ensure a well-balanced treatment of exact exchange in the GS and ES potential energy surfaces. Too large fractions of exchange in the excited state can lead to negative excitation energies. In general, a relatively small amount of exchange $a_{\rm{\tiny{EX}}}=0.2 / 0.1$ is therefore recommended.
 +
 +The most important keywords of the subsection ''sTDA'' are:
 +  * ''FRACTION'': fraction of exact exchange $a_{\rm{\tiny{EX}}}$
 +  * ''MATAGA_NISHIMOTO_CEXP'': keyword to modify the parameter $\alpha$ of $\gamma^{\rm{\tiny{J}}}$
 +  * ''MATAGA_NISHIMOTO_XEXP'': keyword to modify the parameter $\beta$ of $\gamma^{\rm{\tiny{K}}}$
 +  * ''DO_EWALD'': keyword to switch on Ewald summation for the Coulomb contributions, required when treating periodic systems. (Exact exchange is treated within the minimum image convention and does not require any adjustment.)
 +
 +<code cp2k>
 +
 +&PROPERTIES
 +  &TDDFPT                          
 +   KERNEL sTDA     ! switches on the semi-empirical kernel sTDA                   
 +   &sTDA 
 +     FRACTION 0.2  ! it is crucial to adjust the fraction of exact exchange                  
 +   &END sTDA                       
 +   NSTATES 10                      
 +   MAX_ITER   100                  
 +   CONVERGENCE [eV] 1.0e-7         
 +   RKS_TRIPLETS F                  
 +  &END TDDFPT
 +&END PROPERTIES
 +  
 +</code> 
 +
 +In the output, it can be checked that the sTDA kernel was switched on and the correct parameters for $\alpha$, $\beta$ and $a_{\rm{\tiny{EX}}}$ were chosen:
 +
 +<code cp2k>
 +
 + sTDA| HFX Fraction                                                       0.2000
 + sTDA| Mataga-Nishimoto exponent (C)                                      1.5160
 + sTDA| Mataga-Nishimoto exponent (X)                                      0.5660
 + sTDA| TD matrix filter                                           0.10000000E-09
 + 
 + -------------------------------------------------------------------------------
 + -                      sTDA Kernel: Create Matrix SQRT(S)                     -
 + -------------------------------------------------------------------------------
 +
 +</code> 
 + 
 +The subsequent printout is then identical to the printout for GGA and hybrid functional kernels:
 +
 +<code cp2k>
 +
 + R-TDDFPT states of multiplicity 1
 + Transition dipoles calculated using velocity formulation
 +
 +         State    Excitation        Transition dipole (a.u.)        Oscillator
 +         number   energy (eV)       x           y           z     strength (a.u.)
 +         ------------------------------------------------------------------------
 + TDDFPT|      1       4.88982   8.2927E-09  8.6085E-09 -2.2495E-08   7.77366E-17
 + TDDFPT|      2       9.05888   4.5649E-09 -4.9483E-09 -8.2739E-08   1.52940E-15
 + TDDFPT|      3       9.23700   1.4100E-02 -1.2802E-08  3.2962E-08   4.49914E-05
 + TDDFPT|      4       9.52129   1.1297E-08 -4.7865E-08  8.3174E-01   1.61372E-01
 + TDDFPT|      5      10.01790  -9.4987E-09 -4.8324E-01 -5.0387E-08   5.73140E-02
 + TDDFPT|      6      10.86728  -1.6175E-08 -3.6799E-01 -7.5353E-09   3.60529E-02
 + TDDFPT|      7      10.97656  -3.7898E-01 -1.0287E-08 -3.3941E-09   3.86231E-02
 + TDDFPT|      8      11.21831   4.5569E-08 -4.2529E-09  1.9115E-08   6.76117E-16
 + TDDFPT|      9      11.40817   2.7777E-08 -6.8033E-01 -3.7411E-08   1.29362E-01
 + TDDFPT|     10      11.82525  -1.2576E-09  1.3562E-08  3.8592E-01   4.31486E-02
 +
 +</code> 
 +
 +It should be noted that it is possible to combine sTDA with the semi-empirical ground-state reference GFN1-xTB. However, it is then recommended to adjust all parameters [[https://aip.scitation.org/doi/10.1063/1.4959605]] and to apply corrections to shift the virtual KS orbital eigenvalues. A shift can be applied by adding the keyword ''EV_SHIFT'' (for open-shell systems ''EOS_SHIFT'').
 +
 +==== Periodic systems ====
 +
 +Adjustments are required when switching to periodic boundary conditions. As already mentioned above, oscillator strengths should be calculated using the ''BERRY'' or ''VELOCITY'' formulation. When choosing the ''sTDA'' kernel, ''DO_EWALD'' has to be switched to true to activate Ewald summation for Coulomb contributions. When computing ES gradients using ''KERNEL FULL'' in combination with hybrid functionals and ADMM, only the ADMM2 method relying on basis projection is implemented in combination with the default for the exchange functional for the first-order GGA correction term. 
 +
 +<code cp2k>
 +
 +    &AUXILIARY_DENSITY_MATRIX_METHOD
 +      METHOD BASIS_PROJECTION
 +      ADMM_PURIFICATION_METHOD NONE
 +      EXCH_SCALING_MODEL NONE
 +      EXCH_CORRECTION_FUNC PBEX
 +    &END
 +
 +</code> 
 +
 +==== Natural transition orbitals ====
  
-===== FAQ =====+Natural transition orbitals (NTOs) are printed when choosing ''PRINT_LEVEL medium'' in the ''GLOBAL'' section or when adding the keyword ''NTO_ANALYSIS'' to the ''PRINT'' subsection of the ''TDDFPT'' section. For this purpose it is required to generate unoccupied orbitals and the keyword ''LUMO'' enables to adjust the number of unoccupied orbitals. It is possible to print the NTOs as ''CUBE_FILES'' or in Molden format, with the latter being activated with the keyword ''MOS_MOLDEN''.  
  
  
howto/tddft.txt · Last modified: 2024/02/24 10:01 by oschuett