10.11. grad — Analytical nuclear gradients

The grad module provides gradients for mean-field and correlated methods. For example, the RHF gradient can be computed by:

from pyscf import gto, scf
mol = gto.M(
    atom = [
        ['O' , 0. , 0.     , 0],
        ['H' , 0. , -0.757 , 0.587],
        ['H' , 0. ,  0.757 , 0.587]],
        basis = '631g')
mf = scf.RHF(mol)
mf.kernel()
g = mf.nuc_grad_method()
g.kernel()

10.11.1. Examples

Relevant examples examples/grad/01-scf_grad.py examples/grad/02-dft_grad.py examples/grad/03-mp2_grad.py examples/grad/04-cisd_grad.py examples/grad/05-ccsd_grad.py examples/grad/06-tddft_gradients.py examples/grad/10-excited_state_cisd_grad.py examples/grad/11-excited_state_casci_grad.py examples/grad/16-scan_force.py

10.11.2. Program reference

10.11.2.1. Analytical nuclear gradients

Simple usage:

>>> from pyscf import gto, scf, grad
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz')
>>> mf = scf.RHF(mol).run()
>>> grad.RHF(mf).kernel()

10.11.2.2. casci

CASCI analytical nuclear gradients

Ref. J. Comput. Chem., 5, 589

pyscf.grad.casci.Grad

alias of pyscf.grad.casci.Gradients

class pyscf.grad.casci.Gradients(mc)[source]

Non-relativistic restricted Hartree-Fock gradients

as_scanner(state=None)

Generating a nuclear gradients scanner/solver (for geometry optimizer).

