13. tddft — Time dependent density functional theory

13.1. TDHF

pyscf.tdscf.rhf.CIS

alias of TDA

pyscf.tdscf.rhf.RPA

alias of TDHF

class pyscf.tdscf.rhf.TDA(mf)[source]

Tamm-Dancoff approximation

Attributes:
conv_tol
: float
Diagonalization convergence tolerance. Default is 1e-9.
nstates
: int
Number of TD states to be computed. Default is 3.

Saved results:

converged
: bool
Diagonalization converged or not
e
: 1D array
excitation energy for each excited state.
xy
: A list of two 2D arrays
The two 2D arrays are Excitation coefficients X (shape [nocc,nvir]) and de-excitation coefficients Y (shape [nocc,nvir]) for each excited state. (X,Y) are normalized to 1/2 in RHF/RKS methods and normalized to 1 for UHF/UKS methods. In the TDA calculation, Y = 0.
as_scanner(td)

Generating a scanner/solver for TDA/TDHF/TDDFT PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total TDA/TDHF/TDDFT energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the TDA/TDDFT and the underlying SCF objects (conv_tol, max_memory etc) are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, ...) during calculation.

Examples:

>>> from pyscf import gto, scf, tdscf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> td_scanner = tdscf.TDHF(scf.RHF(mol)).as_scanner()
>>> de = td_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))

[ 0.34460866 0.34460866 0.7131453 ] >>> de = td_scanner(gto.M(atom=’H 0 0 0; F 0 0 1.5’)) [ 0.14844013 0.14844013 0.47641829]

e_tot

Excited state energies

gen_vind(mf)[source]

Compute Ax

get_ab(mf=None)[source]

A and B matrices for TDDFT response function.

A[i,a,j,b] = delta_{ab}delta_{ij}(E_a - E_i) + (ia||bj) B[i,a,j,b] = (ia||jb)

get_nto(tdobj, state=1, threshold=0.3, verbose=None)

Natural transition orbital analysis.

The natural transition density matrix between ground state and excited state \(Tia = \langle \Psi_{ex} | i a^\dagger | \Psi_0 \rangle\) can be transformed to diagonal form through SVD \(T = O \sqrt{\lambda} V^\dagger\). O and V are occupied and virtual natural transition orbitals. The diagonal elements \(\lambda\) are the weights of the occupied-virtual orbital pair in the excitation.

Ref: Martin, R. L., JCP, 118, 4775-4777

Note in the TDHF/TDDFT calculations, the excitation part (X) is interpreted as the CIS coefficients and normalized to 1. The de-excitation part (Y) is ignored.

Args:
state
: int
Excited state ID. state = 1 means the first excited state. If state < 0, state ID is counted from the last excited state.
Kwargs:
threshold
: float
Above which the NTO coefficients will be printed in the output.
Returns:
A list (weights, NTOs). NTOs are natural orbitals represented in AO basis. The first N_occ NTOs are occupied NTOs and the rest are virtual NTOs.
kernel(x0=None, nstates=None)[source]

TDA diagonalization solver

transition_dipole(tdobj, xy=None)

Transition dipole moments in the length gauge

transition_magnetic_dipole(tdobj, xy=None)

Transition magnetic dipole moments (imaginary part only)

transition_magnetic_quadrupole(tdobj, xy=None)

Transition magnetic quadrupole moments (imaginary part only)

transition_octupole(tdobj, xy=None)

Transition octupole moments in the length gauge

transition_quadrupole(tdobj, xy=None)

Transition quadrupole moments in the length gauge

transition_velocity_dipole(tdobj, xy=None)

Transition dipole moments in the velocity gauge (imaginary part only)

class pyscf.tdscf.rhf.TDHF(mf)[source]

Time-dependent Hartree-Fock

Attributes:
conv_tol
: float
Diagonalization convergence tolerance. Default is 1e-9.
nstates
: int
Number of TD states to be computed. Default is 3.

Saved results:

converged
: bool
Diagonalization converged or not
e
: 1D array
excitation energy for each excited state.
xy
: A list of two 2D arrays
The two 2D arrays are Excitation coefficients X (shape [nocc,nvir]) and de-excitation coefficients Y (shape [nocc,nvir]) for each excited state. (X,Y) are normalized to 1/2 in RHF/RKS methods and normalized to 1 for UHF/UKS methods. In the TDA calculation, Y = 0.
gen_vind(mf)[source]

Generate function to compute

[ A B][X] [-B -A][Y]

kernel(x0=None, nstates=None)[source]

