pyscf.tdscf package#
Submodules#
pyscf.tdscf.dhf module#
TDA and TDHF for no-pair DKS Hamiltonian
- class pyscf.tdscf.dhf.TDA(mf, frozen=None)[source]#
Bases:
TDBase
Tamm-Dancoff approximation
- Attributes:
- conv_tolfloat
Diagonalization convergence tolerance. Default is 1e-9.
- nstatesint
Number of TD states to be computed. Default is 3.
Saved results:
- convergedbool
Diagonalization converged or not
- e1D array
excitation energy for each excited state.
- xyA 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.
- singlet = None#
- class pyscf.tdscf.dhf.TDBase(mf, frozen=None)[source]#
Bases:
TDBase
- analyze(verbose=None)#
- 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) + (ai||jb) B[i,a,j,b] = (ai||bj)
- get_nto(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:
tdobj : TDA, or TDHF, or TDDFT object
- stateint
Excited state ID. state = 1 means the first excited state. If state < 0, state ID is counted from the last excited state.
- Kwargs:
- thresholdfloat
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.
- class pyscf.tdscf.dhf.TDHF(mf, frozen=None)[source]#
Bases:
TDBase
- get_precond(hdiag)#
- singlet = None#
- pyscf.tdscf.dhf.gen_tda_hop(mf, fock_ao=None, with_nlc=True)#
A x
- pyscf.tdscf.dhf.gen_tdhf_operation(mf, fock_ao=None, with_nlc=True)[source]#
Generate function to compute
[ A B ][X] [-B* -A*][Y]
pyscf.tdscf.dks module#
pyscf.tdscf.ghf module#
- class pyscf.tdscf.ghf.TDA(mf, frozen=None)[source]#
Bases:
TDBase
Tamm-Dancoff approximation
- Attributes:
- conv_tolfloat
Diagonalization convergence tolerance. Default is 1e-9.
- nstatesint
Number of TD states to be computed. Default is 3.
Saved results:
- convergedbool
Diagonalization converged or not
- e1D array
excitation energy for each excited state.
- xyA 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.
- singlet = None#
- class pyscf.tdscf.ghf.TDBase(mf, frozen=None)[source]#
Bases:
TDBase
- analyze(verbose=None)#
- 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) + (ai||jb) B[i,a,j,b] = (ai||bj)
- get_nto(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:
tdobj : TDA, or TDHF, or TDDFT object
- stateint
Excited state ID. state = 1 means the first excited state. If state < 0, state ID is counted from the last excited state.
- Kwargs:
- thresholdfloat
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.ghf.gen_tda_hop(mf, fock_ao=None, wfnsym=None, with_nlc=True)#
A x
- Kwargs:
- wfnsymint or str
Point group symmetry irrep symbol or ID for excited CIS wavefunction.
- pyscf.tdscf.ghf.gen_tda_operation(mf, fock_ao=None, wfnsym=None, with_nlc=True)[source]#
A x
- Kwargs:
- wfnsymint or str
Point group symmetry irrep symbol or ID for excited CIS wavefunction.
- pyscf.tdscf.ghf.gen_tdhf_operation(mf, fock_ao=None, wfnsym=None, with_nlc=True)[source]#
Generate function to compute
[ A B ][X] [-B* -A*][Y]
pyscf.tdscf.gks module#
- class pyscf.tdscf.gks.CasidaTDDFT(mf, frozen=None)[source]#
-
Solve the Casida TDDFT formula (A-B)(A+B)(X+Y) = (X+Y)w^2
- get_init_guess(mf, nstates=None, wfnsym=None, return_symmetry=False)#
- pyscf.tdscf.gks.TDDFTNoHybrid#
alias of
CasidaTDDFT
pyscf.tdscf.rhf module#
- class pyscf.tdscf.rhf.TDA(mf, frozen=None)[source]#
Bases:
TDBase
Tamm-Dancoff approximation
- Attributes:
- conv_tolfloat
Diagonalization convergence tolerance. Default is 1e-9.
- nstatesint
Number of TD states to be computed. Default is 3.
Saved results:
- convergedbool
Diagonalization converged or not
- e1D array
excitation energy for each excited state.
- xyA 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.
- get_init_guess(mf, nstates=None, wfnsym=None, return_symmetry=False)[source]#
Generate initial guess for TDA
- Kwargs:
- nstatesint
The number of initial guess vectors.
- wfnsymint or str
The irrep label or ID of the wavefunction.
- return_symmetrybool
Whether to return symmetry labels for initial guess vectors.
- to_gpu(out=None)#
Convert a method to its corresponding GPU variant, and recursively converts all attributes of a method to cupy objects or gpu4pyscf objects.
- class pyscf.tdscf.rhf.TDBase(mf, frozen=None)[source]#
Bases:
StreamObject
- analyze(verbose=None)#
- as_scanner()#
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]
- check_sanity()[source]#
Check input of class/object attributes, check whether a class method is overwritten. It does not check the attributes which are prefixed with “_”. The return value of method set is the object itself. This allows a series of functions/methods to be executed in pipe.
- conv_tol = 1e-05#
- deg_eia_thresh = 0.001#
- property e_tot#
Excited state energies
- exclude_nlc = True#
- get_ab(mf=None, frozen=None)[source]#
A and B matrices for TDDFT response function.
A[i,a,j,b] = delta_{ab}delta_{ij}(E_a - E_i) + (ai||jb) B[i,a,j,b] = (ai||bj)
Ref: Chem Phys Lett, 256, 454
- get_frozen_mask()#
Get boolean mask for the restricted reference orbitals.
In the returned boolean (mask) array of frozen orbital indices, the element is False if it corresponds to the frozen orbital.
See mp2.get_frozen_mask
- get_nto(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:
tdobj : TDA, or TDHF, or TDDFT object
- stateint
Excited state ID. state = 1 means the first excited state. If state < 0, state ID is counted from the last excited state.
- Kwargs:
- thresholdfloat
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.
- level_shift = 0#
- lindep = 1e-12#
- max_cycle = 100#
- property nroots#
- nstates = 3#
- oscillator_strength(e=None, xy=None, gauge='length', order=0)#
- positive_eig_threshold = 0.001#
- singlet = True#
- transition_dipole(xy=None)#
Transition dipole moments in the length gauge
- transition_magnetic_dipole(xy=None)#
Transition magnetic dipole moments (imaginary part only)
- transition_magnetic_quadrupole(xy=None)#
Transition magnetic quadrupole moments (imaginary part only)
- transition_octupole(xy=None)#
Transition octupole moments in the length gauge
- transition_quadrupole(xy=None)#
Transition quadrupole moments in the length gauge
- transition_velocity_dipole(xy=None)#
Transition dipole moments in the velocity gauge (imaginary part only)
- transition_velocity_octupole(xy=None)#
Transition octupole moments in the velocity gauge (imaginary part only)
- transition_velocity_quadrupole(xy=None)#
Transition quadrupole moments in the velocity gauge (imaginary part only)
- class pyscf.tdscf.rhf.TDHF(mf, frozen=None)[source]#
Bases:
TDBase
Time-dependent Hartree-Fock
- Attributes:
- conv_tolfloat
Diagonalization convergence tolerance. Default is 1e-4.
- nstatesint
Number of TD states to be computed. Default is 3.
Saved results:
- convergedbool
Diagonalization converged or not
- e1D array
excitation energy for each excited state.
- xyA 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.
- to_gpu(out=None)#
Convert a method to its corresponding GPU variant, and recursively converts all attributes of a method to cupy objects or gpu4pyscf objects.
- class pyscf.tdscf.rhf.TD_Scanner(td)[source]#
Bases:
SinglePointScanner
- 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, with_nlc=True)#
Generate function to compute A x
- Kwargs:
- wfnsymint or str
Point group symmetry irrep symbol or ID for excited CIS wavefunction.
- with_nlcboolean
Whether to skip the NLC contribution
- pyscf.tdscf.rhf.gen_tda_operation(mf, fock_ao=None, singlet=True, wfnsym=None, with_nlc=True)[source]#
Generate function to compute A x
- Kwargs:
- wfnsymint or str
Point group symmetry irrep symbol or ID for excited CIS wavefunction.
- with_nlcboolean
Whether to skip the NLC contribution
- pyscf.tdscf.rhf.gen_tdhf_operation(mf, fock_ao=None, singlet=True, wfnsym=None, with_nlc=True)[source]#
Generate function to compute
[ A B ][X] [-B* -A*][Y]
- pyscf.tdscf.rhf.get_ab(mf, frozen=None, 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) + (ai||jb) B[i,a,j,b] = (ai||bj)
Ref: Chem Phys Lett, 256, 454
- pyscf.tdscf.rhf.get_frozen_mask(td)[source]#
Get boolean mask for the restricted reference orbitals.
In the returned boolean (mask) array of frozen orbital indices, the element is False if it corresponds to the frozen orbital.
See mp2.get_frozen_mask
- 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:
tdobj : TDA, or TDHF, or TDDFT object
- stateint
Excited state ID. state = 1 means the first excited state. If state < 0, state ID is counted from the last excited state.
- Kwargs:
- thresholdfloat
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.rks module#
- class pyscf.tdscf.rks.CasidaTDDFT(mf, frozen=None)[source]#
-
Solve the Casida TDDFT formula (A-B)(A+B)(X+Y) = (X+Y)w^2
- get_init_guess(mf, nstates=None, wfnsym=None, return_symmetry=False)#
Generate initial guess for TDA
- Kwargs:
- nstatesint
The number of initial guess vectors.
- wfnsymint or str
The irrep label or ID of the wavefunction.
- return_symmetrybool
Whether to return symmetry labels for initial guess vectors.
- get_precond(hdiag)#
- pyscf.tdscf.rks.TDDFTNoHybrid#
alias of
CasidaTDDFT
- class pyscf.tdscf.rks.dRPA(mf, frozen=None)[source]#
Bases:
CasidaTDDFT
pyscf.tdscf.uhf module#
- class pyscf.tdscf.uhf.TDA(mf, frozen=None)[source]#
Bases:
TDBase
Tamm-Dancoff approximation
- Attributes:
- conv_tolfloat
Diagonalization convergence tolerance. Default is 1e-9.
- nstatesint
Number of TD states to be computed. Default is 3.
Saved results:
- convergedbool
Diagonalization converged or not
- e1D array
excitation energy for each excited state.
- xyA 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.
- singlet = None#
- to_gpu(out=None)#
Convert a method to its corresponding GPU variant, and recursively converts all attributes of a method to cupy objects or gpu4pyscf objects.
- class pyscf.tdscf.uhf.TDBase(mf, frozen=None)[source]#
Bases:
TDBase
- analyze(verbose=None)#
- get_ab(mf=None, frozen=None)[source]#
A and B matrices for TDDFT response function.
A[i,a,j,b] = delta_{ab}delta_{ij}(E_a - E_i) + (ai||jb) B[i,a,j,b] = (ai||bj)
Spin symmetry is considered in the returned A, B lists. List A has three items: (A_aaaa, A_aabb, A_bbbb). A_bbaa = A_aabb.transpose(2,3,0,1). B has three items: (B_aaaa, B_aabb, B_bbbb). B_bbaa = B_aabb.transpose(2,3,0,1).
- get_frozen_mask()#
Get boolean mask for the unrestricted reference orbitals.
In the returned boolean (mask) array of frozen orbital indices, the element is False if it corresponds to the frozen orbital.
- get_nto(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:
- stateint
Excited state ID. state = 1 means the first excited state. If state < 0, state ID is counted from the last excited state.
- Kwargs:
- thresholdfloat
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.
- class pyscf.tdscf.uhf.TDHF(mf, frozen=None)[source]#
Bases:
TDBase
- singlet = None#
- to_gpu(out=None)#
Convert a method to its corresponding GPU variant, and recursively converts all attributes of a method to cupy objects or gpu4pyscf objects.
- pyscf.tdscf.uhf.gen_tda_hop(mf, fock_ao=None, wfnsym=None, with_nlc=True)#
A x
- Kwargs:
- wfnsymint or str
Point group symmetry irrep symbol or ID for excited CIS wavefunction.
- pyscf.tdscf.uhf.gen_tda_operation(mf, fock_ao=None, wfnsym=None, with_nlc=True)[source]#
A x
- Kwargs:
- wfnsymint or str
Point group symmetry irrep symbol or ID for excited CIS wavefunction.
- pyscf.tdscf.uhf.gen_tdhf_operation(mf, fock_ao=None, wfnsym=None, with_nlc=True)[source]#
Generate function to compute
[ A B ][X] [-B* -A*][Y]
- pyscf.tdscf.uhf.get_ab(mf, frozen=None, 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) + (ai||jb) B[i,a,j,b] = (ai||bj)
Spin symmetry is considered in the returned A, B lists. List A has three items: (A_aaaa, A_aabb, A_bbbb). A_bbaa = A_aabb.transpose(2,3,0,1). B has three items: (B_aaaa, B_aabb, B_bbbb). B_bbaa = B_aabb.transpose(2,3,0,1).
- pyscf.tdscf.uhf.get_frozen_mask(td)[source]#
Get boolean mask for the unrestricted reference orbitals.
In the returned boolean (mask) array of frozen orbital indices, the element is False if it corresponds to the frozen orbital.
- pyscf.tdscf.uhf.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:
- stateint
Excited state ID. state = 1 means the first excited state. If state < 0, state ID is counted from the last excited state.
- Kwargs:
- thresholdfloat
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.uks module#
- class pyscf.tdscf.uks.CasidaTDDFT(mf, frozen=None)[source]#
-
Solve the Casida TDDFT formula (A-B)(A+B)(X+Y) = (X+Y)w^2
- get_init_guess(mf, nstates=None, wfnsym=None, return_symmetry=False)#
- get_precond(hdiag)#
- pyscf.tdscf.uks.TDDFTNoHybrid#
alias of
CasidaTDDFT
- class pyscf.tdscf.uks.dRPA(mf, frozen=None)[source]#
Bases:
CasidaTDDFT
Module contents#
- pyscf.tdscf.TD(mf, frozen=None)#