The returned solver is a function. This function requires one argument “mol” as input and returns energy and first order nuclear derivatives.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the nuc-grad object and SCF object (DIIS, 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, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1.1', verbose=0)
>>> mc_grad_scanner = mcscf.CASCI(scf.RHF(mol), 4, 4).nuc_grad_method().as_scanner()
>>> etot, grad = mc_grad_scanner(gto.M(atom='N 0 0 0; N 0 0 1.1'))
>>> etot, grad = mc_grad_scanner(gto.M(atom='N 0 0 0; N 0 0 1.5'))
kernel(mo_coeff=None, ci=None, atmlst=None, state=None, verbose=None)[source]

Kernel function is the main driver of a method. Every method should define the kernel function as the entry of the calculation. Note the return value of kernel function is not strictly defined. It can be anything related to the method (such as the energy, the wave-function, the DFT mesh grids etc.).

pyscf.grad.casci.as_scanner(mcscf_grad, state=None)[source]

Generating a nuclear gradients scanner/solver (for geometry optimizer).

The returned solver is a function. This function requires one argument “mol” as input and returns energy and first order nuclear derivatives.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the nuc-grad object and SCF object (DIIS, 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, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1.1', verbose=0)
>>> mc_grad_scanner = mcscf.CASCI(scf.RHF(mol), 4, 4).nuc_grad_method().as_scanner()
>>> etot, grad = mc_grad_scanner(gto.M(atom='N 0 0 0; N 0 0 1.1'))
>>> etot, grad = mc_grad_scanner(gto.M(atom='N 0 0 0; N 0 0 1.5'))

10.11.2.3. casscf

CASSCF analytical nuclear gradients

Ref. J. Comput. Chem., 5, 589

pyscf.grad.casscf.Grad

alias of pyscf.grad.casscf.Gradients

class pyscf.grad.casscf.Gradients(mc)[source]

Non-relativistic restricted Hartree-Fock gradients

as_scanner()

Generating a nuclear gradients scanner/solver (for geometry optimizer).

The returned solver is a function. This function requires one argument “mol” as input and returns energy and first order nuclear derivatives.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the nuc-grad object and SCF object (DIIS, 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, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1.1', verbose=0)
>>> mc_grad_scanner = mcscf.CASSCF(scf.RHF(mol), 4, 4).nuc_grad_method().as_scanner()
>>> etot, grad = mc_grad_scanner(gto.M(atom='N 0 0 0; N 0 0 1.1'))
>>> etot, grad = mc_grad_scanner(gto.M(atom='N 0 0 0; N 0 0 1.5'))
kernel(mo_coeff=None, ci=None, atmlst=None, verbose=None)[source]

Kernel function is the main driver of a method. Every method should define the kernel function as the entry of the calculation. Note the return value of kernel function is not strictly defined. It can be anything related to the method (such as the energy, the wave-function, the DFT mesh grids etc.).

pyscf.grad.casscf.as_scanner(mcscf_grad)[source]

Generating a nuclear gradients scanner/solver (for geometry optimizer).

The returned solver is a function. This function requires one argument “mol” as input and returns energy and first order nuclear derivatives.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the nuc-grad object and SCF object (DIIS, 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, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1.1', verbose=0)
>>> mc_grad_scanner = mcscf.CASSCF(scf.RHF(mol), 4, 4).nuc_grad_method().as_scanner()
>>> etot, grad = mc_grad_scanner(gto.M(atom='N 0 0 0; N 0 0 1.1'))
>>> etot, grad = mc_grad_scanner(gto.M(atom='N 0 0 0; N 0 0 1.5'))

10.11.2.4. ccsd

CCSD analytical nuclear gradients

pyscf.grad.ccsd.Grad

alias of pyscf.grad.ccsd.Gradients

class pyscf.grad.ccsd.Gradients(method)[source]
as_scanner()

Generating a nuclear gradients scanner/solver (for geometry optimizer).

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

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the CCSD 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, cc
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> cc_scanner = cc.CCSD(scf.RHF(mol)).nuc_grad_method().as_scanner()
>>> e_tot, grad = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
kernel(t1=None, t2=None, l1=None, l2=None, eris=None, atmlst=None, verbose=None)[source]

Kernel function is the main driver of a method. Every method should define the kernel function as the entry of the calculation. Note the return value of kernel function is not strictly defined. It can be anything related to the method (such as the energy, the wave-function, the DFT mesh grids etc.).

pyscf.grad.ccsd.as_scanner(grad_cc)[source]

Generating a nuclear gradients scanner/solver (for geometry optimizer).

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

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the CCSD 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, cc
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> cc_scanner = cc.CCSD(scf.RHF(mol)).nuc_grad_method().as_scanner()
>>> e_tot, grad = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))

10.11.2.5. ccsd_slow

RCCSD

Ref: JCP 90, 1752 (1989); DOI:10.1063/1.456069

10.11.2.6. ccsd_t

class pyscf.grad.ccsd_t.Gradients(method)[source]

10.11.2.7. cisd

CISD analytical nuclear gradients

pyscf.grad.cisd.Grad

alias of pyscf.grad.cisd.Gradients

class pyscf.grad.cisd.Gradients(myci)[source]
as_scanner(state=0)

Generating a nuclear gradients scanner/solver (for geometry optimizer).

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

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the CISD 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, ci
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> ci_scanner = ci.CISD(scf.RHF(mol)).nuc_grad_method().as_scanner()
>>> e_tot, grad = ci_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = ci_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
kernel(civec=None, eris=None, atmlst=None, state=None, verbose=None)[source]

Kernel function is the main driver of a method. Every method should define the kernel function as the entry of the calculation. Note the return value of kernel function is not strictly defined. It can be anything related to the method (such as the energy, the wave-function, the DFT mesh grids etc.).

pyscf.grad.cisd.as_scanner(grad_ci, state=0)[source]

Generating a nuclear gradients scanner/solver (for geometry optimizer).

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

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the CISD 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, ci
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> ci_scanner = ci.CISD(scf.RHF(mol)).nuc_grad_method().as_scanner()
>>> e_tot, grad = ci_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = ci_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))

10.11.2.8. dhf

Relativistic Dirac-Hartree-Fock

pyscf.grad.dhf.Grad

alias of pyscf.grad.dhf.Gradients

class pyscf.grad.dhf.Gradients(scf_method)[source]

Unrestricted Dirac-Hartree-Fock gradients

as_scanner()

Generating a nuclear gradients scanner/solver (for geometry optimizer).