TDHF diagonalization with non-Hermitian eigenvalue solver

pyscf.tdscf.rhf.TDRHF

alias of TDHF

pyscf.tdscf.rhf.as_scanner(td)[source]

Generating a scanner/solver for TDA/TDHF/TDDFT PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total TDA/TDHF/TDDFT energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the TDA/TDDFT and the underlying SCF objects (conv_tol, max_memory etc) are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, ...) during calculation.

Examples:

>>> from pyscf import gto, scf, tdscf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> td_scanner = tdscf.TDHF(scf.RHF(mol)).as_scanner()
>>> de = td_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))

[ 0.34460866 0.34460866 0.7131453 ] >>> de = td_scanner(gto.M(atom=’H 0 0 0; F 0 0 1.5’)) [ 0.14844013 0.14844013 0.47641829]

pyscf.tdscf.rhf.gen_tda_hop(mf, fock_ao=None, singlet=True, wfnsym=None)

Generate function to compute (A+B)x

Kwargs:
wfnsym
: int or str
Point group symmetry irrep symbol or ID for excited CIS wavefunction.
pyscf.tdscf.rhf.gen_tda_operation(mf, fock_ao=None, singlet=True, wfnsym=None)[source]

Generate function to compute (A+B)x

Kwargs:
wfnsym
: int or str
Point group symmetry irrep symbol or ID for excited CIS wavefunction.
pyscf.tdscf.rhf.gen_tdhf_operation(mf, fock_ao=None, singlet=True, wfnsym=None)[source]

Generate function to compute

[ A B][X] [-B -A][Y]

pyscf.tdscf.rhf.get_ab(mf, mo_energy=None, mo_coeff=None, mo_occ=None)[source]

A and B matrices for TDDFT response function.

A[i,a,j,b] = delta_{ab}delta_{ij}(E_a - E_i) + (ia||bj) B[i,a,j,b] = (ia||jb)

pyscf.tdscf.rhf.get_nto(tdobj, state=1, threshold=0.3, verbose=None)[source]

Natural transition orbital analysis.

The natural transition density matrix between ground state and excited state \(Tia = \langle \Psi_{ex} | i a^\dagger | \Psi_0 \rangle\) can be transformed to diagonal form through SVD \(T = O \sqrt{\lambda} V^\dagger\). O and V are occupied and virtual natural transition orbitals. The diagonal elements \(\lambda\) are the weights of the occupied-virtual orbital pair in the excitation.

Ref: Martin, R. L., JCP, 118, 4775-4777

Note in the TDHF/TDDFT calculations, the excitation part (X) is interpreted as the CIS coefficients and normalized to 1. The de-excitation part (Y) is ignored.

Args:
state
: int
Excited state ID. state = 1 means the first excited state. If state < 0, state ID is counted from the last excited state.
Kwargs:
threshold
: float
Above which the NTO coefficients will be printed in the output.
Returns:
A list (weights, NTOs). NTOs are natural orbitals represented in AO basis. The first N_occ NTOs are occupied NTOs and the rest are virtual NTOs.
pyscf.tdscf.rhf.transition_dipole(tdobj, xy=None)[source]

Transition dipole moments in the length gauge

pyscf.tdscf.rhf.transition_magnetic_dipole(tdobj, xy=None)[source]

Transition magnetic dipole moments (imaginary part only)

pyscf.tdscf.rhf.transition_magnetic_quadrupole(tdobj, xy=None)[source]

Transition magnetic quadrupole moments (imaginary part only)

pyscf.tdscf.rhf.transition_octupole(tdobj, xy=None)[source]

Transition octupole moments in the length gauge

pyscf.tdscf.rhf.transition_quadrupole(tdobj, xy=None)[source]

Transition quadrupole moments in the length gauge

pyscf.tdscf.rhf.transition_velocity_dipole(tdobj, xy=None)[source]

Transition dipole moments in the velocity gauge (imaginary part only)

pyscf.tdscf.rhf.transition_velocity_octupole(tdobj, xy=None)[source]

Transition octupole moments in the velocity gauge (imaginary part only)

pyscf.tdscf.rhf.transition_velocity_quadrupole(tdobj, xy=None)[source]

Transition quadrupole moments in the velocity gauge (imaginary part only)

13.2. TDDFT

class pyscf.tdscf.rks.TDDFTNoHybrid(mf)[source]

Solve (A-B)(A+B)(X+Y) = (X+Y)w^2

kernel(x0=None, nstates=None)[source]

TDDFT diagonalization solver