The returned solver is a function. This function requires one argument “mol” as input and returns energy and first order nuclear derivatives.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the nuc-grad object and SCF object (DIIS, 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, grad
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> hf_scanner = scf.RHF(mol).apply(grad.RHF).as_scanner()
>>> e_tot, grad = hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
extra_force(atom_id, envs)[source]

Hook for extra contributions in analytical gradients.

Contributions like the response of auxiliary basis in density fitting method, the grid response in DFT numerical integration can be put in this function.

kernel(mo_energy=None, mo_coeff=None, mo_occ=None, atmlst=None)[source]

Kernel function is the main driver of a method. Every method should define the kernel function as the entry of the calculation. Note the return value of kernel function is not strictly defined. It can be anything related to the method (such as the energy, the wave-function, the DFT mesh grids etc.).

class pyscf.grad.dhf.GradientsBasics(method)[source]

Basic nuclear gradient functions for 4C relativistic methods

pyscf.grad.dhf.get_coulomb_hf(mol, dm, level='SSSS')[source]

Dirac-Hartree-Fock Coulomb repulsion

pyscf.grad.dhf.get_veff(mol, dm, level='SSSS')

Dirac-Hartree-Fock Coulomb repulsion

10.11.2.9. mp2

MP2 analytical nuclear gradients

pyscf.grad.mp2.Grad

alias of pyscf.grad.mp2.Gradients

class pyscf.grad.mp2.Gradients(method)[source]
as_scanner()

Generating a nuclear gradients scanner/solver (for geometry optimizer).

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

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the MP2 and the underlying SCF objects (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, mp
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> mp2_scanner = mp.MP2(scf.RHF(mol)).nuc_grad_method().as_scanner()
>>> e_tot, grad = mp2_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = mp2_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
kernel(t2=None, atmlst=None, verbose=None)[source]

Kernel function is the main driver of a method. Every method should define the kernel function as the entry of the calculation. Note the return value of kernel function is not strictly defined. It can be anything related to the method (such as the energy, the wave-function, the DFT mesh grids etc.).

pyscf.grad.mp2.as_scanner(grad_mp)[source]

Generating a nuclear gradients scanner/solver (for geometry optimizer).

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

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the MP2 and the underlying SCF objects (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, mp
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> mp2_scanner = mp.MP2(scf.RHF(mol)).nuc_grad_method().as_scanner()
>>> e_tot, grad = mp2_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = mp2_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))

10.11.2.10. rhf

Non-relativistic Hartree-Fock analytical nuclear gradients

pyscf.grad.rhf.Grad

alias of pyscf.grad.rhf.Gradients

class pyscf.grad.rhf.Gradients(method)[source]

Non-relativistic restricted Hartree-Fock gradients

as_scanner()

Generating a nuclear gradients scanner/solver (for geometry optimizer).

The returned solver is a function. This function requires one argument “mol” as input and returns energy and first order nuclear derivatives.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the nuc-grad object and SCF object (DIIS, 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, grad
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> hf_scanner = scf.RHF(mol).apply(grad.RHF).as_scanner()
>>> e_tot, grad = hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
extra_force(atom_id, envs)[source]

Hook for extra contributions in analytical gradients.

Contributions like the response of auxiliary basis in density fitting method, the grid response in DFT numerical integration can be put in this function.

grad_elec(mo_energy=None, mo_coeff=None, mo_occ=None, atmlst=None)

Electronic part of RHF/RKS gradients

Args:

mf_grad : grad.rhf.Gradients or grad.rks.Gradients object

kernel(mo_energy=None, mo_coeff=None, mo_occ=None, atmlst=None)[source]

Kernel function is the main driver of a method. Every method should define the kernel function as the entry of the calculation. Note the return value of kernel function is not strictly defined. It can be anything related to the method (such as the energy, the wave-function, the DFT mesh grids etc.).

class pyscf.grad.rhf.GradientsBasics(method)[source]

Basic nuclear gradient functions for non-relativistic methods

as_scanner()[source]

Generate Gradients Scanner

get_jk(mol=None, dm=None, hermi=0)[source]

J = ((-nabla i) j| kl) D_lk K = ((-nabla i) j| kl) D_jk

grad(*args, **kwargs)

An alias to method kernel Function Signature: grad(self) —————————————-

kernel()[source]

Kernel function is the main driver of a method. Every method should define the kernel function as the entry of the calculation. Note the return value of kernel function is not strictly defined. It can be anything related to the method (such as the energy, the wave-function, the DFT mesh grids etc.).

optimizer(solver='geometric')[source]

Geometry optimization solver

Kwargs:

solver (string) : geometry optimization solver, can be “geomeTRIC” (default) or “berny”.

symmetrize(de, atmlst=None)[source]

Symmetrize the gradients wrt the point group symmetry of the molecule.

pyscf.grad.rhf.as_scanner(mf_grad)[source]

Generating a nuclear gradients scanner/solver (for geometry optimizer).

The returned solver is a function. This function requires one argument “mol” as input and returns energy and first order nuclear derivatives.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the nuc-grad object and SCF object (DIIS, 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, grad
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> hf_scanner = scf.RHF(mol).apply(grad.RHF).as_scanner()
>>> e_tot, grad = hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
pyscf.grad.rhf.get_hcore(mol)[source]

Part of the nuclear gradients of core Hamiltonian

pyscf.grad.rhf.get_jk(mol, dm)[source]

J = ((-nabla i) j| kl) D_lk K = ((-nabla i) j| kl) D_jk

pyscf.grad.rhf.get_veff(mf_grad, mol, dm)[source]

NR Hartree-Fock Coulomb repulsion

pyscf.grad.rhf.grad_elec(mf_grad, mo_energy=None, mo_coeff=None, mo_occ=None, atmlst=None)[source]

Electronic part of RHF/RKS gradients

Args:

mf_grad : grad.rhf.Gradients or grad.rks.Gradients object

pyscf.grad.rhf.grad_nuc(mol, atmlst=None)[source]

Derivatives of nuclear repulsion energy wrt nuclear coordinates

pyscf.grad.rhf.make_rdm1e(mo_energy, mo_coeff, mo_occ)[source]

Energy weighted density matrix

pyscf.grad.rhf.symmetrize(mol, de, atmlst=None)[source]

Symmetrize the gradients wrt the point group symmetry of the molecule.

10.11.2.11. rks

Non-relativistic RKS analytical nuclear gradients

pyscf.grad.rks.Grad

alias of pyscf.grad.rks.Gradients

class pyscf.grad.rks.Gradients(mf)[source]
extra_force(atom_id, envs)[source]

Hook for extra contributions in analytical gradients.

Contributions like the response of auxiliary basis in density fitting method, the grid response in DFT numerical integration can be put in this function.

get_veff(mol=None, dm=None)

First order derivative of DFT effective potential matrix (wrt electron coordinates)

Args:

ks_grad : grad.uhf.Gradients or grad.uks.Gradients object

pyscf.grad.rks.get_veff(ks_grad, mol=None, dm=None)[source]

First order derivative of DFT effective potential matrix (wrt electron coordinates)

Args:

ks_grad : grad.uhf.Gradients or grad.uks.Gradients object

pyscf.grad.rks.get_vxc_full_response(ni, mol, grids, xc_code, dms, relativity=0, hermi=1, max_memory=2000, verbose=None)[source]

Full response including the response of the grids

10.11.2.12. rohf

Non-relativistic ROHF analytical nuclear gradients

pyscf.grad.rohf.Grad

alias of pyscf.grad.rohf.Gradients

class pyscf.grad.rohf.Gradients(mf)[source]

Non-relativistic ROHF gradients

10.11.2.13. roks

Non-relativistic ROKS analytical nuclear gradients

pyscf.grad.roks.Grad

alias of pyscf.grad.roks.Gradients

class pyscf.grad.roks.Gradients(mf)[source]

Non-relativistic ROHF gradients

10.11.2.14. tdrhf

pyscf.grad.tdrhf.Grad

alias of pyscf.grad.tdrhf.Gradients

class pyscf.grad.tdrhf.Gradients(td)[source]
as_scanner(state=1)

Generating a nuclear gradients scanner/solver (for geometry optimizer).

The returned solver is a function. This function requires one argument “mol” as input and returns energy and first order nuclear derivatives.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the nuc-grad object and SCF object (DIIS, 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, grad
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> td_grad_scanner = scf.RHF(mol).apply(tdscf.TDA).nuc_grad_method().as_scanner()
>>> e_tot, grad = td_grad_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = td_grad_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
grad_elec(xy, singlet, atmlst=None)[source]

Electronic part of TDA, TDHF nuclear gradients

Args:

td_grad : grad.tdrhf.Gradients or grad.tdrks.Gradients object.

x_ya two-element list of numpy arrays

TDDFT X and Y amplitudes. If Y is set to 0, this function computes TDA energy gradients.

kernel(xy=None, state=None, singlet=None, atmlst=None)[source]
Args:
stateint

Excited state ID. state = 1 means the first excited state.

pyscf.grad.tdrhf.as_scanner(td_grad, state=1)[source]

Generating a nuclear gradients scanner/solver (for geometry optimizer).

The returned solver is a function. This function requires one argument “mol” as input and returns energy and first order nuclear derivatives.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the nuc-grad object and SCF object (DIIS, 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, grad
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> td_grad_scanner = scf.RHF(mol).apply(tdscf.TDA).nuc_grad_method().as_scanner()
>>> e_tot, grad = td_grad_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = td_grad_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
pyscf.grad.tdrhf.grad_elec(td_grad, x_y, singlet=True, atmlst=None, max_memory=2000, verbose=4)[source]

Electronic part of TDA, TDHF nuclear gradients

Args:

td_grad : grad.tdrhf.Gradients or grad.tdrks.Gradients object.

x_ya two-element list of numpy arrays

TDDFT X and Y amplitudes. If Y is set to 0, this function computes TDA energy gradients.

10.11.2.15. tdrhf_slow

10.11.2.16. tdrks

pyscf.grad.tdrks.Grad

alias of pyscf.grad.tdrks.Gradients

class pyscf.grad.tdrks.Gradients(td)[source]
grad_elec(xy, singlet, atmlst=None)[source]

Electronic part of TDA, TDDFT nuclear gradients

Args:

td_grad : grad.tdrhf.Gradients or grad.tdrks.Gradients object.

x_ya two-element list of numpy arrays

TDDFT X and Y amplitudes. If Y is set to 0, this function computes TDA energy gradients.

pyscf.grad.tdrks.grad_elec(td_grad, x_y, singlet=True, atmlst=None, max_memory=2000, verbose=4)[source]

Electronic part of TDA, TDDFT nuclear gradients

Args:

td_grad : grad.tdrhf.Gradients or grad.tdrks.Gradients object.

x_ya two-element list of numpy arrays

TDDFT X and Y amplitudes. If Y is set to 0, this function computes TDA energy gradients.

10.11.2.17. tduhf

pyscf.grad.tduhf.Grad

alias of pyscf.grad.tduhf.Gradients

class pyscf.grad.tduhf.Gradients(td)[source]
grad_elec(xy, singlet=None, atmlst=None)[source]

Electronic part of TDA, TDHF nuclear gradients

Args:

td_grad : grad.tduhf.Gradients or grad.tduks.Gradients object.

x_ya two-element list of numpy arrays

TDDFT X and Y amplitudes. If Y is set to 0, this function computes TDA energy gradients.

pyscf.grad.tduhf.grad_elec(td_grad, x_y, atmlst=None, max_memory=2000, verbose=4)[source]

Electronic part of TDA, TDHF nuclear gradients

Args:

td_grad : grad.tduhf.Gradients or grad.tduks.Gradients object.

x_ya two-element list of numpy arrays

TDDFT X and Y amplitudes. If Y is set to 0, this function computes TDA energy gradients.

10.11.2.18. tduks

pyscf.grad.tduks.Grad

alias of pyscf.grad.tduks.Gradients

class pyscf.grad.tduks.Gradients(td)[source]
grad_elec(xy, singlet=None, atmlst=None)[source]

Electronic part of TDA, TDDFT nuclear gradients

Args:

td_grad : grad.tdrhf.Gradients or grad.tdrks.Gradients object.

x_ya two-element list of numpy arrays

TDDFT X and Y amplitudes. If Y is set to 0, this function computes TDA energy gradients.

pyscf.grad.tduks.grad_elec(td_grad, x_y, atmlst=None, max_memory=2000, verbose=4)[source]

Electronic part of TDA, TDDFT nuclear gradients

Args:

td_grad : grad.tdrhf.Gradients or grad.tdrks.Gradients object.

x_ya two-element list of numpy arrays

TDDFT X and Y amplitudes. If Y is set to 0, this function computes TDA energy gradients.

10.11.2.19. uccsd

UCCSD analytical nuclear gradients

pyscf.grad.uccsd.Grad

alias of pyscf.grad.uccsd.Gradients

class pyscf.grad.uccsd.Gradients(method)[source]

10.11.2.20. uccsd_t

class pyscf.grad.uccsd_t.Gradients(method)[source]

10.11.2.21. ucisd

UCISD analytical nuclear gradients

pyscf.grad.ucisd.Grad

alias of pyscf.grad.ucisd.Gradients

class pyscf.grad.ucisd.Gradients(myci)[source]

10.11.2.22. uhf

Non-relativistic unrestricted Hartree-Fock analytical nuclear gradients

pyscf.grad.uhf.Grad

alias of pyscf.grad.uhf.Gradients

class pyscf.grad.uhf.Gradients(method)[source]

Non-relativistic unrestricted Hartree-Fock gradients

grad_elec(mo_energy=None, mo_coeff=None, mo_occ=None, atmlst=None)

Electronic part of UHF/UKS gradients

Args:

mf_grad : grad.uhf.Gradients or grad.uks.Gradients object

pyscf.grad.uhf.get_veff(mf_grad, mol, dm)[source]

First order derivative of HF potential matrix (wrt electron coordinates)

Args:

mf_grad : grad.uhf.Gradients or grad.uks.Gradients object

pyscf.grad.uhf.grad_elec(mf_grad, mo_energy=None, mo_coeff=None, mo_occ=None, atmlst=None)[source]

Electronic part of UHF/UKS gradients

Args:

mf_grad : grad.uhf.Gradients or grad.uks.Gradients object

pyscf.grad.uhf.make_rdm1e(mo_energy, mo_coeff, mo_occ)[source]

Energy weighted density matrix

10.11.2.23. uks

10.11.2.24. ump2