pyscf.scf package#

Submodules#

pyscf.scf.addons module#

pyscf.scf.addons.canonical_orth_(S, thr=1e-07)[source]#

Löwdin’s canonical orthogonalization

pyscf.scf.addons.convert_to_ghf(mf, out=None, remove_df=False)[source]#

Convert the given mean-field object to the generalized HF/KS object

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mf object. If mf is an second order SCF (SOSCF) object, the SOSCF layer will be discarded. Its underlying SCF object mf._scf will be converted.

Args:

mf : SCF object

Kwargs
remove_dfbool

Whether to convert the DF-SCF object to the normal SCF object. This conversion is not applied by default.

Returns:

An generalized SCF object

pyscf.scf.addons.convert_to_rhf(mf, out=None, remove_df=False)[source]#

Convert the given mean-field object to the restricted HF/KS object

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mf object. If mf is an second order SCF (SOSCF) object, the SOSCF layer will be discarded. Its underlying SCF object mf._scf will be converted.

Args:

mf : SCF object

Kwargs
remove_dfbool

Whether to convert the DF-SCF object to the normal SCF object. This conversion is not applied by default.

Returns:

An unrestricted SCF object

pyscf.scf.addons.convert_to_uhf(mf, out=None, remove_df=False)[source]#

Convert the given mean-field object to the unrestricted HF/KS object

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mf object. If mf is an second order SCF (SOSCF) object, the SOSCF layer will be discarded. Its underlying SCF object mf._scf will be converted.

Args:

mf : SCF object

Kwargs
remove_dfbool

Whether to convert the DF-SCF object to the normal SCF object. This conversion is not applied by default.

Returns:

An unrestricted SCF object

pyscf.scf.addons.dynamic_level_shift(mf, factor=1.0)#

Dynamically change the level shift in each SCF cycle. The level shift value is set to (HF energy change * factor)

pyscf.scf.addons.dynamic_level_shift_(mf, factor=1.0)[source]#

Dynamically change the level shift in each SCF cycle. The level shift value is set to (HF energy change * factor)

pyscf.scf.addons.dynamic_occ(mf, tol=0.001)#

Dynamically adjust the occupancy to avoid degeneracy between HOMO and LUMO

pyscf.scf.addons.dynamic_occ_(mf, tol=0.001)[source]#

Dynamically adjust the occupancy to avoid degeneracy between HOMO and LUMO

pyscf.scf.addons.dynamic_sz_(mf)#

For UHF, allowing the Sz value being changed during SCF iteration. Determine occupation of alpha and beta electrons based on energy spectrum

pyscf.scf.addons.fast_newton(mf, mo_coeff=None, mo_occ=None, dm0=None, auxbasis=None, dual_basis=None, **newton_kwargs)[source]#

This is a wrap function which combines several operations. This function first setup the initial guess from density fitting calculation then use for Newton solver and call Newton solver.

Newton solver attributes [max_cycle_inner, max_stepsize, ah_start_tol, ah_conv_tol, ah_grad_trust_region, …] can be passed through **newton_kwargs.

pyscf.scf.addons.float_occ(mf)#

For UHF, allowing the Sz value being changed during SCF iteration. Determine occupation of alpha and beta electrons based on energy spectrum

pyscf.scf.addons.float_occ_(mf)[source]#

For UHF, allowing the Sz value being changed during SCF iteration. Determine occupation of alpha and beta electrons based on energy spectrum

pyscf.scf.addons.follow_state(mf, occorb=None)#
pyscf.scf.addons.follow_state_(mf, occorb=None)[source]#
pyscf.scf.addons.frac_occ(mf, tol=0.001)#

Addons for SCF methods to assign fractional occupancy for degenerated occupied HOMOs.

Examples:

>>> mf = gto.M(atom='O 0 0 0; O 0 0 1', verbose=4).RHF()
>>> mf = scf.addons.frac_occ(mf)
>>> mf.run()
pyscf.scf.addons.frac_occ_(mf, tol=0.001)[source]#

Addons for SCF methods to assign fractional occupancy for degenerated occupied HOMOs.

Examples:

>>> mf = gto.M(atom='O 0 0 0; O 0 0 1', verbose=4).RHF()
>>> mf = scf.addons.frac_occ(mf)
>>> mf.run()
pyscf.scf.addons.get_ghf_orbspin(mo_energy, mo_occ, is_rhf=None)[source]#

Spin of each GHF orbital when the GHF orbitals are converted from RHF/UHF orbitals

For RHF orbitals, the orbspin corresponds to first occupied orbitals then unoccupied orbitals. In the occupied orbital space, if degenerated, first alpha then beta, last the (open-shell) singly occupied (alpha) orbitals. In the unoccupied orbital space, first the (open-shell) unoccupied (beta) orbitals if applicable, then alpha and beta orbitals

For UHF orbitals, the orbspin corresponds to first occupied orbitals then unoccupied orbitals.

pyscf.scf.addons.mom_occ(mf, occorb, setocc)#

Use maximum overlap method to determine occupation number for each orbital in every iteration.

pyscf.scf.addons.mom_occ_(mf, occorb, setocc)[source]#

Use maximum overlap method to determine occupation number for each orbital in every iteration.

pyscf.scf.addons.partial_cholesky_orth_(S, canthr=1e-07, cholthr=1e-09)[source]#

Partial Cholesky orthogonalization for curing overcompleteness.

References:

Susi Lehtola, Curing basis set overcompleteness with pivoted Cholesky decompositions, J. Chem. Phys. 151, 241102 (2019), doi:10.1063/1.5139948.

Susi Lehtola, Accurate reproduction of strongly repulsive interatomic potentials, Phys. Rev. A 101, 032504 (2020), doi:10.1103/PhysRevA.101.032504.

pyscf.scf.addons.project_dm_nr2nr(mol1, dm1, mol2)[source]#

Project density matrix representation from basis set 1 (mol1) to basis set 2 (mol2).

\[ \begin{align}\begin{aligned}|AO2\rangle DM_AO2 \langle AO2|\\= |AO2\rangle P DM_AO1 P \langle AO2|\\DM_AO2 = P DM_AO1 P\\P = S_{AO2}^{-1}\langle AO2|AO1\rangle\end{aligned}\end{align} \]

There are three relevant functions: project_dm_nr2nr() is the projection for non-relativistic (scalar) basis. project_dm_nr2r() projects from non-relativistic to relativistic basis. project_dm_r2r() is the projection between relativistic (spinor) basis.

pyscf.scf.addons.project_dm_nr2r(mol1, dm1, mol2)[source]#

Project density matrix representation from basis set 1 (mol1) to basis set 2 (mol2).

\[ \begin{align}\begin{aligned}|AO2\rangle DM_AO2 \langle AO2|\\= |AO2\rangle P DM_AO1 P \langle AO2|\\DM_AO2 = P DM_AO1 P\\P = S_{AO2}^{-1}\langle AO2|AO1\rangle\end{aligned}\end{align} \]

There are three relevant functions: project_dm_nr2nr() is the projection for non-relativistic (scalar) basis. project_dm_nr2r() projects from non-relativistic to relativistic basis. project_dm_r2r() is the projection between relativistic (spinor) basis.

pyscf.scf.addons.project_dm_r2r(mol1, dm1, mol2)[source]#

Project density matrix representation from basis set 1 (mol1) to basis set 2 (mol2).

\[ \begin{align}\begin{aligned}|AO2\rangle DM_AO2 \langle AO2|\\= |AO2\rangle P DM_AO1 P \langle AO2|\\DM_AO2 = P DM_AO1 P\\P = S_{AO2}^{-1}\langle AO2|AO1\rangle\end{aligned}\end{align} \]

There are three relevant functions: project_dm_nr2nr() is the projection for non-relativistic (scalar) basis. project_dm_nr2r() projects from non-relativistic to relativistic basis. project_dm_r2r() is the projection between relativistic (spinor) basis.

pyscf.scf.addons.project_mo_nr2nr(mol1, mo1, mol2)[source]#

Project orbital coefficients from basis set 1 (C1 for mol1) to basis set 2 (C2 for mol2).

\[ \begin{align}\begin{aligned}|\psi1\rangle = |AO1\rangle C1\\|\psi2\rangle = P |\psi1\rangle = |AO2\rangle S^{-1}\langle AO2| AO1\rangle> C1 = |AO2\rangle> C2\\C2 = S^{-1}\langle AO2|AO1\rangle C1\end{aligned}\end{align} \]

There are three relevant functions: project_mo_nr2nr() is the projection for non-relativistic (scalar) basis. project_mo_nr2r() projects from non-relativistic to relativistic basis. project_mo_r2r() is the projection between relativistic (spinor) basis.

pyscf.scf.addons.project_mo_nr2r(mol1, mo1, mol2)[source]#

Project orbital coefficients from basis set 1 (C1 for mol1) to basis set 2 (C2 for mol2).

\[ \begin{align}\begin{aligned}|\psi1\rangle = |AO1\rangle C1\\|\psi2\rangle = P |\psi1\rangle = |AO2\rangle S^{-1}\langle AO2| AO1\rangle> C1 = |AO2\rangle> C2\\C2 = S^{-1}\langle AO2|AO1\rangle C1\end{aligned}\end{align} \]

There are three relevant functions: project_mo_nr2nr() is the projection for non-relativistic (scalar) basis. project_mo_nr2r() projects from non-relativistic to relativistic basis. project_mo_r2r() is the projection between relativistic (spinor) basis.

pyscf.scf.addons.project_mo_r2r(mol1, mo1, mol2)[source]#

Project orbital coefficients from basis set 1 (C1 for mol1) to basis set 2 (C2 for mol2).

\[ \begin{align}\begin{aligned}|\psi1\rangle = |AO1\rangle C1\\|\psi2\rangle = P |\psi1\rangle = |AO2\rangle S^{-1}\langle AO2| AO1\rangle> C1 = |AO2\rangle> C2\\C2 = S^{-1}\langle AO2|AO1\rangle C1\end{aligned}\end{align} \]

There are three relevant functions: project_mo_nr2nr() is the projection for non-relativistic (scalar) basis. project_mo_nr2r() projects from non-relativistic to relativistic basis. project_mo_r2r() is the projection between relativistic (spinor) basis.

pyscf.scf.addons.remove_linear_dep(mf, threshold=1e-08, lindep=1e-10, cholesky_threshold=1e-10, force_pivoted_cholesky=False)#
Args:
thresholdfloat

The threshold under which the eigenvalues of the overlap matrix are discarded to avoid numerical instability.

lindepfloat

The threshold that triggers the special treatment of the linear dependence issue.

pyscf.scf.addons.remove_linear_dep_(mf, threshold=1e-08, lindep=1e-10, cholesky_threshold=1e-10, force_pivoted_cholesky=False)[source]#
Args:
thresholdfloat

The threshold under which the eigenvalues of the overlap matrix are discarded to avoid numerical instability.

lindepfloat

The threshold that triggers the special treatment of the linear dependence issue.

pyscf.scf.addons.smearing(mf, sigma=None, method='fermi', mu0=None, fix_spin=False)[source]#

Fermi-Dirac or Gaussian smearing

pyscf.scf.addons.smearing_(mf, *args, **kwargs)[source]#

pyscf.scf.atom_hf module#

class pyscf.scf.atom_hf.AtomHF1e(mol)[source]#

Bases: HF1e, AtomSphAverageRHF

eig(f, s)#

Solver for generalized eigenvalue problem

\[HC = SCE\]
class pyscf.scf.atom_hf.AtomSphAverageRHF(mol)[source]#

Bases: RHF

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.

dump_flags(verbose=None)[source]#
eig(f, s)[source]#

Solver for generalized eigenvalue problem

\[HC = SCE\]
get_grad(mo_coeff, mo_occ, fock=None)[source]#

RHF orbital gradients

Args:
mo_coeff2D ndarray

Obital coefficients

mo_occ1D ndarray

Orbital occupancy

fock_ao2D ndarray

Fock matrix in AO representation

Returns:

Gradients in MO representation. It’s a num_occ*num_vir vector.

get_occ(mo_energy=None, mo_coeff=None)[source]#

spherically averaged fractional occupancy

scf(*args, **kwargs)[source]#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
pyscf.scf.atom_hf.AtomSphericAverageRHF#

alias of AtomSphAverageRHF

pyscf.scf.atom_hf.frac_occ(symb, l, atomic_configuration=[[0, 0, 0, 0], [1, 0, 0, 0], [2, 0, 0, 0], [3, 0, 0, 0], [4, 0, 0, 0], [4, 1, 0, 0], [4, 2, 0, 0], [4, 3, 0, 0], [4, 4, 0, 0], [4, 5, 0, 0], [4, 6, 0, 0], [5, 6, 0, 0], [6, 6, 0, 0], [6, 7, 0, 0], [6, 8, 0, 0], [6, 9, 0, 0], [6, 10, 0, 0], [6, 11, 0, 0], [6, 12, 0, 0], [7, 12, 0, 0], [8, 12, 0, 0], [8, 13, 0, 0], [8, 12, 2, 0], [8, 12, 3, 0], [8, 12, 4, 0], [6, 12, 7, 0], [6, 12, 8, 0], [6, 12, 9, 0], [6, 12, 10, 0], [7, 12, 10, 0], [8, 12, 10, 0], [8, 13, 10, 0], [8, 14, 10, 0], [8, 15, 10, 0], [8, 16, 10, 0], [8, 17, 10, 0], [8, 18, 10, 0], [9, 18, 10, 0], [10, 18, 10, 0], [10, 19, 10, 0], [10, 18, 12, 0], [10, 18, 13, 0], [8, 18, 16, 0], [8, 18, 17, 0], [8, 18, 18, 0], [8, 18, 19, 0], [8, 18, 20, 0], [9, 18, 20, 0], [10, 18, 20, 0], [10, 19, 20, 0], [10, 20, 20, 0], [10, 21, 20, 0], [10, 22, 20, 0], [10, 23, 20, 0], [10, 24, 20, 0], [11, 24, 20, 0], [12, 24, 20, 0], [12, 24, 21, 0], [12, 24, 22, 0], [12, 24, 21, 2], [12, 24, 20, 4], [12, 24, 20, 5], [12, 24, 20, 6], [12, 24, 20, 7], [11, 24, 20, 9], [10, 24, 20, 11], [10, 24, 20, 12], [10, 24, 20, 13], [10, 24, 20, 14], [11, 24, 20, 14], [12, 24, 20, 14], [12, 25, 20, 14], [12, 24, 22, 14], [12, 24, 23, 14], [10, 24, 26, 14], [10, 24, 27, 14], [10, 24, 28, 14], [10, 24, 29, 14], [10, 24, 30, 14], [11, 24, 30, 14], [12, 24, 30, 14], [12, 25, 30, 14], [12, 26, 30, 14], [12, 27, 30, 14], [12, 28, 30, 14], [12, 29, 30, 14], [12, 30, 30, 14], [13, 30, 30, 14], [14, 30, 30, 14], [14, 30, 31, 14], [14, 30, 32, 14], [14, 30, 30, 17], [14, 30, 30, 18], [14, 30, 30, 19], [13, 30, 30, 21], [12, 30, 30, 23], [12, 30, 30, 24], [12, 30, 30, 25], [12, 30, 30, 26], [12, 30, 30, 27], [12, 30, 30, 28], [13, 30, 30, 28], [14, 30, 30, 28], [14, 30, 31, 28], [14, 30, 32, 28], [14, 30, 33, 28], [12, 30, 36, 28], [12, 30, 37, 28], [12, 30, 38, 28], [12, 30, 39, 28], [12, 30, 40, 28], [13, 30, 40, 28], [14, 30, 40, 28], [14, 31, 40, 28], [14, 32, 40, 28], [14, 33, 40, 28], [14, 34, 40, 28], [14, 35, 40, 28], [14, 36, 40, 28]])[source]#
pyscf.scf.atom_hf.get_atm_nrhf(mol, atomic_configuration=[[0, 0, 0, 0], [1, 0, 0, 0], [2, 0, 0, 0], [3, 0, 0, 0], [4, 0, 0, 0], [4, 1, 0, 0], [4, 2, 0, 0], [4, 3, 0, 0], [4, 4, 0, 0], [4, 5, 0, 0], [4, 6, 0, 0], [5, 6, 0, 0], [6, 6, 0, 0], [6, 7, 0, 0], [6, 8, 0, 0], [6, 9, 0, 0], [6, 10, 0, 0], [6, 11, 0, 0], [6, 12, 0, 0], [7, 12, 0, 0], [8, 12, 0, 0], [8, 13, 0, 0], [8, 12, 2, 0], [8, 12, 3, 0], [8, 12, 4, 0], [6, 12, 7, 0], [6, 12, 8, 0], [6, 12, 9, 0], [6, 12, 10, 0], [7, 12, 10, 0], [8, 12, 10, 0], [8, 13, 10, 0], [8, 14, 10, 0], [8, 15, 10, 0], [8, 16, 10, 0], [8, 17, 10, 0], [8, 18, 10, 0], [9, 18, 10, 0], [10, 18, 10, 0], [10, 19, 10, 0], [10, 18, 12, 0], [10, 18, 13, 0], [8, 18, 16, 0], [8, 18, 17, 0], [8, 18, 18, 0], [8, 18, 19, 0], [8, 18, 20, 0], [9, 18, 20, 0], [10, 18, 20, 0], [10, 19, 20, 0], [10, 20, 20, 0], [10, 21, 20, 0], [10, 22, 20, 0], [10, 23, 20, 0], [10, 24, 20, 0], [11, 24, 20, 0], [12, 24, 20, 0], [12, 24, 21, 0], [12, 24, 22, 0], [12, 24, 21, 2], [12, 24, 20, 4], [12, 24, 20, 5], [12, 24, 20, 6], [12, 24, 20, 7], [11, 24, 20, 9], [10, 24, 20, 11], [10, 24, 20, 12], [10, 24, 20, 13], [10, 24, 20, 14], [11, 24, 20, 14], [12, 24, 20, 14], [12, 25, 20, 14], [12, 24, 22, 14], [12, 24, 23, 14], [10, 24, 26, 14], [10, 24, 27, 14], [10, 24, 28, 14], [10, 24, 29, 14], [10, 24, 30, 14], [11, 24, 30, 14], [12, 24, 30, 14], [12, 25, 30, 14], [12, 26, 30, 14], [12, 27, 30, 14], [12, 28, 30, 14], [12, 29, 30, 14], [12, 30, 30, 14], [13, 30, 30, 14], [14, 30, 30, 14], [14, 30, 31, 14], [14, 30, 32, 14], [14, 30, 30, 17], [14, 30, 30, 18], [14, 30, 30, 19], [13, 30, 30, 21], [12, 30, 30, 23], [12, 30, 30, 24], [12, 30, 30, 25], [12, 30, 30, 26], [12, 30, 30, 27], [12, 30, 30, 28], [13, 30, 30, 28], [14, 30, 30, 28], [14, 30, 31, 28], [14, 30, 32, 28], [14, 30, 33, 28], [12, 30, 36, 28], [12, 30, 37, 28], [12, 30, 38, 28], [12, 30, 39, 28], [12, 30, 40, 28], [13, 30, 40, 28], [14, 30, 40, 28], [14, 31, 40, 28], [14, 32, 40, 28], [14, 33, 40, 28], [14, 34, 40, 28], [14, 35, 40, 28], [14, 36, 40, 28]])[source]#

pyscf.scf.atom_hf_pp module#

class pyscf.scf.atom_hf_pp.AtomHF1ePP(mol)[source]#

Bases: HF1e, AtomSCFPP

eig(f, s)#

Solver for generalized eigenvalue problem

\[HC = SCE\]
get_hcore(mol=None)#
class pyscf.scf.atom_hf_pp.AtomSCFPP(mol)[source]#

Bases: AtomSphAverageRHF

get_hcore(mol=None)[source]#
pyscf.scf.atom_hf_pp.get_pp_loc(mol)[source]#
pyscf.scf.atom_hf_pp.get_pp_loc_part1_rs(mol, coords)[source]#
pyscf.scf.atom_hf_pp.get_pp_loc_part2(mol)[source]#
pyscf.scf.atom_hf_pp.get_pp_nl(mol)[source]#

pyscf.scf.atom_ks module#

class pyscf.scf.atom_ks.AtomSphAverageRKS(mol, *args, **kwargs)[source]#

Bases: RKS, AtomSphAverageRHF

eig(f, s)#

Solver for generalized eigenvalue problem

\[HC = SCE\]
get_grad(mo_coeff, mo_occ, fock=None)#

RHF orbital gradients

Args:
mo_coeff2D ndarray

Obital coefficients

mo_occ1D ndarray

Orbital occupancy

fock_ao2D ndarray

Fock matrix in AO representation

Returns:

Gradients in MO representation. It’s a num_occ*num_vir vector.

get_occ(mo_energy=None, mo_coeff=None)#

spherically averaged fractional occupancy

pyscf.scf.atom_ks.AtomSphericAverageRKS#

alias of AtomSphAverageRKS

pyscf.scf.atom_ks.get_atm_nrks(mol, atomic_configuration=[[0, 0, 0, 0], [1, 0, 0, 0], [2, 0, 0, 0], [3, 0, 0, 0], [4, 0, 0, 0], [4, 1, 0, 0], [4, 2, 0, 0], [4, 3, 0, 0], [4, 4, 0, 0], [4, 5, 0, 0], [4, 6, 0, 0], [5, 6, 0, 0], [6, 6, 0, 0], [6, 7, 0, 0], [6, 8, 0, 0], [6, 9, 0, 0], [6, 10, 0, 0], [6, 11, 0, 0], [6, 12, 0, 0], [7, 12, 0, 0], [8, 12, 0, 0], [8, 12, 1, 0], [8, 12, 2, 0], [8, 12, 3, 0], [8, 12, 4, 0], [7, 12, 6, 0], [7, 12, 7, 0], [7, 12, 8, 0], [7, 12, 9, 0], [7, 12, 10, 0], [8, 12, 10, 0], [8, 13, 10, 0], [8, 14, 10, 0], [8, 15, 10, 0], [8, 16, 10, 0], [8, 17, 10, 0], [8, 18, 10, 0], [9, 18, 10, 0], [10, 18, 10, 0], [10, 18, 11, 0], [10, 18, 12, 0], [10, 18, 13, 0], [9, 18, 15, 0], [9, 18, 16, 0], [8, 18, 18, 0], [8, 18, 19, 0], [8, 18, 20, 0], [9, 18, 20, 0], [10, 18, 20, 0], [10, 19, 20, 0], [10, 20, 20, 0], [10, 21, 20, 0], [10, 22, 20, 0], [10, 23, 20, 0], [10, 24, 20, 0], [11, 24, 20, 0], [12, 24, 20, 0], [12, 24, 20, 1], [12, 24, 20, 2], [12, 24, 20, 3], [12, 24, 20, 4], [12, 24, 20, 5], [12, 24, 20, 6], [12, 24, 20, 7], [12, 24, 20, 8], [12, 24, 20, 9], [12, 24, 20, 10], [12, 24, 20, 11], [12, 24, 20, 12], [12, 24, 20, 13], [12, 24, 20, 14], [12, 24, 21, 14], [12, 24, 22, 14], [12, 24, 23, 14], [11, 24, 25, 14], [11, 24, 26, 14], [10, 24, 28, 14], [10, 24, 29, 14], [10, 24, 30, 14], [11, 24, 30, 14], [12, 24, 30, 14], [12, 25, 30, 14], [12, 26, 30, 14], [12, 27, 30, 14], [12, 28, 30, 14], [12, 29, 30, 14], [12, 30, 30, 14], [13, 30, 30, 14], [14, 30, 30, 14], [14, 30, 30, 15], [14, 30, 30, 16], [14, 30, 30, 17], [13, 30, 30, 19], [13, 30, 30, 20], [13, 30, 30, 21], [13, 30, 30, 22], [13, 30, 30, 23], [13, 30, 30, 24], [13, 30, 30, 25], [13, 30, 30, 26], [13, 30, 30, 27], [13, 30, 30, 28], [14, 30, 30, 28], [14, 30, 31, 28], [14, 30, 32, 28], [13, 30, 34, 28], [12, 30, 36, 28], [12, 30, 37, 28], [12, 30, 38, 28], [12, 30, 39, 28], [12, 30, 40, 28], [13, 30, 40, 28], [14, 30, 40, 28], [14, 31, 40, 28], [14, 32, 40, 28], [14, 33, 40, 28], [14, 34, 40, 28], [14, 35, 40, 28], [14, 36, 40, 28]], xc='slater', grid=(100, 434))[source]#

pyscf.scf.chkfile module#

pyscf.scf.chkfile.dump_scf(mol, chkfile, e_tot, mo_energy, mo_coeff, mo_occ, overwrite_mol=True)[source]#

save temporary results

pyscf.scf.chkfile.load_scf(chkfile)[source]#

pyscf.scf.cphf module#

Restricted coupled perturbed Hartree-Fock solver

pyscf.scf.cphf.kernel(fvind, mo_energy, mo_occ, h1, s1=None, max_cycle=50, tol=1e-09, hermi=False, verbose=2, level_shift=0)#
Args:
fvindfunction

Given density matrix, compute (ij|kl)D_{lk}*2 - (ij|kl)D_{jk}

Kwargs:
hermiboolean

Whether the matrix defined by fvind is Hermitian or not.

level_shiftfloat

Add to diagonal terms to slightly improve the convergence speed of Krylov solver

pyscf.scf.cphf.solve(fvind, mo_energy, mo_occ, h1, s1=None, max_cycle=50, tol=1e-09, hermi=False, verbose=2, level_shift=0)[source]#
Args:
fvindfunction

Given density matrix, compute (ij|kl)D_{lk}*2 - (ij|kl)D_{jk}

Kwargs:
hermiboolean

Whether the matrix defined by fvind is Hermitian or not.

level_shiftfloat

Add to diagonal terms to slightly improve the convergence speed of Krylov solver

pyscf.scf.cphf.solve_nos1(fvind, mo_energy, mo_occ, h1, max_cycle=50, tol=1e-09, hermi=False, verbose=2, level_shift=0)[source]#

For field independent basis. First order overlap matrix is zero

Kwargs:
level_shiftfloat

Add to diagonal terms to slightly improve the convergence speed of Krylov solver

pyscf.scf.cphf.solve_withs1(fvind, mo_energy, mo_occ, h1, s1, max_cycle=50, tol=1e-09, hermi=False, verbose=2, level_shift=0)[source]#

For field dependent basis. First order overlap matrix is non-zero. The first order orbitals are set to C^1_{ij} = -1/2 S1 e1 = h1 - s1*e0 + (e0_j-e0_i)*c1 + vhf[c1]

Kwargs:
level_shiftfloat

Add to diagonal terms to slightly improve the convergence speed of Krylov solver

Returns:

First order orbital coefficients (in MO basis) and first order orbital energy matrix

pyscf.scf.dhf module#

Dirac Hartree-Fock

class pyscf.scf.dhf.DHF(mol)[source]#

Bases: SCF

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
Attributes for Dirac-Hartree-Fock
with_ssssbool or string, for Dirac-Hartree-Fock only

If False, ignore small component integrals (SS|SS). Default is True.

with_gauntbool, for Dirac-Hartree-Fock only

Default is False.

with_breitbool, for Dirac-Hartree-Fock only

Gaunt + gauge term. Default is False.

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> e0 = mf.scf()
>>> mf = scf.DHF(mol)
>>> e1 = mf.scf()
>>> print('Relativistic effects = %.12f' % (e1-e0))
Relativistic effects = -0.000008854205
Gradients(*args, **kwargs)#

Unrestricted Dirac-Hartree-Fock gradients

TDA(*args, **kwargs)#

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.

TDHF(*args, **kwargs)#
analyze(verbose=None)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis; Diople moment.

build(mol=None)[source]#
conv_tol = 1e-08#
dip_moment(mol=None, dm=None, unit='Debye', verbose=3, **kwargs)[source]#

Dipole moment calculation

\[\begin{split}\mu_x = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|x|\mu) + \sum_A Q_A X_A\\ \mu_y = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|y|\mu) + \sum_A Q_A Y_A\\ \mu_z = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|z|\mu) + \sum_A Q_A Z_A\end{split}\]

where \(\mu_x, \mu_y, \mu_z\) are the x, y and z components of dipole moment

Args:

mol: an instance of Mole dm : a 2D ndarrays density matrices

Return:

A list: the dipole moment on x, y and z component

dump_flags(verbose=None)[source]#
energy_elec(dm=None, h1e=None, vhf=None)#

Electronic part of Dirac-Hartree-Fock energy

Args:

mf : an instance of SCF class

Kwargs:
dm2D ndarray

one-particle density matrix

h1e2D ndarray

Core hamiltonian

vhf2D ndarray

HF potential

Returns:

Hartree-Fock electronic energy and the Coulomb energy

gen_response(mo_coeff=None, mo_occ=None, with_j=True, hermi=0, max_memory=None)#

Generate a function to compute the product of DHF response function and DHF density matrices.

get_grad(mo_coeff, mo_occ, fock=None)[source]#

RHF orbital gradients

Args:
mo_coeff2D ndarray

Obital coefficients

mo_occ1D ndarray

Orbital occupancy

fock_ao2D ndarray

Fock matrix in AO representation

Returns:

Gradients in MO representation. It’s a num_occ*num_vir vector.

get_hcore(mol=None)[source]#
get_jk(mol=None, dm=None, hermi=1, with_j=True, with_k=True, omega=None)[source]#

Compute J, K matrices for all input density matrices

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
hermiint

Whether J, K matrix is hermitian

0 : not hermitian and not symmetric
1 : hermitian or symmetric
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

with_jboolean

Whether to compute J matrices

with_kboolean

Whether to compute K matrices

omegafloat

Parameter of range-separated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-separated Coulomb operator.

Returns:

Depending on the given dm, the function returns one J and one K matrix, or a list of J matrices and a list of K matrices, corresponding to the input density matrices.

Examples:

>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dms = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> j, k = scf.hf.get_jk(mol, dms, hermi=0)
>>> print(j.shape)
(3, 2, 2)
get_occ(mo_energy=None, mo_coeff=None)[source]#

Label the occupancies for each orbital

Kwargs:
mo_energy1D ndarray

Obital energies

mo_coeff2D ndarray

Obital coefficients

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> energy = numpy.array([-10., -1., 1, -2., 0, -3])
>>> mf.get_occ(energy)
array([2, 2, 0, 2, 2, 2])
get_ovlp(mol=None)[source]#
get_veff(mol=None, dm=None, dm_last=0, vhf_last=0, hermi=1)[source]#

Dirac-Coulomb

init_direct_scf(mol=None)[source]#
init_guess_by_atom(mol=None)[source]#

Generate initial guess density matrix from superposition of atomic HF density matrix. The atomic HF is occupancy averaged RHF

Returns:

Density matrix, 2D ndarray

init_guess_by_chkfile(chkfile=None, project=None)[source]#

Read the HF results from checkpoint file, then project it to the basis defined by mol

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

Returns:

Density matrix, 2D ndarray

init_guess_by_huckel(mol=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

Returns:

Density matrix, 2D ndarray

init_guess_by_minao(mol=None)[source]#

Initial guess in terms of the overlap to minimal basis.

init_guess_by_mod_huckel(mol=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

In contrast to init_guess_by_huckel, this routine employs the updated GWH rule from doi:10.1021/ja00480a005 to form the guess.

Returns:

Density matrix, 2D ndarray

init_guess_by_sap(mol=None, **kwargs)[source]#

Generate initial guess density matrix from a superposition of atomic potentials (SAP), doi:10.1021/acs.jctc.8b01089. This is the Gaussian fit implementation, see doi:10.1063/5.0004046.

Args:
molMoleBase object

the molecule object for which the initial guess is evaluated

sap_basisdict

SAP basis in internal format (python dictionary)

Returns:
dm0ndarray

SAP initial guess density matrix

make_rdm1(mo_coeff=None, mo_occ=None, **kwargs)#

One-particle density matrix in AO representation

Args:
mo_coeff2D ndarray

Orbital coefficients. Each column is one orbital.

mo_occ1D ndarray

Occupancy

Returns:

One-particle density matrix, 2D ndarray

mulliken_pop(mol=None, dm=None, s=None, verbose=5)[source]#

Mulliken population analysis

\[M_{ij} = D_{ij} S_{ji}\]

Mulliken charges

\[\delta_i = \sum_j M_{ij}\]
nuc_grad_method()[source]#

Hook to create object for analytical nuclear gradients.

reset(mol=None)[source]#

Reset mol and clean up relevant attributes for scanner mode

scf(dm0=None)[source]#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
sfx2c1e()[source]#
ssss_approx = 'Visscher'#
stability(internal=None, external=None, verbose=None, return_status=False, **kwargs)[source]#

DHF/DKS stability analysis.

See also pyscf.scf.stability.rhf_stability function.

Kwargs:
return_status: bool

Whether to return stable_i and stable_e

Returns:

If return_status is False (default), the return value includes two set of orbitals, which are more close to the stable condition. The first corresponds to the internal stability and the second corresponds to the external stability.

Else, another two boolean variables (indicating current status: stable or unstable) are returned. The first corresponds to the internal stability and the second corresponds to the external stability.

to_dhf()[source]#
to_dks(xc='HF')[source]#

Convert the input mean-field object to a DKS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_ghf()[source]#

Convert the input mean-field object to a GHF object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_gks(xc=None)[source]#

Convert the input mean-field object to a GKS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

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.

to_ks(xc='HF')#

Convert the input mean-field object to a DKS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_rhf()[source]#

Convert the input mean-field object to a RHF/ROHF object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_rks(xc=None)[source]#

Convert the input mean-field object to a RKS/ROKS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_uhf()[source]#

Convert the input mean-field object to a UHF object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_uks(xc=None)[source]#

Convert the input mean-field object to a UKS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

with_breit = False#
with_gaunt = False#
with_ssss = True#
x2c()#
x2c1e()[source]#
class pyscf.scf.dhf.HF1e(mol)[source]#

Bases: DHF

scf(*args)#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
class pyscf.scf.dhf.RDHF(mol)[source]#

Bases: DHF

Kramers restricted Dirac-Hartree-Fock

TDA(*args, **kwargs)#

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.

TDHF(*args, **kwargs)#
to_dks(xc='HF')[source]#

Convert the input mean-field object to a DKS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

x2c()#
x2c1e()[source]#
pyscf.scf.dhf.RHF#

alias of RDHF

pyscf.scf.dhf.UDHF#

alias of DHF

pyscf.scf.dhf.UHF#

alias of DHF

pyscf.scf.dhf.analyze(mf, verbose=5, **kwargs)[source]#
pyscf.scf.dhf.dip_moment(mol, dm, unit='Debye', verbose=3, **kwargs)[source]#

Dipole moment calculation

\[\begin{split}\mu_x = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|x|\mu) + \sum_A Q_A X_A\\ \mu_y = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|y|\mu) + \sum_A Q_A Y_A\\ \mu_z = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|z|\mu) + \sum_A Q_A Z_A\end{split}\]

where \(\mu_x, \mu_y, \mu_z\) are the x, y and z components of dipole moment

Args:

mol: an instance of Mole dm : a 2D ndarrays density matrices

Return:

A list: the dipole moment on x, y and z component

pyscf.scf.dhf.energy_elec(mf, dm=None, h1e=None, vhf=None)[source]#

Electronic part of Dirac-Hartree-Fock energy

Args:

mf : an instance of SCF class

Kwargs:
dm2D ndarray

one-particle density matrix

h1e2D ndarray

Core hamiltonian

vhf2D ndarray

HF potential

Returns:

Hartree-Fock electronic energy and the Coulomb energy

pyscf.scf.dhf.get_grad(mo_coeff, mo_occ, fock_ao)[source]#

DHF Gradients

pyscf.scf.dhf.get_hcore(mol)[source]#
pyscf.scf.dhf.get_init_guess(mol, key='minao', **kwargs)[source]#

Generate density matrix for initial guess

Kwargs:
keystr

One of ‘minao’, ‘atom’, ‘huckel’, ‘mod_huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’.

pyscf.scf.dhf.get_jk(mol, dm, hermi=1, coulomb_allow='SSSS', opt_llll=None, opt_ssll=None, opt_ssss=None, omega=None, verbose=None)#
pyscf.scf.dhf.get_jk_coulomb(mol, dm, hermi=1, coulomb_allow='SSSS', opt_llll=None, opt_ssll=None, opt_ssss=None, omega=None, verbose=None)[source]#
pyscf.scf.dhf.get_ovlp(mol)[source]#
pyscf.scf.dhf.init_guess_by_1e(mol)[source]#

Initial guess from one electron system.

pyscf.scf.dhf.init_guess_by_atom(mol)[source]#

Initial guess from atom calculation.

pyscf.scf.dhf.init_guess_by_chkfile(mol, chkfile_name, project=None)[source]#

Read SCF chkfile and make the density matrix for 4C-DHF initial guess.

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

pyscf.scf.dhf.init_guess_by_huckel(mol)[source]#

Initial guess from on-the-fly Huckel, doi:10.1021/acs.jctc.8b01089.

pyscf.scf.dhf.init_guess_by_minao(mol)[source]#

Initial guess in terms of the overlap to minimal basis.

pyscf.scf.dhf.init_guess_by_mod_huckel(mol)[source]#

Initial guess from on-the-fly Huckel, doi:10.1021/acs.jctc.8b01089, employing the updated GWH rule from doi:10.1021/ja00480a005.

pyscf.scf.dhf.init_guess_by_sap(mol, mf)[source]#

Generate initial guess density matrix from a superposition of atomic potentials (SAP), doi:10.1021/acs.jctc.8b01089. This is the Gaussian fit implementation, see doi:10.1063/5.0004046.

Args:
molMoleBase object

the molecule object for which the initial guess is evaluated

sap_basisdict

SAP basis in internal format (python dictionary)

Returns:
dm0ndarray

SAP initial guess density matrix

pyscf.scf.dhf.kernel(mf, conv_tol=1e-09, conv_tol_grad=None, dump_chk=True, dm0=None, callback=None, conv_check=True)[source]#

the modified SCF kernel for Dirac-Hartree-Fock. In this kernel, the SCF is carried out in three steps. First the 2-electron part is approximated by large component integrals (LL|LL); Next, (SS|LL) the interaction between large and small components are added; Finally, converge the SCF with the small component contributions (SS|SS)

pyscf.scf.dhf.mulliken_pop(mol, dm, s=None, verbose=5)[source]#

Mulliken population analysis

\[M_{ij} = D_{ij} S_{ji}\]

Mulliken charges

\[\delta_i = \sum_j M_{ij}\]
pyscf.scf.dhf.time_reversal_matrix(mol, mat)[source]#

T(A_ij) = A[T(i),T(j)]^*

pyscf.scf.diis module#

DIIS

class pyscf.scf.diis.ADIIS(dev=None, filename=None, incore=False)[source]#

Bases: DIIS

Ref: JCP 132, 054109 (2010); DOI:10.1063/1.3304922

update(s, d, f, mf, h1e, vhf, *args, **kwargs)[source]#

Extrapolate vector

  • If xerr the error vector is given, this function will push the target

vector and error vector in the DIIS subspace, and use the error vector to extrapolate the vector and return the extrapolated vector. * If xerr is None, this function will take the difference between the current given vector and the last given vector as the error vector to extrapolate the vector.

class pyscf.scf.diis.CDIIS(mf=None, filename=None, Corth=None)[source]#

Bases: DIIS

get_num_vec()[source]#
update(s, d, f, *args, **kwargs)[source]#

Extrapolate vector

  • If xerr the error vector is given, this function will push the target

vector and error vector in the DIIS subspace, and use the error vector to extrapolate the vector and return the extrapolated vector. * If xerr is None, this function will take the difference between the current given vector and the last given vector as the error vector to extrapolate the vector.

pyscf.scf.diis.DIIS#

alias of CDIIS

class pyscf.scf.diis.EDIIS(dev=None, filename=None, incore=False)[source]#

Bases: DIIS

SCF-EDIIS Ref: JCP 116, 8255 (2002); DOI:10.1063/1.1470195

update(s, d, f, mf, h1e, vhf, *args, **kwargs)[source]#

Extrapolate vector

  • If xerr the error vector is given, this function will push the target

vector and error vector in the DIIS subspace, and use the error vector to extrapolate the vector and return the extrapolated vector. * If xerr is None, this function will take the difference between the current given vector and the last given vector as the error vector to extrapolate the vector.

pyscf.scf.diis.SCFDIIS#

alias of CDIIS

pyscf.scf.diis.SCF_DIIS#

alias of CDIIS

pyscf.scf.diis.adiis_minimize(ds, fs, idnewest)[source]#
pyscf.scf.diis.ediis_minimize(es, ds, fs)[source]#
pyscf.scf.diis.get_err_vec(s, d, f, Corth=None)[source]#
pyscf.scf.diis.get_err_vec_orig(s, d, f)[source]#

error vector = SDF - FDS

pyscf.scf.diis.get_err_vec_orth(s, d, f, Corth)[source]#

error vector in orthonormal basis = C.T.conj() (SDF - FDS) C

pyscf.scf.dispersion module#

dispersion correction for HF and DFT

pyscf.scf.dispersion.check_disp(mf, disp=None)[source]#

Check whether to apply dispersion correction based on the xc attribute. If dispersion is allowed, return the DFTD3 disp version, such as d3bj, d3zero, d4.

pyscf.scf.dispersion.get_dispersion(mf, disp=None, with_3body=None, verbose=None)[source]#
pyscf.scf.dispersion.parse_dft(xc_code)[source]#

Extract (xc, nlc, disp) from xc_code

pyscf.scf.dispersion.parse_disp(dft_method)[source]#

Decode the disp parameters based on the xc code. Returns xc_code_for_dftd3, disp_version, with_3body

Example: b3lyp-d3bj2b -> (b3lyp, d3bj, False)

wb97x-d3bj -> (wb97x, d3bj, False)

pyscf.scf.ghf module#

Non-relativistic generalized Hartree-Fock

class pyscf.scf.ghf.GHF(mol)[source]#

Bases: SCF

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
Attributes for GHF method

GHF orbital coefficients are 2D array. Let nao be the number of spatial AOs, mo_coeff[:nao] are the coefficients of AO with alpha spin; mo_coeff[nao:nao*2] are the coefficients of AO with beta spin.

CISD(*args, **kwargs)#
MP2(*args, **kwargs)#
TDA(*args, **kwargs)#

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.

TDHF(*args, **kwargs)#
analyze(verbose=None, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis; Diople moment.

convert_from_(mf)[source]#

Create GHF object based on the RHF/UHF object

det_ovlp(mo1, mo2, occ1, occ2, ovlp=None)[source]#

Calculate the overlap between two different determinants. It is the product of single values of molecular orbital overlap matrix.

Return:
A list:

the product of single values: float x_a: \(\mathbf{U} \mathbf{\Lambda}^{-1} \mathbf{V}^\dagger\) They are used to calculate asymmetric density matrix

dip_moment(mol=None, dm=None, unit_symbol='Debye', origin=None, verbose=3)[source]#

Dipole moment calculation

\[\begin{split}\mu_x = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|x|\mu) + \sum_A Q_A X_A\\ \mu_y = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|y|\mu) + \sum_A Q_A Y_A\\ \mu_z = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|z|\mu) + \sum_A Q_A Z_A\end{split}\]

where \(\mu_x, \mu_y, \mu_z\) are the x, y and z components of dipole moment

Args:

mol: an instance of Mole dm : a 2D ndarrays density matrices origin : optional; length 3 list, tuple, or 1D array

Location of the origin. By default, the point (0, 0, 0) is used.

Return:

A list: the dipole moment on x, y and z component

gen_response(mo_coeff=None, mo_occ=None, with_j=True, hermi=0, max_memory=None)#

Generate a function to compute the product of GHF response function and GHF density matrices.

get_grad(mo_coeff, mo_occ, fock=None)[source]#

RHF orbital gradients

Args:
mo_coeff2D ndarray

Obital coefficients

mo_occ1D ndarray

Orbital occupancy

fock_ao2D ndarray

Fock matrix in AO representation

Returns:

Gradients in MO representation. It’s a num_occ*num_vir vector.

get_hcore(mol=None)[source]#
get_init_guess(mol=None, key='minao', **kwargs)#
get_jk(mol=None, dm=None, hermi=0, with_j=True, with_k=True, omega=None)[source]#

Compute J, K matrices for all input density matrices

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
hermiint

Whether J, K matrix is hermitian

0 : not hermitian and not symmetric
1 : hermitian or symmetric
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

with_jboolean

Whether to compute J matrices

with_kboolean

Whether to compute K matrices

omegafloat

Parameter of range-separated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-separated Coulomb operator.

Returns:

Depending on the given dm, the function returns one J and one K matrix, or a list of J matrices and a list of K matrices, corresponding to the input density matrices.

Examples:

>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dms = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> j, k = scf.hf.get_jk(mol, dms, hermi=0)
>>> print(j.shape)
(3, 2, 2)
get_occ(mo_energy=None, mo_coeff=None)#

Label the occupancies for each orbital

Kwargs:
mo_energy1D ndarray

Obital energies

mo_coeff2D ndarray

Obital coefficients

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> energy = numpy.array([-10., -1., 1, -2., 0, -3])
>>> mf.get_occ(energy)
array([2, 2, 0, 2, 2, 2])
get_ovlp(mol=None)[source]#
get_veff(mol=None, dm=None, dm_last=0, vhf_last=0, hermi=1)[source]#

Hartree-Fock potential matrix for the given density matrix

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
dm_lastndarray or a list of ndarrays or 0

The density matrix baseline. If not 0, this function computes the increment of HF potential w.r.t. the reference HF potential matrix.

vhf_lastndarray or a list of ndarrays or 0

The reference HF potential matrix.

hermiint

Whether J, K matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

Returns:

matrix Vhf = 2*J - K. Vhf can be a list matrices, corresponding to the input density matrices.

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dm0 = numpy.random.random((mol.nao_nr(),mol.nao_nr()))
>>> vhf0 = scf.hf.get_veff(mol, dm0, hermi=0)
>>> dm1 = numpy.random.random((mol.nao_nr(),mol.nao_nr()))
>>> vhf1 = scf.hf.get_veff(mol, dm1, hermi=0)
>>> vhf2 = scf.hf.get_veff(mol, dm1, dm_last=dm0, vhf_last=vhf0, hermi=0)
>>> numpy.allclose(vhf1, vhf2)
True
init_guess_by_atom(mol=None)[source]#

Generate initial guess density matrix from superposition of atomic HF density matrix. The atomic HF is occupancy averaged RHF

Returns:

Density matrix, 2D ndarray

init_guess_by_chkfile(chkfile=None, project=None)[source]#

Read the HF results from checkpoint file, then project it to the basis defined by mol

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

Returns:

Density matrix, 2D ndarray

init_guess_by_huckel(mol=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

Returns:

Density matrix, 2D ndarray

init_guess_by_minao(mol=None)[source]#

Generate initial guess density matrix based on ANO basis, then project the density matrix to the basis set defined by mol

Returns:

Density matrix, 2D ndarray

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> scf.hf.init_guess_by_minao(mol)
array([[ 0.94758917,  0.09227308],
       [ 0.09227308,  0.94758917]])
init_guess_by_mod_huckel(mol=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

In contrast to init_guess_by_huckel, this routine employs the updated GWH rule from doi:10.1021/ja00480a005 to form the guess.

Returns:

Density matrix, 2D ndarray

init_guess_by_sap(mol=None, **kwargs)[source]#

Generate initial guess density matrix from a superposition of atomic potentials (SAP), doi:10.1021/acs.jctc.8b01089. This is the Gaussian fit implementation, see doi:10.1063/5.0004046.

Args:
molMoleBase object

the molecule object for which the initial guess is evaluated

sap_basisdict

SAP basis in internal format (python dictionary)

Returns:
dm0ndarray

SAP initial guess density matrix

mulliken_meta(mol=None, dm=None, verbose=5, pre_orth_method='ANO', s=None)[source]#

Mulliken population analysis, based on meta-Lowdin AOs. In the meta-lowdin, the AOs are grouped in three sets: core, valence and Rydberg, the orthogonalization are carried out within each subsets.

Args:

mol : an instance of Mole

dmndarray or 2-item list of ndarray

Density matrix. ROHF dm is a 2-item list of 2D array

Kwargs:

verbose : int or instance of lib.logger.Logger

pre_orth_methodstr

Pre-orthogonalization, which localized GTOs for each atom. To obtain the occupied and unoccupied atomic shells, there are three methods

‘ano’ : Project GTOs to ANO basis
‘minao’ : Project GTOs to MINAO basis
‘scf’ : Symmetry-averaged fractional occupation atomic RHF
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

mulliken_pop(mol=None, dm=None, s=None, verbose=5)[source]#

Mulliken population analysis

\[M_{ij} = D_{ij} S_{ji}\]

Mulliken charges

\[\delta_i = \sum_j M_{ij}\]
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

nuc_grad_method()[source]#

Hook to create object for analytical nuclear gradients.

spin_square(mo_coeff=None, s=None)[source]#

Spin of the GHF wavefunction

\[S^2 = \frac{1}{2}(S_+ S_- + S_- S_+) + S_z^2\]

where \(S_+ = \sum_i S_{i+}\) is effective for all beta occupied orbitals; \(S_- = \sum_i S_{i-}\) is effective for all alpha occupied orbitals.

  1. There are two possibilities for \(S_+ S_-\)
    1. same electron \(S_+ S_- = \sum_i s_{i+} s_{i-}\),

    \[\sum_i \langle UHF|s_{i+} s_{i-}|UHF\rangle = \sum_{pq}\langle p|s_+s_-|q\rangle \gamma_{qp} = n_\alpha\]

    2) different electrons \(S_+ S_- = \sum s_{i+} s_{j-}, (i\neq j)\). There are in total \(n(n-1)\) terms. As a two-particle operator,

    \[\langle S_+ S_- \rangle =\sum_{ij}(\langle i^\alpha|i^\beta\rangle \langle j^\beta|j^\alpha\rangle - \langle i^\alpha|j^\beta\rangle \langle j^\beta|i^\alpha\rangle)\]
  2. Similarly, for \(S_- S_+\)
    1. same electron

    \[\sum_i \langle s_{i-} s_{i+}\rangle = n_\beta\]
    1. different electrons

    \[\langle S_- S_+ \rangle =\sum_{ij}(\langle i^\beta|i^\alpha\rangle \langle j^\alpha|j^\beta\rangle - \langle i^\beta|j^\alpha\rangle \langle j^\alpha|i^\beta\rangle)\]
  3. For \(S_z^2\)
    1. same electron

    \[\langle s_z^2\rangle = \frac{1}{4}(n_\alpha + n_\beta)\]
    1. different electrons

    \[\begin{split}&\sum_{ij}(\langle ij|s_{z1}s_{z2}|ij\rangle -\langle ij|s_{z1}s_{z2}|ji\rangle) \\ &=\frac{1}{4}\sum_{ij}(\langle i^\alpha|i^\alpha\rangle \langle j^\alpha|j^\alpha\rangle - \langle i^\alpha|i^\alpha\rangle \langle j^\beta|j^\beta\rangle - \langle i^\beta|i^\beta\rangle \langle j^\alpha|j^\alpha\rangle + \langle i^\beta|i^\beta\rangle \langle j^\beta|j^\beta\rangle) \\ &-\frac{1}{4}\sum_{ij}(\langle i^\alpha|j^\alpha\rangle \langle j^\alpha|i^\alpha\rangle - \langle i^\alpha|j^\alpha\rangle \langle j^\beta|i^\beta\rangle - \langle i^\beta|j^\beta\rangle \langle j^\alpha|i^\alpha\rangle + \langle i^\beta|j^\beta\rangle\langle j^\beta|i^\beta\rangle) \\ &=\frac{1}{4}\sum_{ij}|\langle i^\alpha|i^\alpha\rangle - \langle i^\beta|i^\beta\rangle|^2 -\frac{1}{4}\sum_{ij}|\langle i^\alpha|j^\alpha\rangle - \langle i^\beta|j^\beta\rangle|^2 \\ &=\frac{1}{4}(n_\alpha - n_\beta)^2 -\frac{1}{4}\sum_{ij}|\langle i^\alpha|j^\alpha\rangle - \langle i^\beta|j^\beta\rangle|^2\end{split}\]
Args:
moa list of 2 ndarrays

Occupied alpha and occupied beta orbitals

Kwargs:
sndarray

AO overlap

Returns:

A list of two floats. The first is the expectation value of S^2. The second is the corresponding 2S+1

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', charge=1, spin=1, verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.kernel()
-75.623975516256706
>>> mo = (mf.mo_coeff[0][:,mf.mo_occ[0]>0], mf.mo_coeff[1][:,mf.mo_occ[1]>0])
>>> print('S^2 = %.7f, 2S+1 = %.7f' % spin_square(mo, mol.intor('int1e_ovlp_sph')))
S^2 = 0.7570150, 2S+1 = 2.0070027
stability(internal=None, external=None, verbose=None, return_status=False, **kwargs)[source]#
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.

to_ks(xc='HF')[source]#

Convert to GKS object.

with_soc = None#
x2c()#

X2C with spin-orbit coupling effects.

Starting from PySCF 2.1, this function (mol.GHF().x2c()) produces an X2C calculation in spherical GTO bases. The results in theory are equivalent to those obtained from the mol.X2C() method, which is computed in the spinor GTO bases. This function called the spin-free X2C1E method in the older versions.

Please note the difference from other SCF methods, such as RHF and UHF. In those methods, the .x2c() method produces a scalar relativistic calculation using the X2C Hamiltonian.

x2c1e()[source]#

X2C with spin-orbit coupling effects.

Starting from PySCF 2.1, this function (mol.GHF().x2c()) produces an X2C calculation in spherical GTO bases. The results in theory are equivalent to those obtained from the mol.X2C() method, which is computed in the spinor GTO bases. This function called the spin-free X2C1E method in the older versions.

Please note the difference from other SCF methods, such as RHF and UHF. In those methods, the .x2c() method produces a scalar relativistic calculation using the X2C Hamiltonian.

class pyscf.scf.ghf.HF1e(mol)[source]#

Bases: GHF

scf(*args)#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
pyscf.scf.ghf.analyze(mf, verbose=5, with_meta_lowdin=True, origin=None, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis; Dipole moment

pyscf.scf.ghf.det_ovlp(mo1, mo2, occ1, occ2, ovlp)[source]#

Calculate the overlap between two different determinants. It is the product of single values of molecular orbital overlap matrix.

Return:
A list:

the product of single values: float x_a: \(\mathbf{U} \mathbf{\Lambda}^{-1} \mathbf{V}^\dagger\) They are used to calculate asymmetric density matrix

pyscf.scf.ghf.dip_moment(mol, dm, unit_symbol='Debye', origin=None, verbose=3)[source]#
pyscf.scf.ghf.get_jk(mol, dm, hermi=0, with_j=True, with_k=True, jkbuild=<function get_jk>, omega=None)[source]#

Compute J, K matrices for all input density matrices

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
hermiint

Whether J, K matrix is hermitian

0 : not hermitian and not symmetric
1 : hermitian or symmetric
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

with_jboolean

Whether to compute J matrices

with_kboolean

Whether to compute K matrices

omegafloat

Parameter of range-separated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-separated Coulomb operator.

Returns:

Depending on the given dm, the function returns one J and one K matrix, or a list of J matrices and a list of K matrices, corresponding to the input density matrices.

Examples:

>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dms = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> j, k = scf.hf.get_jk(mol, dms, hermi=0)
>>> print(j.shape)
(3, 2, 2)
pyscf.scf.ghf.get_occ(mf, mo_energy=None, mo_coeff=None)[source]#
pyscf.scf.ghf.guess_orbspin(mo_coeff)[source]#

Guess the orbital spin (alpha 0, beta 1, unknown -1) based on the orbital coefficients

pyscf.scf.ghf.init_guess_by_chkfile(mol, chkfile_name, project=None)[source]#

Read SCF chkfile and make the density matrix for GHF initial guess.

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

pyscf.scf.ghf.mulliken_meta(mol, dm_ao, verbose=5, pre_orth_method='ANO', s=None)[source]#

Mulliken population analysis, based on meta-Lowdin AOs.

pyscf.scf.ghf.mulliken_pop(mol, dm, s=None, verbose=5)[source]#

Mulliken population analysis

pyscf.scf.ghf.spin_square(mo, s=1)[source]#

Spin of the GHF wavefunction

\[S^2 = \frac{1}{2}(S_+ S_- + S_- S_+) + S_z^2\]

where \(S_+ = \sum_i S_{i+}\) is effective for all beta occupied orbitals; \(S_- = \sum_i S_{i-}\) is effective for all alpha occupied orbitals.

  1. There are two possibilities for \(S_+ S_-\)
    1. same electron \(S_+ S_- = \sum_i s_{i+} s_{i-}\),

    \[\sum_i \langle UHF|s_{i+} s_{i-}|UHF\rangle = \sum_{pq}\langle p|s_+s_-|q\rangle \gamma_{qp} = n_\alpha\]

    2) different electrons \(S_+ S_- = \sum s_{i+} s_{j-}, (i\neq j)\). There are in total \(n(n-1)\) terms. As a two-particle operator,

    \[\langle S_+ S_- \rangle =\sum_{ij}(\langle i^\alpha|i^\beta\rangle \langle j^\beta|j^\alpha\rangle - \langle i^\alpha|j^\beta\rangle \langle j^\beta|i^\alpha\rangle)\]
  2. Similarly, for \(S_- S_+\)
    1. same electron

    \[\sum_i \langle s_{i-} s_{i+}\rangle = n_\beta\]
    1. different electrons

    \[\langle S_- S_+ \rangle =\sum_{ij}(\langle i^\beta|i^\alpha\rangle \langle j^\alpha|j^\beta\rangle - \langle i^\beta|j^\alpha\rangle \langle j^\alpha|i^\beta\rangle)\]
  3. For \(S_z^2\)
    1. same electron

    \[\langle s_z^2\rangle = \frac{1}{4}(n_\alpha + n_\beta)\]
    1. different electrons

    \[\begin{split}&\sum_{ij}(\langle ij|s_{z1}s_{z2}|ij\rangle -\langle ij|s_{z1}s_{z2}|ji\rangle) \\ &=\frac{1}{4}\sum_{ij}(\langle i^\alpha|i^\alpha\rangle \langle j^\alpha|j^\alpha\rangle - \langle i^\alpha|i^\alpha\rangle \langle j^\beta|j^\beta\rangle - \langle i^\beta|i^\beta\rangle \langle j^\alpha|j^\alpha\rangle + \langle i^\beta|i^\beta\rangle \langle j^\beta|j^\beta\rangle) \\ &-\frac{1}{4}\sum_{ij}(\langle i^\alpha|j^\alpha\rangle \langle j^\alpha|i^\alpha\rangle - \langle i^\alpha|j^\alpha\rangle \langle j^\beta|i^\beta\rangle - \langle i^\beta|j^\beta\rangle \langle j^\alpha|i^\alpha\rangle + \langle i^\beta|j^\beta\rangle\langle j^\beta|i^\beta\rangle) \\ &=\frac{1}{4}\sum_{ij}|\langle i^\alpha|i^\alpha\rangle - \langle i^\beta|i^\beta\rangle|^2 -\frac{1}{4}\sum_{ij}|\langle i^\alpha|j^\alpha\rangle - \langle i^\beta|j^\beta\rangle|^2 \\ &=\frac{1}{4}(n_\alpha - n_\beta)^2 -\frac{1}{4}\sum_{ij}|\langle i^\alpha|j^\alpha\rangle - \langle i^\beta|j^\beta\rangle|^2\end{split}\]
Args:
moa list of 2 ndarrays

Occupied alpha and occupied beta orbitals

Kwargs:
sndarray

AO overlap

Returns:

A list of two floats. The first is the expectation value of S^2. The second is the corresponding 2S+1

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', charge=1, spin=1, verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.kernel()
-75.623975516256706
>>> mo = (mf.mo_coeff[0][:,mf.mo_occ[0]>0], mf.mo_coeff[1][:,mf.mo_occ[1]>0])
>>> print('S^2 = %.7f, 2S+1 = %.7f' % spin_square(mo, mol.intor('int1e_ovlp_sph')))
S^2 = 0.7570150, 2S+1 = 2.0070027

pyscf.scf.ghf_symm module#

Non-relativistic generalized Hartree-Fock with point group symmetry.

pyscf.scf.ghf_symm.GHF#

alias of SymAdaptedGHF

class pyscf.scf.ghf_symm.HF1e(mol)[source]#

Bases: SymAdaptedGHF

scf(*args)#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
class pyscf.scf.ghf_symm.SymAdaptedGHF(mol)[source]#

Bases: GHF

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
Attributes for GHF method

GHF orbital coefficients are 2D array. Let nao be the number of spatial AOs, mo_coeff[:nao] are the coefficients of AO with alpha spin; mo_coeff[nao:nao*2] are the coefficients of AO with beta spin.

Attributes for symmetry allowed GHF:
irrep_nelecdict

Specify the number of electrons for particular irrep {‘ir_name’:int, …}. For the irreps not listed in these dicts, the program will choose the occupancy based on the orbital energies.

analyze(verbose=None, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis; Diople moment.

build(mol=None)[source]#
canonicalize(mo_coeff, mo_occ, fock=None)#

Canonicalization diagonalizes the UHF Fock matrix in occupied, virtual subspaces separatedly (without change occupancy).

dump_flags(verbose=None)[source]#
eig(h, s, symm_orb=None, irrep_id=None)[source]#

Solver for generalized eigenvalue problem

\[HC = SCE\]
get_grad(mo_coeff, mo_occ, fock=None)[source]#

RHF orbital gradients

Args:
mo_coeff2D ndarray

Obital coefficients

mo_occ1D ndarray

Orbital occupancy

fock_ao2D ndarray

Fock matrix in AO representation

Returns:

Gradients in MO representation. It’s a num_occ*num_vir vector.

get_irrep_nelec(mol=None, mo_coeff=None, mo_occ=None, s=None)[source]#

Electron numbers for each irreducible representation.

Args:
molan instance of Mole

To provide irrep_id, and spin-adapted basis

mo_coeff2D ndarray

Regular orbital coefficients, without grouping for irreps

mo_occ1D ndarray

Regular occupancy, without grouping for irreps

Returns:
irrep_nelecdict

The number of electrons for each irrep {‘ir_name’:int,…}.

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', symmetry=True, verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
-76.016789472074251
>>> scf.hf_symm.get_irrep_nelec(mol, mf.mo_coeff, mf.mo_occ)
{'A1': 6, 'A2': 0, 'B1': 2, 'B2': 2}
get_occ(mo_energy=None, mo_coeff=None)[source]#

We assumed mo_energy are grouped by symmetry irreps, (see function self.eig). The orbitals are sorted after SCF.

get_orbsym(mo_coeff=None, s=None)[source]#
property orbsym#
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.scf.ghf_symm.analyze(mf, verbose=5, with_meta_lowdin=True, **kwargs)[source]#
pyscf.scf.ghf_symm.canonicalize(mf, mo_coeff, mo_occ, fock=None)[source]#

Canonicalization diagonalizes the UHF Fock matrix in occupied, virtual subspaces separatedly (without change occupancy).

pyscf.scf.ghf_symm.get_orbsym(mol, mo_coeff, s=None, check=False, symm_orb=None, irrep_id=None)[source]#

pyscf.scf.hf module#

Hartree-Fock

pyscf.scf.hf.Kgwh(Ei, Ej, updated_rule=False)[source]#

Computes the generalized Wolfsberg-Helmholtz parameter

class pyscf.scf.hf.RHF(mol)[source]#

Bases: SCF

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
CASCI(*args, **kwargs)#
CASSCF(*args, **kwargs)#

CASCI/CASSCF

Args:
mf_or_molSCF object or Mole object

SCF or Mole to define the problem size.

ncasint

Number of active orbitals.

nelecasint or a pair of int

Number of electrons in active space.

Kwargs:
ncoreint

Number of doubly occupied core orbitals. If not presented, this parameter can be automatically determined.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose.

max_memoryfloat or int

Allowed memory in MB. Default value equals to Mole.max_memory.

ncasint

Active space size.

nelecastuple of int

Active (nelec_alpha, nelec_beta)

ncoreint or tuple of int

Core electron number. In UHF-CASSCF, it’s a tuple to indicate the different core electron numbers.

natorbbool

Whether to transform natural orbitals in active space. Note: when CASCI/CASSCF are combined with DMRG solver or selected CI solver, enabling this parameter may slightly change the total energy. False by default.

canonicalizationbool

Whether to canonicalize orbitals in core and external space against the general Fock matrix. The orbitals in active space are NOT transformed by default. To get the natural orbitals in active space, the attribute .natorb needs to be enabled. True by default.

sorting_mo_energybool

Whether to sort the orbitals based on the diagonal elements of the general Fock matrix. Default is False.

fcisolveran instance of FCISolver

The pyscf.fci module provides several FCISolver for different scenario. Generally, fci.direct_spin1.FCISolver can be used for all RHF-CASSCF. However, a proper FCISolver can provide better performance and better numerical stability. One can either use fci.solver() function to pick the FCISolver by the program or manually assigen the FCISolver to this attribute, e.g.

>>> from pyscf import fci
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> mc.fcisolver = fci.solver(mol, singlet=True)
>>> mc.fcisolver = fci.direct_spin1.FCISolver(mol)

You can control FCISolver by setting e.g.:

>>> mc.fcisolver.max_cycle = 30
>>> mc.fcisolver.conv_tol = 1e-7

For more details of the parameter for FCISolver, See fci.

Saved results

e_totfloat

Total MCSCF energy (electronic energy plus nuclear repulsion)

e_casfloat

CAS space FCI energy

cindarray

CAS space FCI coefficients

mo_coeffndarray

When canonicalization is specified, the orbitals are canonical orbitals which make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is specified, the active segment of the mo_coeff is natural orbitals.

mo_energyndarray

Diagonal elements of general Fock matrix (in mo_coeff representation).

mo_occndarray

Occupation numbers of natural orbitals if natorb is specified.

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASCI(mf, 6, 6)
>>> mc.kernel()[0]
-108.980200816243354

Extra attributes for CASSCF:

conv_tolfloat

Converge threshold. Default is 1e-7

conv_tol_gradfloat

Converge threshold for CI gradients and orbital rotation gradients. If not specified, it is set to sqrt(conv_tol).

max_stepsizefloat

The step size for orbital rotation. Small step (0.005 - 0.05) is prefered. Default is 0.02.

max_cycle_macroint

Max number of macro iterations. Default is 50.

max_cycle_microint

Max number of micro iterations in each macro iteration. Depending on systems, increasing this value might reduce the total macro iterations. Generally, 2 - 5 steps should be enough. Default is 4.

small_rot_tolfloat

Threshold for orbital rotation to be considered small. If the largest orbital rotation is smaller than this value, the CI solver will restart from the previous iteration if supported. Default is 0.01

ah_level_shiftfloat, for AH solver.

Level shift for the Davidson diagonalization in AH solver. Default is 1e-8.

ah_conv_tolfloat, for AH solver.

converge threshold for AH solver. Default is 1e-12.

ah_max_cyclefloat, for AH solver.

Max number of iterations allowd in AH solver. Default is 30.

ah_lindepfloat, for AH solver.

Linear dependence threshold for AH solver. Default is 1e-14.

ah_start_tolflat, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 2.5.

ah_start_cycleint, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 3.

ah_conv_tol, ah_max_cycle, ah_lindep, ah_start_tol and ah_start_cycle can affect the accuracy and performance of CASSCF solver. Lower ah_conv_tol and ah_lindep might improve the accuracy of CASSCF optimization, but decrease the performance.

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.conv_tol = 1e-10
>>> mc.ah_conv_tol = 1e-5
>>> mc.kernel()[0]
-109.044401898486001
>>> mc.ah_conv_tol = 1e-10
>>> mc.kernel()[0]
-109.044401887945668
chkfilestr

Checkpoint file to save the intermediate orbitals during the CASSCF optimization. Default is the checkpoint file of mean field object.

ci_response_spaceint

subspace size to solve the CI vector response. Default is 3.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

scale_restorationfloat

When a step of orbital rotation moves out of trust region, the orbital optimization will be restored to previous state and the step size of the orbital rotation needs to be reduced. scale_restoration controls how much to scale down the step size.

Saved results

e_totfloat

Total MCSCF energy (electronic energy plus nuclear repulsion)

e_casfloat

CAS space FCI energy

cindarray

CAS space FCI coefficients

mo_coeffndarray

Optimized CASSCF orbitals coefficients. When canonicalization is specified, the returned orbitals make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is enabled, the active segment of mo_coeff is transformed to natural orbitals.

mo_energyndarray

Diagonal elements of general Fock matrix (in mo_coeff representation).

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.kernel()[0]
-109.044401882238134
CISD(*args, **kwargs)#
DFMP2(*args, **kwargs)#
Gradients(*args, **kwargs)#

Non-relativistic restricted Hartree-Fock gradients

MP2(*args, **kwargs)#

restricted MP2 with canonical HF and non-canonical HF reference

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default value equals to Mole.max_memory

conv_tolfloat

For non-canonical MP2, converge threshold for MP2 correlation energy. Default value is 1e-7.

conv_tol_normtfloat

For non-canonical MP2, converge threshold for norm(t2). Default value is 1e-5.

max_cycleint

For non-canonical MP2, max number of MP2 iterations. Default value is 50.

diis_spaceint

For non-canonical MP2, DIIS space size in MP2 iterations. Default is 6.

level_shiftfloat

A shift on virtual orbital energies to stabilize the MP2 iterations.

frozenint or list

If integer is given, the inner-most orbitals are excluded from MP2 amplitudes. Given the orbital indices (0-based) in a list, both occupied and virtual orbitals can be frozen in MP2 calculation.

>>> mol = gto.M(atom = 'H 0 0 0; F 0 0 1.1', basis = 'ccpvdz')
>>> mf = scf.RHF(mol).run()
>>> # freeze 2 core orbitals
>>> pt = mp.MP2(mf).set(frozen = 2).run()
>>> # auto-generate the number of core orbitals to be frozen (1 in this case)
>>> pt = mp.MP2(mf).set_frozen().run()
>>> # freeze 2 core orbitals and 3 high lying unoccupied orbitals
>>> pt.set(frozen = [0,1,16,17,18]).run()

Saved results

e_corrfloat

MP2 correlation correction

e_corr_ss/osfloat

Same-spin and opposite-spin component of the MP2 correlation energy

e_totfloat

Total MP2 energy (HF + correlation)

t2 :

T amplitudes t2[i,j,a,b] (i,j in occ, a,b in virt)

TDA(*args, **kwargs)#

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.

TDHF(*args, **kwargs)#

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.

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.

convert_from_(mf)[source]#

Convert the input mean-field object to RHF/ROHF

gen_response(mo_coeff=None, mo_occ=None, singlet=None, hermi=0, max_memory=None)#

Generate a function to compute the product of RHF response function and RHF density matrices.

Kwargs:
singlet (None or boolean)If singlet is None, response function for

orbital hessian or CPHF will be generated. If singlet is boolean, it is used in TDDFT response kernel.

get_init_guess(mol=None, key='minao', **kwargs)[source]#
get_jk(mol=None, dm=None, hermi=1, with_j=True, with_k=True, omega=None)[source]#

Compute J, K matrices for all input density matrices

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
hermiint

Whether J, K matrix is hermitian

0 : not hermitian and not symmetric
1 : hermitian or symmetric
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

with_jboolean

Whether to compute J matrices

with_kboolean

Whether to compute K matrices

omegafloat

Parameter of range-separated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-separated Coulomb operator.

Returns:

Depending on the given dm, the function returns one J and one K matrix, or a list of J matrices and a list of K matrices, corresponding to the input density matrices.

Examples:

>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dms = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> j, k = scf.hf.get_jk(mol, dms, hermi=0)
>>> print(j.shape)
(3, 2, 2)
get_veff(mol=None, dm=None, dm_last=0, vhf_last=0, hermi=1)[source]#

Hartree-Fock potential matrix for the given density matrix

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
dm_lastndarray or a list of ndarrays or 0

The density matrix baseline. If not 0, this function computes the increment of HF potential w.r.t. the reference HF potential matrix.

vhf_lastndarray or a list of ndarrays or 0

The reference HF potential matrix.

hermiint

Whether J, K matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

Returns:

matrix Vhf = 2*J - K. Vhf can be a list matrices, corresponding to the input density matrices.

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dm0 = numpy.random.random((mol.nao_nr(),mol.nao_nr()))
>>> vhf0 = scf.hf.get_veff(mol, dm0, hermi=0)
>>> dm1 = numpy.random.random((mol.nao_nr(),mol.nao_nr()))
>>> vhf1 = scf.hf.get_veff(mol, dm1, hermi=0)
>>> vhf2 = scf.hf.get_veff(mol, dm1, dm_last=dm0, vhf_last=vhf0, hermi=0)
>>> numpy.allclose(vhf1, vhf2)
True
nuc_grad_method()[source]#

Hook to create object for analytical nuclear gradients.

spin_square(mo_coeff=None, s=None)[source]#

Spin square and multiplicity of RHF determinant

stability(internal=True, external=False, verbose=None, return_status=False, **kwargs)[source]#

RHF/RKS stability analysis.

See also pyscf.scf.stability.rhf_stability function.

Kwargs:
internalbool

Internal stability, within the RHF optimization space.

externalbool

External stability. Including the RHF -> UHF and real -> complex stability analysis.

return_status: bool

Whether to return stable_i and stable_e

Returns:

If return_status is False (default), the return value includes two set of orbitals, which are more close to the stable condition. The first corresponds to the internal stability and the second corresponds to the external stability.

Else, another two boolean variables (indicating current status: stable or unstable) are returned. The first corresponds to the internal stability and the second corresponds to the external stability.

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.

to_ks(xc='HF')[source]#

Convert to RKS object.

class pyscf.scf.hf.SCF(mol)[source]#

Bases: StreamObject

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
CCSD(frozen=None, mo_coeff=None, mo_occ=None)#

restricted CCSD

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default value equals to Mole.max_memory

conv_tolfloat

converge threshold. Default is 1e-7.

conv_tol_normtfloat

converge threshold for norm(t1,t2). Default is 1e-5.

max_cycleint

max number of iterations. Default is 50.

diis_spaceint

DIIS space size. Default is 6.

diis_start_cycleint

The step to start DIIS. Default is 0.

iterative_dampingfloat

The self consistent damping parameter.

directbool

AO-direct CCSD. Default is False.

async_iobool

Allow for asynchronous function execution. Default is True.

incore_completebool

Avoid all I/O (also for DIIS). Default is False.

level_shiftfloat

A shift on virtual orbital energies to stabilize the CCSD iteration

frozenint or list

If integer is given, the inner-most orbitals are frozen from CC amplitudes. Given the orbital indices (0-based) in a list, both occupied and virtual orbitals can be frozen in CC calculation.

>>> mol = gto.M(atom = 'H 0 0 0; F 0 0 1.1', basis = 'ccpvdz')
>>> mf = scf.RHF(mol).run()
>>> # freeze 2 core orbitals
>>> mycc = cc.CCSD(mf).set(frozen = 2).run()
>>> # auto-generate the number of core orbitals to be frozen (1 in this case)
>>> mycc = cc.CCSD(mf).set_frozen().run()
>>> # freeze 2 core orbitals and 3 high lying unoccupied orbitals
>>> mycc.set(frozen = [0,1,16,17,18]).run()
callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

Saved results:

convergedbool

Whether the CCSD iteration converged

e_corrfloat

CCSD correlation correction

e_totfloat

Total CCSD energy (HF + correlation)

t1, t2 :

T amplitudes t1[i,a], t2[i,j,a,b] (i,j in occ, a,b in virt)

l1, l2 :

Lambda amplitudes l1[i,a], l2[i,j,a,b] (i,j in occ, a,b in virt)

cyclesint

The number of iteration cycles performed

DIIS#

alias of CDIIS

QCISD(frozen=None, mo_coeff=None, mo_occ=None)#

restricted QCISD

QMMM(atoms_or_coords, charges, radii=None, unit=None)#

Embedding the one-electron (non-relativistic) potential generated by MM point charges into QM Hamiltonian.

The total energy includes the regular QM energy, the interaction between the nuclei in QM region and the MM charges, and the static Coulomb interaction between the electron density and the MM charges. It does not include the static Coulomb interactions of the MM point charges, the MM energy, the vdw interaction or other bonding/non-bonding effects between QM region and MM particles.

Args:

scf_method : a HF or DFT object

atoms_or_coords2D array, shape (N,3)

MM particle coordinates

charges1D array

MM particle charges

Kwargs:
radii1D array

The Gaussian charge distribution radii of MM atoms.

unitstr

Bohr, AU, Ang (case insensitive). Default is the same to mol.unit

Returns:

Same method object as the input scf_method with modified 1e Hamiltonian

Note:

1. if MM charge and X2C correction are used together, function mm_charge needs to be applied after X2C decoration (.x2c method), eg mf = mm_charge(scf.RHF(mol).x2c()), [(0.5,0.6,0.8)], [-0.5]). 2. Once mm_charge function is applied on the SCF object, it affects all the post-HF calculations eg MP2, CCSD, MCSCF etc

Examples:

>>> mol = gto.M(atom='H 0 0 0; F 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = mm_charge(dft.RKS(mol), [(0.5,0.6,0.8)], [-0.3])
>>> mf.kernel()
-101.940495711284
analyze(verbose=None, with_meta_lowdin=True, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis; Diople moment.

apply(fn, *args, **kwargs)[source]#

Apply the fn to rest arguments: return fn(*args, **kwargs). The return value of method set is the object itself. This allows a series of functions/methods to be executed in pipe.

as_scanner()#

Generating a scanner/solver for HF PES.

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

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the 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
>>> hf_scanner = scf.RHF(gto.Mole().set(verbose=0)).as_scanner()
>>> hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
-98.552190448277955
>>> hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
-98.414750424294368
build(mol=None)[source]#
callback = None#
canonicalize(mo_coeff, mo_occ, fock=None)#

Canonicalization diagonalizes the Fock matrix within occupied, open, virtual subspaces separatedly (without change occupancy).

check_convergence = None#
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_check = True#
conv_tol = 1e-09#
conv_tol_cpscf = 1e-08#
conv_tol_grad = None#
convert_from_(mf)[source]#

Convert the abinput mean-field object to the associated KS object.

damp = 0#
density_fit(auxbasis=None, with_df=None, only_dfj=False)[source]#
diis = True#
diis_damp = 0#
diis_file = None#
diis_space = 8#
diis_space_rollback = 0#
diis_start_cycle = 1#
dip_moment(mol=None, dm=None, unit='Debye', origin=None, verbose=3, **kwargs)[source]#

Dipole moment calculation

\[\begin{split}\mu_x = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|x|\mu) + \sum_A Q_A X_A\\ \mu_y = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|y|\mu) + \sum_A Q_A Y_A\\ \mu_z = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|z|\mu) + \sum_A Q_A Z_A\end{split}\]

where \(\mu_x, \mu_y, \mu_z\) are the x, y and z components of dipole moment

Args:

mol: an instance of Mole dm : a 2D ndarrays density matrices origin : optional; length 3 list, tuple, or 1D array

Location of the origin. By default, the point (0, 0, 0) is used.

Return:

A list: the dipole moment on x, y and z component

direct_scf = True#
direct_scf_tol = 1e-13#
disp = None#
do_disp(disp=None)#

Check whether to apply dispersion correction based on the xc attribute. If dispersion is allowed, return the DFTD3 disp version, such as d3bj, d3zero, d4.

dump_chk(envs_or_file)[source]#

Serialize the SCF object and save it to the specified chkfile.

Args:
envs_or_file:

If this argument is a file path, the serialized SCF object is saved to the file specified by this argument. If this attribute is a dict (created by locals()), the necessary variables are saved to the file specified by the attribute mf.chkfile.

dump_flags(verbose=None)[source]#
dump_scf_summary(verbose=5)#
eig(h, s)[source]#

Solver for generalized eigenvalue problem

\[HC = SCE\]
energy_elec(dm=None, h1e=None, vhf=None)#

Electronic part of Hartree-Fock energy, for given core hamiltonian and HF potential

… math:

E = \sum_{ij}h_{ij} \gamma_{ji}
  + \frac{1}{2}\sum_{ijkl} \gamma_{ji}\gamma_{lk} \langle ik||jl\rangle

Note this function has side effects which cause mf.scf_summary updated.

Args:

mf : an instance of SCF class

Kwargs:
dm2D ndarray

one-particle density matrix

h1e2D ndarray

Core hamiltonian

vhf2D ndarray

HF potential

Returns:

Hartree-Fock electronic energy and the Coulomb energy

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> dm = mf.make_rdm1()
>>> scf.hf.energy_elec(mf, dm)
(-1.5176090667746334, 0.60917167853723675)
>>> mf.energy_elec(dm)
(-1.5176090667746334, 0.60917167853723675)
energy_nuc()[source]#
energy_tot(dm=None, h1e=None, vhf=None)#

Total Hartree-Fock energy, electronic part plus nuclear repulsion See scf.hf.energy_elec() for the electron part

Note this function has side effects which cause mf.scf_summary updated.

from_chk(chkfile=None, project=None)[source]#

Read the HF results from checkpoint file, then project it to the basis defined by mol

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

Returns:

Density matrix, 2D ndarray

get_dispersion(disp=None, with_3body=None, verbose=None)#
get_fock(h1e=None, s1e=None, vhf=None, dm=None, cycle=-1, diis=None, diis_start_cycle=None, level_shift_factor=None, damp_factor=None, fock_last=None)#

F = h^{core} + V^{HF}

Special treatment (damping, DIIS, or level shift) will be applied to the Fock matrix if diis and cycle is specified (The two parameters are passed to get_fock function during the SCF iteration)

Kwargs:
h1e2D ndarray

Core hamiltonian

s1e2D ndarray

Overlap matrix, for DIIS

vhf2D ndarray

HF potential matrix

dm2D ndarray

Density matrix, for DIIS

cycleint

Then present SCF iteration step, for DIIS

diisan object of SCF.DIIS class

DIIS object to hold intermediate Fock and error vectors

diis_start_cycleint

The step to start DIIS. Default is 0.

level_shift_factorfloat or int

Level shift (in AU) for virtual space. Default is 0.

get_grad(mo_coeff, mo_occ, fock=None)[source]#

RHF orbital gradients

Args:
mo_coeff2D ndarray

Obital coefficients

mo_occ1D ndarray

Orbital occupancy

fock_ao2D ndarray

Fock matrix in AO representation

Returns:

Gradients in MO representation. It’s a num_occ*num_vir vector.

get_hcore(mol=None)[source]#
get_init_guess(mol=None, key='minao', **kwargs)[source]#
get_j(mol=None, dm=None, hermi=1, omega=None)[source]#

Compute J matrices for all input density matrices

get_jk(mol=None, dm=None, hermi=1, with_j=True, with_k=True, omega=None)[source]#

Compute J, K matrices for all input density matrices

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
hermiint

Whether J, K matrix is hermitian

0 : not hermitian and not symmetric
1 : hermitian or symmetric
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

with_jboolean

Whether to compute J matrices

with_kboolean

Whether to compute K matrices

omegafloat

Parameter of range-separated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-separated Coulomb operator.

Returns:

Depending on the given dm, the function returns one J and one K matrix, or a list of J matrices and a list of K matrices, corresponding to the input density matrices.

Examples:

>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dms = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> j, k = scf.hf.get_jk(mol, dms, hermi=0)
>>> print(j.shape)
(3, 2, 2)
get_k(mol=None, dm=None, hermi=1, omega=None)[source]#

Compute K matrices for all input density matrices

get_occ(mo_energy=None, mo_coeff=None)#

Label the occupancies for each orbital

Kwargs:
mo_energy1D ndarray

Obital energies

mo_coeff2D ndarray

Obital coefficients

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> energy = numpy.array([-10., -1., 1, -2., 0, -3])
>>> mf.get_occ(energy)
array([2, 2, 0, 2, 2, 2])
get_ovlp(mol=None)[source]#
get_veff(mol=None, dm=None, dm_last=0, vhf_last=0, hermi=1)[source]#

Hartree-Fock potential matrix for the given density matrix

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
dm_lastndarray or a list of ndarrays or 0

The density matrix baseline. If not 0, this function computes the increment of HF potential w.r.t. the reference HF potential matrix.

vhf_lastndarray or a list of ndarrays or 0

The reference HF potential matrix.

hermiint

Whether J, K matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

Returns:

matrix Vhf = 2*J - K. Vhf can be a list matrices, corresponding to the input density matrices.

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dm0 = numpy.random.random((mol.nao_nr(),mol.nao_nr()))
>>> vhf0 = scf.hf.get_veff(mol, dm0, hermi=0)
>>> dm1 = numpy.random.random((mol.nao_nr(),mol.nao_nr()))
>>> vhf1 = scf.hf.get_veff(mol, dm1, hermi=0)
>>> vhf2 = scf.hf.get_veff(mol, dm1, dm_last=dm0, vhf_last=vhf0, hermi=0)
>>> numpy.allclose(vhf1, vhf2)
True
init_direct_scf(mol=None)[source]#
init_guess = 'minao'#
init_guess_by_1e(mol=None)[source]#

Generate initial guess density matrix from core hamiltonian

Returns:

Density matrix, 2D ndarray

init_guess_by_atom(mol=None)[source]#

Generate initial guess density matrix from superposition of atomic HF density matrix. The atomic HF is occupancy averaged RHF

Returns:

Density matrix, 2D ndarray

init_guess_by_chkfile(chkfile=None, project=None)[source]#

Read the HF results from checkpoint file, then project it to the basis defined by mol

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

Returns:

Density matrix, 2D ndarray

init_guess_by_huckel(mol=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

Returns:

Density matrix, 2D ndarray

init_guess_by_minao(mol=None)[source]#

Generate initial guess density matrix based on ANO basis, then project the density matrix to the basis set defined by mol

Returns:

Density matrix, 2D ndarray

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> scf.hf.init_guess_by_minao(mol)
array([[ 0.94758917,  0.09227308],
       [ 0.09227308,  0.94758917]])
init_guess_by_mod_huckel(updated_rule, mol=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

In contrast to init_guess_by_huckel, this routine employs the updated GWH rule from doi:10.1021/ja00480a005 to form the guess.

Returns:

Density matrix, 2D ndarray

init_guess_by_sap(mol=None, **kwargs)[source]#

Generate initial guess density matrix from a superposition of atomic potentials (SAP), doi:10.1021/acs.jctc.8b01089. This is the Gaussian fit implementation, see doi:10.1063/5.0004046.

Args:
molMoleBase object

the molecule object for which the initial guess is evaluated

sap_basisdict

SAP basis in internal format (python dictionary)

Returns:
dm0ndarray

SAP initial guess density matrix

istype(type_code)[source]#

Checks if the object is an instance of the class specified by the type_code. type_code can be a class or a str. If the type_code is a class, it is equivalent to the Python built-in function isinstance. If the type_code is a str, it checks the type_code against the names of the object and all its parent classes.

kernel(dm0=None, **kwargs)#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
level_shift = 0#
make_rdm1(mo_coeff=None, mo_occ=None, **kwargs)#

One-particle density matrix in AO representation

Args:
mo_coeff2D ndarray

Orbital coefficients. Each column is one orbital.

mo_occ1D ndarray

Occupancy

Returns:

One-particle density matrix, 2D ndarray

make_rdm2(mo_coeff=None, mo_occ=None, **kwargs)#

Two-particle density matrix in AO representation

Args:
mo_coeff2D ndarray

Orbital coefficients. Each column is one orbital.

mo_occ1D ndarray

Occupancy

Returns:

Two-particle density matrix, 4D ndarray

max_cycle = 50#
mulliken_meta(mol=None, dm=None, verbose=5, pre_orth_method='ANO', s=None)[source]#

Mulliken population analysis, based on meta-Lowdin AOs. In the meta-lowdin, the AOs are grouped in three sets: core, valence and Rydberg, the orthogonalization are carried out within each subsets.

Args:

mol : an instance of Mole

dmndarray or 2-item list of ndarray

Density matrix. ROHF dm is a 2-item list of 2D array

Kwargs:

verbose : int or instance of lib.logger.Logger

pre_orth_methodstr

Pre-orthogonalization, which localized GTOs for each atom. To obtain the occupied and unoccupied atomic shells, there are three methods

‘ano’ : Project GTOs to ANO basis
‘minao’ : Project GTOs to MINAO basis
‘scf’ : Symmetry-averaged fractional occupation atomic RHF
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

mulliken_pop(mol=None, dm=None, s=None, verbose=5)[source]#

Mulliken population analysis

\[M_{ij} = D_{ij} S_{ji}\]

Mulliken charges

\[\delta_i = \sum_j M_{ij}\]
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

mulliken_pop_meta_lowdin_ao(*args, **kwargs)#

Mulliken population analysis, based on meta-Lowdin AOs. In the meta-lowdin, the AOs are grouped in three sets: core, valence and Rydberg, the orthogonalization are carried out within each subsets.

Args:

mol : an instance of Mole

dmndarray or 2-item list of ndarray

Density matrix. ROHF dm is a 2-item list of 2D array

Kwargs:

verbose : int or instance of lib.logger.Logger

pre_orth_methodstr

Pre-orthogonalization, which localized GTOs for each atom. To obtain the occupied and unoccupied atomic shells, there are three methods

‘ano’ : Project GTOs to ANO basis
‘minao’ : Project GTOs to MINAO basis
‘scf’ : Symmetry-averaged fractional occupation atomic RHF
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

newton()[source]#

Create an SOSCF object based on the mean-field object

nuc_grad_method()[source]#

Hook to create object for analytical nuclear gradients.

property opt#
pop(*args, **kwargs)[source]#

Mulliken population analysis, based on meta-Lowdin AOs. In the meta-lowdin, the AOs are grouped in three sets: core, valence and Rydberg, the orthogonalization are carried out within each subsets.

Args:

mol : an instance of Mole

dmndarray or 2-item list of ndarray

Density matrix. ROHF dm is a 2-item list of 2D array

Kwargs:

verbose : int or instance of lib.logger.Logger

pre_orth_methodstr

Pre-orthogonalization, which localized GTOs for each atom. To obtain the occupied and unoccupied atomic shells, there are three methods

‘ano’ : Project GTOs to ANO basis
‘minao’ : Project GTOs to MINAO basis
‘scf’ : Symmetry-averaged fractional occupation atomic RHF
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

quad_moment(mol=None, dm=None, unit='DebyeAngstrom', origin=None, verbose=3, **kwargs)[source]#

Calculates traceless quadrupole moment tensor.

The traceless quadrupole tensor is given by

\[\begin{split}Q_{ij} &= - \frac{1}{2} \sum_{\mu \nu} P_{\mu \nu} \left[ 3 (\nu | r_i r_j | \mu) - \delta_{ij} (\nu | r^2 | \mu) \right] \\ &+ \frac{1}{2} \sum_A Q_A \left( R_{iA} R_{jA} - \delta_{ij} \|\mathbf{R}_A\|^2 \right).\end{split}\]

If the molecule has a dipole, the quadrupole moment depends on the location of the origin. By default, the origin is taken to be (0, 0, 0), but it can be set manually via the keyword argument origin.

Args:

mol: an instance of Mole dm : a 2D ndarrays density matrices origin : optional; length 3 list, tuple, or 1D array

Location of the origin. By default, it is (0, 0, 0).

Return:

Traceless quadrupole tensor, 2D ndarray.

remove_soscf()[source]#

Remove the SOSCF decorator

reset(mol=None)[source]#

Reset mol and relevant attributes associated to the old mol object

sap_basis = 'sapgrasplarge'#
scf(dm0=None, **kwargs)[source]#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
sfx2c1e()[source]#
stability()[source]#
to_ghf()[source]#

Convert the input mean-field object to a GHF object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_gks(xc='HF')[source]#

Convert the input mean-field object to a GKS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_gpu()[source]#

Converts to the object with GPU support.

to_ks(xc='HF')[source]#

Convert the input mean-field object to the associated KS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_rhf()[source]#

Convert the input mean-field object to a RHF/ROHF object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_rks(xc='HF')[source]#

Convert the input mean-field object to a RKS/ROKS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_uhf()[source]#

Convert the input mean-field object to a UHF object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

to_uks(xc='HF')[source]#

Convert the input mean-field object to a UKS object.

Note this conversion only changes the class of the mean-field object. The total energy and wave-function are the same as them in the input mean-field object.

update(chkfile=None)#

Read attributes from the chkfile then replace the attributes of current object. It’s an alias of function update_from_chk_.

update_(chkfile=None)[source]#

Read attributes from the chkfile then replace the attributes of current object. It’s an alias of function update_from_chk_.

update_from_chk(chkfile=None)#

Read attributes from the chkfile then replace the attributes of current object. It’s an alias of function update_from_chk_.

update_from_chk_(chkfile=None)#

Read attributes from the chkfile then replace the attributes of current object. It’s an alias of function update_from_chk_.

x2c()#
x2c1e()#
class pyscf.scf.hf.SCF_Scanner(mf_obj)[source]#

Bases: SinglePointScanner

pyscf.scf.hf.analyze(mf, verbose=5, with_meta_lowdin=True, origin=None, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis; Diople moment.

pyscf.scf.hf.as_scanner(mf)[source]#

Generating a scanner/solver for HF PES.

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

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the 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
>>> hf_scanner = scf.RHF(gto.Mole().set(verbose=0)).as_scanner()
>>> hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
-98.552190448277955
>>> hf_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
-98.414750424294368
pyscf.scf.hf.canonicalize(mf, mo_coeff, mo_occ, fock=None)[source]#

Canonicalization diagonalizes the Fock matrix within occupied, open, virtual subspaces separatedly (without change occupancy).

pyscf.scf.hf.damping(f, f_prev, factor)[source]#
pyscf.scf.hf.dip_moment(mol, dm, unit='Debye', origin=None, verbose=3, **kwargs)[source]#

Dipole moment calculation

\[\begin{split}\mu_x = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|x|\mu) + \sum_A Q_A X_A\\ \mu_y = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|y|\mu) + \sum_A Q_A Y_A\\ \mu_z = -\sum_{\mu}\sum_{\nu} P_{\mu\nu}(\nu|z|\mu) + \sum_A Q_A Z_A\end{split}\]

where \(\mu_x, \mu_y, \mu_z\) are the x, y and z components of dipole moment

Args:

mol: an instance of Mole dm : a 2D ndarrays density matrices origin : optional; length 3 list, tuple, or 1D array

Location of the origin. By default, the point (0, 0, 0) is used.

Return:

A list: the dipole moment on x, y and z component

pyscf.scf.hf.dot_eri_dm(eri, dm, hermi=0, with_j=True, with_k=True)[source]#

Compute J, K matrices in terms of the given 2-electron integrals and density matrix:

J ~ numpy.einsum(‘pqrs,qp->rs’, eri, dm) K ~ numpy.einsum(‘pqrs,qr->ps’, eri, dm)

Args:
erindarray

8-fold or 4-fold ERIs or complex integral array with N^4 elements (N is the number of orbitals)

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
hermiint

Whether J, K matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
Returns:

Depending on the given dm, the function returns one J and one K matrix, or a list of J matrices and a list of K matrices, corresponding to the input density matrices.

Examples:

>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> eri = _vhf.int2e_sph(mol._atm, mol._bas, mol._env)
>>> dms = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> j, k = scf.hf.dot_eri_dm(eri, dms, hermi=0)
>>> print(j.shape)
(3, 2, 2)
pyscf.scf.hf.dump_scf_summary(mf, verbose=5)[source]#
pyscf.scf.hf.eig(h, s)[source]#

Solver for generalized eigenvalue problem

\[HC = SCE\]
pyscf.scf.hf.energy_elec(mf, dm=None, h1e=None, vhf=None)[source]#

Electronic part of Hartree-Fock energy, for given core hamiltonian and HF potential

… math:

E = \sum_{ij}h_{ij} \gamma_{ji}
  + \frac{1}{2}\sum_{ijkl} \gamma_{ji}\gamma_{lk} \langle ik||jl\rangle

Note this function has side effects which cause mf.scf_summary updated.

Args:

mf : an instance of SCF class

Kwargs:
dm2D ndarray

one-particle density matrix

h1e2D ndarray

Core hamiltonian

vhf2D ndarray

HF potential

Returns:

Hartree-Fock electronic energy and the Coulomb energy

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> dm = mf.make_rdm1()
>>> scf.hf.energy_elec(mf, dm)
(-1.5176090667746334, 0.60917167853723675)
>>> mf.energy_elec(dm)
(-1.5176090667746334, 0.60917167853723675)
pyscf.scf.hf.energy_tot(mf, dm=None, h1e=None, vhf=None)[source]#

Total Hartree-Fock energy, electronic part plus nuclear repulsion See scf.hf.energy_elec() for the electron part

Note this function has side effects which cause mf.scf_summary updated.

pyscf.scf.hf.get_fock(mf, h1e=None, s1e=None, vhf=None, dm=None, cycle=-1, diis=None, diis_start_cycle=None, level_shift_factor=None, damp_factor=None, fock_last=None)[source]#

F = h^{core} + V^{HF}

Special treatment (damping, DIIS, or level shift) will be applied to the Fock matrix if diis and cycle is specified (The two parameters are passed to get_fock function during the SCF iteration)

Kwargs:
h1e2D ndarray

Core hamiltonian

s1e2D ndarray

Overlap matrix, for DIIS

vhf2D ndarray

HF potential matrix

dm2D ndarray

Density matrix, for DIIS

cycleint

Then present SCF iteration step, for DIIS

diisan object of SCF.DIIS class

DIIS object to hold intermediate Fock and error vectors

diis_start_cycleint

The step to start DIIS. Default is 0.

level_shift_factorfloat or int

Level shift (in AU) for virtual space. Default is 0.

pyscf.scf.hf.get_grad(mo_coeff, mo_occ, fock_ao)[source]#

RHF orbital gradients

Args:
mo_coeff2D ndarray

Obital coefficients

mo_occ1D ndarray

Orbital occupancy

fock_ao2D ndarray

Fock matrix in AO representation

Returns:

Gradients in MO representation. It’s a num_occ*num_vir vector.

pyscf.scf.hf.get_hcore(mol)[source]#

Core Hamiltonian

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> scf.hf.get_hcore(mol)
array([[-0.93767904, -0.59316327],
       [-0.59316327, -0.93767904]])
pyscf.scf.hf.get_init_guess(mol, key='minao', **kwargs)[source]#

Generate density matrix for initial guess

Kwargs:
keystr

One of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’.

pyscf.scf.hf.get_jk(mol, dm, hermi=1, vhfopt=None, with_j=True, with_k=True, omega=None)[source]#

Compute J, K matrices for all input density matrices

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
hermiint

Whether J, K matrix is hermitian

0 : not hermitian and not symmetric
1 : hermitian or symmetric
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

with_jboolean

Whether to compute J matrices

with_kboolean

Whether to compute K matrices

omegafloat

Parameter of range-separated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-separated Coulomb operator.

Returns:

Depending on the given dm, the function returns one J and one K matrix, or a list of J matrices and a list of K matrices, corresponding to the input density matrices.

Examples:

>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dms = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> j, k = scf.hf.get_jk(mol, dms, hermi=0)
>>> print(j.shape)
(3, 2, 2)
pyscf.scf.hf.get_occ(mf, mo_energy=None, mo_coeff=None)[source]#

Label the occupancies for each orbital

Kwargs:
mo_energy1D ndarray

Obital energies

mo_coeff2D ndarray

Obital coefficients

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> energy = numpy.array([-10., -1., 1, -2., 0, -3])
>>> mf.get_occ(energy)
array([2, 2, 0, 2, 2, 2])
pyscf.scf.hf.get_ovlp(mol)[source]#

Overlap matrix

pyscf.scf.hf.get_veff(mol, dm, dm_last=None, vhf_last=None, hermi=1, vhfopt=None)[source]#

Hartree-Fock potential matrix for the given density matrix

Args:

mol : an instance of Mole

dmndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
dm_lastndarray or a list of ndarrays or 0

The density matrix baseline. If not 0, this function computes the increment of HF potential w.r.t. the reference HF potential matrix.

vhf_lastndarray or a list of ndarrays or 0

The reference HF potential matrix.

hermiint

Whether J, K matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

Returns:

matrix Vhf = 2*J - K. Vhf can be a list matrices, corresponding to the input density matrices.

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> from pyscf.scf import _vhf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dm0 = numpy.random.random((mol.nao_nr(),mol.nao_nr()))
>>> vhf0 = scf.hf.get_veff(mol, dm0, hermi=0)
>>> dm1 = numpy.random.random((mol.nao_nr(),mol.nao_nr()))
>>> vhf1 = scf.hf.get_veff(mol, dm1, hermi=0)
>>> vhf2 = scf.hf.get_veff(mol, dm1, dm_last=dm0, vhf_last=vhf0, hermi=0)
>>> numpy.allclose(vhf1, vhf2)
True
pyscf.scf.hf.init_guess_by_1e(mol)[source]#

Generate initial guess density matrix from core hamiltonian

Returns:

Density matrix, 2D ndarray

pyscf.scf.hf.init_guess_by_atom(mol)[source]#

Generate initial guess density matrix from superposition of atomic HF density matrix. The atomic HF is occupancy averaged RHF

Returns:

Density matrix, 2D ndarray

pyscf.scf.hf.init_guess_by_chkfile(mol, chkfile_name, project=None)[source]#

Read the HF results from checkpoint file, then project it to the basis defined by mol

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

Returns:

Density matrix, 2D ndarray

pyscf.scf.hf.init_guess_by_huckel(mol)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

Returns:

Density matrix, 2D ndarray

pyscf.scf.hf.init_guess_by_minao(mol)[source]#

Generate initial guess density matrix based on ANO basis, then project the density matrix to the basis set defined by mol

Returns:

Density matrix, 2D ndarray

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> scf.hf.init_guess_by_minao(mol)
array([[ 0.94758917,  0.09227308],
       [ 0.09227308,  0.94758917]])
pyscf.scf.hf.init_guess_by_mod_huckel(mol)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

In contrast to init_guess_by_huckel, this routine employs the updated GWH rule from doi:10.1021/ja00480a005 to form the guess.

Returns:

Density matrix, 2D ndarray

pyscf.scf.hf.init_guess_by_sap(mol, sap_basis, **kwargs)[source]#

Generate initial guess density matrix from a superposition of atomic potentials (SAP), doi:10.1021/acs.jctc.8b01089. This is the Gaussian fit implementation, see doi:10.1063/5.0004046.

Args:
molMoleBase object

the molecule object for which the initial guess is evaluated

sap_basisdict

SAP basis in internal format (python dictionary)

Returns:
dm0ndarray

SAP initial guess density matrix

pyscf.scf.hf.kernel(mf, conv_tol=1e-10, conv_tol_grad=None, dump_chk=True, dm0=None, callback=None, conv_check=True, **kwargs)[source]#

kernel: the SCF driver.

Args:
mfan instance of SCF class

mf object holds all parameters to control SCF. One can modify its member functions to change the behavior of SCF. The member functions which are called in kernel are

mf.get_init_guess
mf.get_hcore
mf.get_ovlp
mf.get_veff
mf.get_fock
mf.get_grad
mf.eig
mf.get_occ
mf.make_rdm1
mf.energy_tot
mf.dump_chk
Kwargs:
conv_tolfloat

converge threshold.

conv_tol_gradfloat

gradients converge threshold.

dump_chkbool

Whether to save SCF intermediate results in the checkpoint file

dm0ndarray

Initial guess density matrix. If not given (the default), the kernel takes the density matrix generated by mf.get_init_guess.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

sap_basisstr

SAP basis name

Returns:

A list : scf_conv, e_tot, mo_energy, mo_coeff, mo_occ

scf_convbool

True means SCF converged

e_totfloat

Hartree-Fock energy of last iteration

mo_energy1D float array

Orbital energies. Depending the eig function provided by mf object, the orbital energies may NOT be sorted.

mo_coeff2D array

Orbital coefficients.

mo_occ1D array

Orbital occupancies. The occupancies may NOT be sorted from large to small.

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> conv, e, mo_e, mo, mo_occ = scf.hf.kernel(scf.hf.SCF(mol), dm0=numpy.eye(mol.nao_nr()))
>>> print('conv = %s, E(HF) = %.12f' % (conv, e))
conv = True, E(HF) = -1.081170784378
pyscf.scf.hf.level_shift(s, d, f, factor)[source]#

Apply level shift \(\Delta\) to virtual orbitals

\begin{align} FC &= SCE \\ F &= F + SC \Lambda C^\dagger S \\ \Lambda_{ij} &= \begin{cases} \delta_{ij}\Delta & i \in \text{virtual} \\ 0 & \text{otherwise} \end{cases} \end{align}
Returns:

New Fock matrix, 2D ndarray

pyscf.scf.hf.make_rdm1(mo_coeff, mo_occ, **kwargs)[source]#

One-particle density matrix in AO representation

Args:
mo_coeff2D ndarray

Orbital coefficients. Each column is one orbital.

mo_occ1D ndarray

Occupancy

Returns:

One-particle density matrix, 2D ndarray

pyscf.scf.hf.make_rdm2(mo_coeff, mo_occ, **kwargs)[source]#

Two-particle density matrix in AO representation

Args:
mo_coeff2D ndarray

Orbital coefficients. Each column is one orbital.

mo_occ1D ndarray

Occupancy

Returns:

Two-particle density matrix, 4D ndarray

pyscf.scf.hf.make_sap(mol, sap_basis)[source]#

Superposition of atomic potentials (SAP) potential matrix

Args:
molMoleBase object

molecule for which SAP is computed

sap_basisdict

SAP basis

Returns:
Vsapndarray

SAP potential matrix

pyscf.scf.hf.mulliken_meta(mol, dm, verbose=5, pre_orth_method='ANO', s=None)[source]#

Mulliken population analysis, based on meta-Lowdin AOs. In the meta-lowdin, the AOs are grouped in three sets: core, valence and Rydberg, the orthogonalization are carried out within each subsets.

Args:

mol : an instance of Mole

dmndarray or 2-item list of ndarray

Density matrix. ROHF dm is a 2-item list of 2D array

Kwargs:

verbose : int or instance of lib.logger.Logger

pre_orth_methodstr

Pre-orthogonalization, which localized GTOs for each atom. To obtain the occupied and unoccupied atomic shells, there are three methods

‘ano’ : Project GTOs to ANO basis
‘minao’ : Project GTOs to MINAO basis
‘scf’ : Symmetry-averaged fractional occupation atomic RHF
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

pyscf.scf.hf.mulliken_pop(mol, dm, s=None, verbose=5)[source]#

Mulliken population analysis

\[M_{ij} = D_{ij} S_{ji}\]

Mulliken charges

\[\delta_i = \sum_j M_{ij}\]
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

pyscf.scf.hf.mulliken_pop_meta_lowdin_ao(mol, dm, verbose=5, pre_orth_method='ANO', s=None)#

Mulliken population analysis, based on meta-Lowdin AOs. In the meta-lowdin, the AOs are grouped in three sets: core, valence and Rydberg, the orthogonalization are carried out within each subsets.

Args:

mol : an instance of Mole

dmndarray or 2-item list of ndarray

Density matrix. ROHF dm is a 2-item list of 2D array

Kwargs:

verbose : int or instance of lib.logger.Logger

pre_orth_methodstr

Pre-orthogonalization, which localized GTOs for each atom. To obtain the occupied and unoccupied atomic shells, there are three methods

‘ano’ : Project GTOs to ANO basis
‘minao’ : Project GTOs to MINAO basis
‘scf’ : Symmetry-averaged fractional occupation atomic RHF
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

pyscf.scf.hf.pack_uniq_var(x, mo_occ)[source]#

Extract the unique variables from the full orbital-gradients (or orbital-rotation) matrix

pyscf.scf.hf.quad_moment(mol, dm, unit='DebyeAngstrom', origin=None, verbose=3, **kwargs)[source]#

Calculates traceless quadrupole moment tensor.

The traceless quadrupole tensor is given by

\[\begin{split}Q_{ij} &= - \frac{1}{2} \sum_{\mu \nu} P_{\mu \nu} \left[ 3 (\nu | r_i r_j | \mu) - \delta_{ij} (\nu | r^2 | \mu) \right] \\ &+ \frac{1}{2} \sum_A Q_A \left( R_{iA} R_{jA} - \delta_{ij} \|\mathbf{R}_A\|^2 \right).\end{split}\]

If the molecule has a dipole, the quadrupole moment depends on the location of the origin. By default, the origin is taken to be (0, 0, 0), but it can be set manually via the keyword argument origin.

Args:

mol: an instance of Mole dm : a 2D ndarrays density matrices origin : optional; length 3 list, tuple, or 1D array

Location of the origin. By default, it is (0, 0, 0).

Return:

Traceless quadrupole tensor, 2D ndarray.

pyscf.scf.hf.uniq_var_indices(mo_occ)[source]#

Indices of the unique variables for the orbital-gradients (or orbital-rotation) matrix.

pyscf.scf.hf.unpack_uniq_var(dx, mo_occ)[source]#

Fill the full orbital-gradients (or orbital-rotation) matrix with the unique variables.

pyscf.scf.hf_symm module#

Non-relativistic restricted Hartree-Fock with point group symmetry.

The symmetry are not handled in a separate data structure. Note that during the SCF iteration, the orbitals are grouped in terms of symmetry irreps. But the orbitals in the result are sorted based on the orbital energies. Function symm.label_orb_symm can be used to detect the symmetry of the molecular orbitals.

class pyscf.scf.hf_symm.HF1e(mol)[source]#

Bases: SymAdaptedROHF

scf(*args)#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
pyscf.scf.hf_symm.RHF#

alias of SymAdaptedRHF

pyscf.scf.hf_symm.ROHF#

alias of SymAdaptedROHF

class pyscf.scf.hf_symm.SymAdaptedRHF(mol)[source]#

Bases: RHF

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
Attributes for symmetry allowed RHF:
irrep_nelecdict

Specify the number of electrons for particular irrep {‘ir_name’:int,…}. For the irreps not listed in this dict, the program will choose the occupancy based on the orbital energies.

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', symmetry=True, verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
-76.016789472074251
>>> mf.get_irrep_nelec()
{'A1': 6, 'A2': 0, 'B1': 2, 'B2': 2}
>>> mf.irrep_nelec = {'A2': 2}
>>> mf.scf()
-72.768201804695622
>>> mf.get_irrep_nelec()
{'A1': 6, 'A2': 2, 'B1': 2, 'B2': 0}
CASCI(*args, **kwargs)#
CASSCF(*args, **kwargs)#

CASCI/CASSCF

Args:
mf_or_molSCF object or Mole object

SCF or Mole to define the problem size.

ncasint

Number of active orbitals.

nelecasint or a pair of int

Number of electrons in active space.

Kwargs:
ncoreint

Number of doubly occupied core orbitals. If not presented, this parameter can be automatically determined.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose.

max_memoryfloat or int

Allowed memory in MB. Default value equals to Mole.max_memory.

ncasint

Active space size.

nelecastuple of int

Active (nelec_alpha, nelec_beta)

ncoreint or tuple of int

Core electron number. In UHF-CASSCF, it’s a tuple to indicate the different core electron numbers.

natorbbool

Whether to transform natural orbitals in active space. Note: when CASCI/CASSCF are combined with DMRG solver or selected CI solver, enabling this parameter may slightly change the total energy. False by default.

canonicalizationbool

Whether to canonicalize orbitals in core and external space against the general Fock matrix. The orbitals in active space are NOT transformed by default. To get the natural orbitals in active space, the attribute .natorb needs to be enabled. True by default.

sorting_mo_energybool

Whether to sort the orbitals based on the diagonal elements of the general Fock matrix. Default is False.

fcisolveran instance of FCISolver

The pyscf.fci module provides several FCISolver for different scenario. Generally, fci.direct_spin1.FCISolver can be used for all RHF-CASSCF. However, a proper FCISolver can provide better performance and better numerical stability. One can either use fci.solver() function to pick the FCISolver by the program or manually assigen the FCISolver to this attribute, e.g.

>>> from pyscf import fci
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> mc.fcisolver = fci.solver(mol, singlet=True)
>>> mc.fcisolver = fci.direct_spin1.FCISolver(mol)

You can control FCISolver by setting e.g.:

>>> mc.fcisolver.max_cycle = 30
>>> mc.fcisolver.conv_tol = 1e-7

For more details of the parameter for FCISolver, See fci.

Saved results

e_totfloat

Total MCSCF energy (electronic energy plus nuclear repulsion)

e_casfloat

CAS space FCI energy

cindarray

CAS space FCI coefficients

mo_coeffndarray

When canonicalization is specified, the orbitals are canonical orbitals which make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is specified, the active segment of the mo_coeff is natural orbitals.

mo_energyndarray

Diagonal elements of general Fock matrix (in mo_coeff representation).

mo_occndarray

Occupation numbers of natural orbitals if natorb is specified.

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASCI(mf, 6, 6)
>>> mc.kernel()[0]
-108.980200816243354

Extra attributes for CASSCF:

conv_tolfloat

Converge threshold. Default is 1e-7

conv_tol_gradfloat

Converge threshold for CI gradients and orbital rotation gradients. If not specified, it is set to sqrt(conv_tol).

max_stepsizefloat

The step size for orbital rotation. Small step (0.005 - 0.05) is prefered. Default is 0.02.

max_cycle_macroint

Max number of macro iterations. Default is 50.

max_cycle_microint

Max number of micro iterations in each macro iteration. Depending on systems, increasing this value might reduce the total macro iterations. Generally, 2 - 5 steps should be enough. Default is 4.

small_rot_tolfloat

Threshold for orbital rotation to be considered small. If the largest orbital rotation is smaller than this value, the CI solver will restart from the previous iteration if supported. Default is 0.01

ah_level_shiftfloat, for AH solver.

Level shift for the Davidson diagonalization in AH solver. Default is 1e-8.

ah_conv_tolfloat, for AH solver.

converge threshold for AH solver. Default is 1e-12.

ah_max_cyclefloat, for AH solver.

Max number of iterations allowd in AH solver. Default is 30.

ah_lindepfloat, for AH solver.

Linear dependence threshold for AH solver. Default is 1e-14.

ah_start_tolflat, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 2.5.

ah_start_cycleint, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 3.

ah_conv_tol, ah_max_cycle, ah_lindep, ah_start_tol and ah_start_cycle can affect the accuracy and performance of CASSCF solver. Lower ah_conv_tol and ah_lindep might improve the accuracy of CASSCF optimization, but decrease the performance.

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.conv_tol = 1e-10
>>> mc.ah_conv_tol = 1e-5
>>> mc.kernel()[0]
-109.044401898486001
>>> mc.ah_conv_tol = 1e-10
>>> mc.kernel()[0]
-109.044401887945668
chkfilestr

Checkpoint file to save the intermediate orbitals during the CASSCF optimization. Default is the checkpoint file of mean field object.

ci_response_spaceint

subspace size to solve the CI vector response. Default is 3.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

scale_restorationfloat

When a step of orbital rotation moves out of trust region, the orbital optimization will be restored to previous state and the step size of the orbital rotation needs to be reduced. scale_restoration controls how much to scale down the step size.

Saved results

e_totfloat

Total MCSCF energy (electronic energy plus nuclear repulsion)

e_casfloat

CAS space FCI energy

cindarray

CAS space FCI coefficients

mo_coeffndarray

Optimized CASSCF orbitals coefficients. When canonicalization is specified, the returned orbitals make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is enabled, the active segment of mo_coeff is transformed to natural orbitals.

mo_energyndarray

Diagonal elements of general Fock matrix (in mo_coeff representation).

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.kernel()[0]
-109.044401882238134
analyze(verbose=None, with_meta_lowdin=True, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Occupancy for each irreps; Mulliken population analysis

build(mol=None)[source]#
canonicalize(mo_coeff, mo_occ, fock=None)#

Canonicalization diagonalizes the Fock matrix in occupied, open, virtual subspaces separatedly (without change occupancy).

eig(h, s, symm_orb=None, irrep_id=None)#

Solve generalized eigenvalue problem, for each irrep. The eigenvalues and eigenvectors are not sorted to ascending order. Instead, they are grouped based on irreps.

get_grad(mo_coeff, mo_occ, fock=None)[source]#

RHF orbital gradients

Args:
mo_coeff2D ndarray

Obital coefficients

mo_occ1D ndarray

Orbital occupancy

fock_ao2D ndarray

Fock matrix in AO representation

Returns:

Gradients in MO representation. It’s a num_occ*num_vir vector.

get_irrep_nelec(mol=None, mo_coeff=None, mo_occ=None, s=None)[source]#

Electron numbers for each irreducible representation.

Args:
molan instance of Mole

To provide irrep_id, and spin-adapted basis

mo_coeff2D ndarray

Regular orbital coefficients, without grouping for irreps

mo_occ1D ndarray

Regular occupancy, without grouping for irreps

Returns:
irrep_nelecdict

The number of electrons for each irrep {‘ir_name’:int,…}.

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', symmetry=True, verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
-76.016789472074251
>>> scf.hf_symm.get_irrep_nelec(mol, mf.mo_coeff, mf.mo_occ)
{'A1': 6, 'A2': 0, 'B1': 2, 'B2': 2}
get_occ(mo_energy=None, mo_coeff=None)[source]#

We assumed mo_energy are grouped by symmetry irreps, (see function self.eig). The orbitals are sorted after SCF.

get_orbsym(mo_coeff=None, s=None)[source]#
get_wfnsym(mo_coeff=None, mo_occ=None, orbsym=None)#
property orbsym#
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.

property wfnsym#
class pyscf.scf.hf_symm.SymAdaptedROHF(mol)[source]#

Bases: ROHF

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
Attributes for symmetry allowed ROHF:
irrep_nelecdict

Specify the number of alpha/beta electrons for particular irrep {‘ir_name’:(int,int), …}. For the irreps not listed in these dicts, the program will choose the occupancy based on the orbital energies.

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', symmetry=True, charge=1, spin=1, verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
-75.619358861084052
>>> mf.get_irrep_nelec()
{'A1': (3, 3), 'A2': (0, 0), 'B1': (1, 1), 'B2': (1, 0)}
>>> mf.irrep_nelec = {'B1': (1, 0)}
>>> mf.scf()
-75.425669486776457
>>> mf.get_irrep_nelec()
{'A1': (3, 3), 'A2': (0, 0), 'B1': (1, 0), 'B2': (1, 1)}
CASCI(*args, **kwargs)#
CASSCF(*args, **kwargs)#

CASCI/CASSCF

Args:
mf_or_molSCF object or Mole object

SCF or Mole to define the problem size.

ncasint

Number of active orbitals.

nelecasint or a pair of int

Number of electrons in active space.

Kwargs:
ncoreint

Number of doubly occupied core orbitals. If not presented, this parameter can be automatically determined.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose.

max_memoryfloat or int

Allowed memory in MB. Default value equals to Mole.max_memory.

ncasint

Active space size.

nelecastuple of int

Active (nelec_alpha, nelec_beta)

ncoreint or tuple of int

Core electron number. In UHF-CASSCF, it’s a tuple to indicate the different core electron numbers.

natorbbool

Whether to transform natural orbitals in active space. Note: when CASCI/CASSCF are combined with DMRG solver or selected CI solver, enabling this parameter may slightly change the total energy. False by default.

canonicalizationbool

Whether to canonicalize orbitals in core and external space against the general Fock matrix. The orbitals in active space are NOT transformed by default. To get the natural orbitals in active space, the attribute .natorb needs to be enabled. True by default.

sorting_mo_energybool

Whether to sort the orbitals based on the diagonal elements of the general Fock matrix. Default is False.

fcisolveran instance of FCISolver

The pyscf.fci module provides several FCISolver for different scenario. Generally, fci.direct_spin1.FCISolver can be used for all RHF-CASSCF. However, a proper FCISolver can provide better performance and better numerical stability. One can either use fci.solver() function to pick the FCISolver by the program or manually assigen the FCISolver to this attribute, e.g.

>>> from pyscf import fci
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> mc.fcisolver = fci.solver(mol, singlet=True)
>>> mc.fcisolver = fci.direct_spin1.FCISolver(mol)

You can control FCISolver by setting e.g.:

>>> mc.fcisolver.max_cycle = 30
>>> mc.fcisolver.conv_tol = 1e-7

For more details of the parameter for FCISolver, See fci.

Saved results

e_totfloat

Total MCSCF energy (electronic energy plus nuclear repulsion)

e_casfloat

CAS space FCI energy

cindarray

CAS space FCI coefficients

mo_coeffndarray

When canonicalization is specified, the orbitals are canonical orbitals which make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is specified, the active segment of the mo_coeff is natural orbitals.

mo_energyndarray

Diagonal elements of general Fock matrix (in mo_coeff representation).

mo_occndarray

Occupation numbers of natural orbitals if natorb is specified.

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASCI(mf, 6, 6)
>>> mc.kernel()[0]
-108.980200816243354

Extra attributes for CASSCF:

conv_tolfloat

Converge threshold. Default is 1e-7

conv_tol_gradfloat

Converge threshold for CI gradients and orbital rotation gradients. If not specified, it is set to sqrt(conv_tol).

max_stepsizefloat

The step size for orbital rotation. Small step (0.005 - 0.05) is prefered. Default is 0.02.

max_cycle_macroint

Max number of macro iterations. Default is 50.

max_cycle_microint

Max number of micro iterations in each macro iteration. Depending on systems, increasing this value might reduce the total macro iterations. Generally, 2 - 5 steps should be enough. Default is 4.

small_rot_tolfloat

Threshold for orbital rotation to be considered small. If the largest orbital rotation is smaller than this value, the CI solver will restart from the previous iteration if supported. Default is 0.01

ah_level_shiftfloat, for AH solver.

Level shift for the Davidson diagonalization in AH solver. Default is 1e-8.

ah_conv_tolfloat, for AH solver.

converge threshold for AH solver. Default is 1e-12.

ah_max_cyclefloat, for AH solver.

Max number of iterations allowd in AH solver. Default is 30.

ah_lindepfloat, for AH solver.

Linear dependence threshold for AH solver. Default is 1e-14.

ah_start_tolflat, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 2.5.

ah_start_cycleint, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 3.

ah_conv_tol, ah_max_cycle, ah_lindep, ah_start_tol and ah_start_cycle can affect the accuracy and performance of CASSCF solver. Lower ah_conv_tol and ah_lindep might improve the accuracy of CASSCF optimization, but decrease the performance.

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.conv_tol = 1e-10
>>> mc.ah_conv_tol = 1e-5
>>> mc.kernel()[0]
-109.044401898486001
>>> mc.ah_conv_tol = 1e-10
>>> mc.kernel()[0]
-109.044401887945668
chkfilestr

Checkpoint file to save the intermediate orbitals during the CASSCF optimization. Default is the checkpoint file of mean field object.

ci_response_spaceint

subspace size to solve the CI vector response. Default is 3.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

scale_restorationfloat

When a step of orbital rotation moves out of trust region, the orbital optimization will be restored to previous state and the step size of the orbital rotation needs to be reduced. scale_restoration controls how much to scale down the step size.

Saved results

e_totfloat

Total MCSCF energy (electronic energy plus nuclear repulsion)

e_casfloat

CAS space FCI energy

cindarray

CAS space FCI coefficients

mo_coeffndarray

Optimized CASSCF orbitals coefficients. When canonicalization is specified, the returned orbitals make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is enabled, the active segment of mo_coeff is transformed to natural orbitals.

mo_energyndarray

Diagonal elements of general Fock matrix (in mo_coeff representation).

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.kernel()[0]
-109.044401882238134
TDA = None#
TDHF = None#
analyze(verbose=None, with_meta_lowdin=True, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis

build(mol=None)[source]#
canonicalize(mo_coeff, mo_occ, fock=None)[source]#

Canonicalization diagonalizes the Fock matrix in occupied, open, virtual subspaces separatedly (without change occupancy).

dump_flags(verbose=None)[source]#
eig(fock, s)[source]#

Solve generalized eigenvalue problem, for each irrep. The eigenvalues and eigenvectors are not sorted to ascending order. Instead, they are grouped based on irreps.

get_grad(mo_coeff, mo_occ, fock=None)[source]#

ROHF gradients is the off-diagonal block [co + cv + ov], where [ cc co cv ] [ oc oo ov ] [ vc vo vv ]

get_irrep_nelec(mol=None, mo_coeff=None, mo_occ=None, s=None)[source]#
get_occ(mo_energy=None, mo_coeff=None)[source]#

Label the occupancies for each orbital. NOTE the occupancies are not assigned based on the orbital energy ordering. The first N orbitals are assigned to be occupied orbitals.

Examples:

>>> mol = gto.M(atom='H 0 0 0; O 0 0 1.1', spin=1)
>>> mf = scf.hf.SCF(mol)
>>> energy = numpy.array([-10., -1., 1, -2., 0, -3])
>>> mf.get_occ(energy)
array([2, 2, 2, 2, 1, 0])
get_orbsym(mo_coeff=None, s=None)#
get_wfnsym(mo_coeff=None, mo_occ=None, orbsym=None)#
make_rdm1(mo_coeff=None, mo_occ=None, **kwargs)[source]#

One-particle density matrix. mo_occ is a 1D array, with occupancy 1 or 2.

property orbsym#
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.

property wfnsym#
pyscf.scf.hf_symm.analyze(mf, verbose=5, with_meta_lowdin=True, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Occupancy for each irreps; Mulliken population analysis

pyscf.scf.hf_symm.canonicalize(mf, mo_coeff, mo_occ, fock=None)[source]#

Canonicalization diagonalizes the Fock matrix in occupied, open, virtual subspaces separatedly (without change occupancy).

pyscf.scf.hf_symm.check_irrep_nelec(mol, irrep_nelec, nelec)[source]#
pyscf.scf.hf_symm.eig(mf, h, s, symm_orb=None, irrep_id=None)[source]#

Solve generalized eigenvalue problem, for each irrep. The eigenvalues and eigenvectors are not sorted to ascending order. Instead, they are grouped based on irreps.

pyscf.scf.hf_symm.get_irrep_nelec(mol, mo_coeff, mo_occ, s=None)[source]#

Electron numbers for each irreducible representation.

Args:
molan instance of Mole

To provide irrep_id, and spin-adapted basis

mo_coeff2D ndarray

Regular orbital coefficients, without grouping for irreps

mo_occ1D ndarray

Regular occupancy, without grouping for irreps

Returns:
irrep_nelecdict

The number of electrons for each irrep {‘ir_name’:int,…}.

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', symmetry=True, verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
-76.016789472074251
>>> scf.hf_symm.get_irrep_nelec(mol, mf.mo_coeff, mf.mo_occ)
{'A1': 6, 'A2': 0, 'B1': 2, 'B2': 2}
pyscf.scf.hf_symm.get_orbsym(mol, mo_coeff, s=None, check=False, symm_orb=None, irrep_id=None)[source]#
pyscf.scf.hf_symm.get_wfnsym(mf, mo_coeff=None, mo_occ=None, orbsym=None)[source]#
pyscf.scf.hf_symm.map_degeneracy(mo_energy, orbsym)[source]#

Find degeneracy correspondence for cylindrical symmetry

pyscf.scf.hf_symm.so2ao_mo_coeff(so, irrep_mo_coeff)[source]#

Transfer the basis of MO coefficients, from symmetry-adapted basis to AO basis

pyscf.scf.jk module#

General JK contraction function for * arbitrary integrals * 4 different molecules * multiple density matrices * arbitrary basis subset for the 4 indices

pyscf.scf.jk.get_jk(mols, dms, scripts=['ijkl,ji->kl'], intor='int2e_sph', aosym='s1', comp=None, hermi=0, shls_slice=None, verbose=2, vhfopt=None)[source]#

Compute J/K matrices for the given density matrix

Args:

mols : an instance of Mole or a list of Mole objects

dmsndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
hermiint

Whether the returned J (K) matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
intorstr

2-electron integral name. See getints() for the complete list of available 2-electron integral names

aosymint or str

Permutation symmetry for the AO integrals

4 or ‘4’ or ‘s4’: 4-fold symmetry (default)
‘2ij’ or ‘s2ij’ : symmetry between i, j in (ij|kl)
‘2kl’ or ‘s2kl’ : symmetry between k, l in (ij|kl)
1 or ‘1’ or ‘s1’: no symmetry
‘a4ij’ : 4-fold symmetry with anti-symmetry between i, j in (ij|kl)
‘a4kl’ : 4-fold symmetry with anti-symmetry between k, l in (ij|kl)
‘a2ij’ : anti-symmetry between i, j in (ij|kl)
‘a2kl’ : anti-symmetry between k, l in (ij|kl)
compint

Components of the integrals, e.g. cint2e_ip_sph has 3 components.

scriptsstring or a list of strings

Contraction description (following numpy.einsum convention) based on letters [ijkl]. Each script will be one-to-one applied to each entry of dms. So it must have the same number of elements as the dms, len(scripts) == len(dms).

shls_slice8-element list

(ish_start, ish_end, jsh_start, jsh_end, ksh_start, ksh_end, lsh_start, lsh_end)

Returns:

Depending on the number of density matrices, the function returns one J/K matrix or a list of J/K matrices (the same number of entries as the input dms). Each JK matrices may be a 2D array or 3D array if the AO integral has multiple components.

Examples:

>>> from pyscf import gto
>>> mol = gto.M(atom='H 0 -.5 0; H 0 .5 0', basis='cc-pvdz')
>>> nao = mol.nao_nr()
>>> dm = numpy.random.random((nao,nao))
>>> # Default, Coulomb matrix
>>> vj = get_jk(mol, dm)
>>> # Coulomb matrix with 8-fold permutation symmetry for AO integrals
>>> vj = get_jk(mol, dm, 'ijkl,ji->kl', aosym='s8')
>>> # Exchange matrix with 8-fold permutation symmetry for AO integrals
>>> vk = get_jk(mol, dm, 'ijkl,jk->il', aosym='s8')
>>> # Compute coulomb and exchange matrices together
>>> vj, vk = get_jk(mol, (dm,dm), ('ijkl,ji->kl','ijkl,li->kj'), aosym='s8')
>>> # Analytical gradients for coulomb matrix
>>> j1 = get_jk(mol, dm, 'ijkl,lk->ij', intor='int2e_ip1_sph', aosym='s2kl', comp=3)
>>> # contraction across two molecules
>>> mol1 = gto.M(atom='He 2 0 0', basis='6-31g')
>>> nao1 = mol1.nao_nr()
>>> dm1 = numpy.random.random((nao1,nao1))
>>> # Coulomb interaction between two molecules, note 4-fold symmetry can be applied
>>> jcross = get_jk((mol1,mol1,mol,mol), dm, scripts='ijkl,lk->ij', aosym='s4')
>>> ecoul = numpy.einsum('ij,ij', jcross, dm1)
>>> # Exchange interaction between two molecules, no symmetry can be used
>>> kcross = get_jk((mol1,mol,mol,mol1), dm, scripts='ijkl,jk->il')
>>> ex = numpy.einsum('ij,ji', kcross, dm1)
>>> # Analytical gradients for coulomb matrix between two molecules
>>> jcros1 = get_jk((mol1,mol1,mol,mol), dm, scripts='ijkl,lk->ij', intor='int2e_ip1_sph', comp=3)
>>> # Analytical gradients for coulomb interaction between 1s density and the other molecule
>>> jpart1 = get_jk((mol1,mol1,mol,mol), dm, scripts='ijkl,lk->ij', intor='int2e_ip1_sph', comp=3,
...                 shls_slice=(0,1,0,1,0,mol.nbas,0,mol.nbas))
pyscf.scf.jk.jk_build(mols, dms, scripts=['ijkl,ji->kl'], intor='int2e_sph', aosym='s1', comp=None, hermi=0, shls_slice=None, verbose=2, vhfopt=None)#

Compute J/K matrices for the given density matrix

Args:

mols : an instance of Mole or a list of Mole objects

dmsndarray or list of ndarrays

A density matrix or a list of density matrices

Kwargs:
hermiint

Whether the returned J (K) matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
intorstr

2-electron integral name. See getints() for the complete list of available 2-electron integral names

aosymint or str

Permutation symmetry for the AO integrals

4 or ‘4’ or ‘s4’: 4-fold symmetry (default)
‘2ij’ or ‘s2ij’ : symmetry between i, j in (ij|kl)
‘2kl’ or ‘s2kl’ : symmetry between k, l in (ij|kl)
1 or ‘1’ or ‘s1’: no symmetry
‘a4ij’ : 4-fold symmetry with anti-symmetry between i, j in (ij|kl)
‘a4kl’ : 4-fold symmetry with anti-symmetry between k, l in (ij|kl)
‘a2ij’ : anti-symmetry between i, j in (ij|kl)
‘a2kl’ : anti-symmetry between k, l in (ij|kl)
compint

Components of the integrals, e.g. cint2e_ip_sph has 3 components.

scriptsstring or a list of strings

Contraction description (following numpy.einsum convention) based on letters [ijkl]. Each script will be one-to-one applied to each entry of dms. So it must have the same number of elements as the dms, len(scripts) == len(dms).

shls_slice8-element list

(ish_start, ish_end, jsh_start, jsh_end, ksh_start, ksh_end, lsh_start, lsh_end)

Returns:

Depending on the number of density matrices, the function returns one J/K matrix or a list of J/K matrices (the same number of entries as the input dms). Each JK matrices may be a 2D array or 3D array if the AO integral has multiple components.

Examples:

>>> from pyscf import gto
>>> mol = gto.M(atom='H 0 -.5 0; H 0 .5 0', basis='cc-pvdz')
>>> nao = mol.nao_nr()
>>> dm = numpy.random.random((nao,nao))
>>> # Default, Coulomb matrix
>>> vj = get_jk(mol, dm)
>>> # Coulomb matrix with 8-fold permutation symmetry for AO integrals
>>> vj = get_jk(mol, dm, 'ijkl,ji->kl', aosym='s8')
>>> # Exchange matrix with 8-fold permutation symmetry for AO integrals
>>> vk = get_jk(mol, dm, 'ijkl,jk->il', aosym='s8')
>>> # Compute coulomb and exchange matrices together
>>> vj, vk = get_jk(mol, (dm,dm), ('ijkl,ji->kl','ijkl,li->kj'), aosym='s8')
>>> # Analytical gradients for coulomb matrix
>>> j1 = get_jk(mol, dm, 'ijkl,lk->ij', intor='int2e_ip1_sph', aosym='s2kl', comp=3)
>>> # contraction across two molecules
>>> mol1 = gto.M(atom='He 2 0 0', basis='6-31g')
>>> nao1 = mol1.nao_nr()
>>> dm1 = numpy.random.random((nao1,nao1))
>>> # Coulomb interaction between two molecules, note 4-fold symmetry can be applied
>>> jcross = get_jk((mol1,mol1,mol,mol), dm, scripts='ijkl,lk->ij', aosym='s4')
>>> ecoul = numpy.einsum('ij,ij', jcross, dm1)
>>> # Exchange interaction between two molecules, no symmetry can be used
>>> kcross = get_jk((mol1,mol,mol,mol1), dm, scripts='ijkl,jk->il')
>>> ex = numpy.einsum('ij,ji', kcross, dm1)
>>> # Analytical gradients for coulomb matrix between two molecules
>>> jcros1 = get_jk((mol1,mol1,mol,mol), dm, scripts='ijkl,lk->ij', intor='int2e_ip1_sph', comp=3)
>>> # Analytical gradients for coulomb interaction between 1s density and the other molecule
>>> jpart1 = get_jk((mol1,mol1,mol,mol), dm, scripts='ijkl,lk->ij', intor='int2e_ip1_sph', comp=3,
...                 shls_slice=(0,1,0,1,0,mol.nbas,0,mol.nbas))

pyscf.scf.rohf module#

Restricted Open-shell Hartree-Fock

class pyscf.scf.rohf.HF1e(mol)[source]#

Bases: ROHF

scf(*args)#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
class pyscf.scf.rohf.ROHF(mol)[source]#

Bases: RHF

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
CASCI(*args, **kwargs)#
CASSCF(*args, **kwargs)#

CASCI/CASSCF

Args:
mf_or_molSCF object or Mole object

SCF or Mole to define the problem size.

ncasint

Number of active orbitals.

nelecasint or a pair of int

Number of electrons in active space.

Kwargs:
ncoreint

Number of doubly occupied core orbitals. If not presented, this parameter can be automatically determined.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose.

max_memoryfloat or int

Allowed memory in MB. Default value equals to Mole.max_memory.

ncasint

Active space size.

nelecastuple of int

Active (nelec_alpha, nelec_beta)

ncoreint or tuple of int

Core electron number. In UHF-CASSCF, it’s a tuple to indicate the different core electron numbers.

natorbbool

Whether to transform natural orbitals in active space. Note: when CASCI/CASSCF are combined with DMRG solver or selected CI solver, enabling this parameter may slightly change the total energy. False by default.

canonicalizationbool

Whether to canonicalize orbitals in core and external space against the general Fock matrix. The orbitals in active space are NOT transformed by default. To get the natural orbitals in active space, the attribute .natorb needs to be enabled. True by default.

sorting_mo_energybool

Whether to sort the orbitals based on the diagonal elements of the general Fock matrix. Default is False.

fcisolveran instance of FCISolver

The pyscf.fci module provides several FCISolver for different scenario. Generally, fci.direct_spin1.FCISolver can be used for all RHF-CASSCF. However, a proper FCISolver can provide better performance and better numerical stability. One can either use fci.solver() function to pick the FCISolver by the program or manually assigen the FCISolver to this attribute, e.g.

>>> from pyscf import fci
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> mc.fcisolver = fci.solver(mol, singlet=True)
>>> mc.fcisolver = fci.direct_spin1.FCISolver(mol)

You can control FCISolver by setting e.g.:

>>> mc.fcisolver.max_cycle = 30
>>> mc.fcisolver.conv_tol = 1e-7

For more details of the parameter for FCISolver, See fci.

Saved results

e_totfloat

Total MCSCF energy (electronic energy plus nuclear repulsion)

e_casfloat

CAS space FCI energy

cindarray

CAS space FCI coefficients

mo_coeffndarray

When canonicalization is specified, the orbitals are canonical orbitals which make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is specified, the active segment of the mo_coeff is natural orbitals.

mo_energyndarray

Diagonal elements of general Fock matrix (in mo_coeff representation).

mo_occndarray

Occupation numbers of natural orbitals if natorb is specified.

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASCI(mf, 6, 6)
>>> mc.kernel()[0]
-108.980200816243354

Extra attributes for CASSCF:

conv_tolfloat

Converge threshold. Default is 1e-7

conv_tol_gradfloat

Converge threshold for CI gradients and orbital rotation gradients. If not specified, it is set to sqrt(conv_tol).

max_stepsizefloat

The step size for orbital rotation. Small step (0.005 - 0.05) is prefered. Default is 0.02.

max_cycle_macroint

Max number of macro iterations. Default is 50.

max_cycle_microint

Max number of micro iterations in each macro iteration. Depending on systems, increasing this value might reduce the total macro iterations. Generally, 2 - 5 steps should be enough. Default is 4.

small_rot_tolfloat

Threshold for orbital rotation to be considered small. If the largest orbital rotation is smaller than this value, the CI solver will restart from the previous iteration if supported. Default is 0.01

ah_level_shiftfloat, for AH solver.

Level shift for the Davidson diagonalization in AH solver. Default is 1e-8.

ah_conv_tolfloat, for AH solver.

converge threshold for AH solver. Default is 1e-12.

ah_max_cyclefloat, for AH solver.

Max number of iterations allowd in AH solver. Default is 30.

ah_lindepfloat, for AH solver.

Linear dependence threshold for AH solver. Default is 1e-14.

ah_start_tolflat, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 2.5.

ah_start_cycleint, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 3.

ah_conv_tol, ah_max_cycle, ah_lindep, ah_start_tol and ah_start_cycle can affect the accuracy and performance of CASSCF solver. Lower ah_conv_tol and ah_lindep might improve the accuracy of CASSCF optimization, but decrease the performance.

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.conv_tol = 1e-10
>>> mc.ah_conv_tol = 1e-5
>>> mc.kernel()[0]
-109.044401898486001
>>> mc.ah_conv_tol = 1e-10
>>> mc.kernel()[0]
-109.044401887945668
chkfilestr

Checkpoint file to save the intermediate orbitals during the CASSCF optimization. Default is the checkpoint file of mean field object.

ci_response_spaceint

subspace size to solve the CI vector response. Default is 3.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

scale_restorationfloat

When a step of orbital rotation moves out of trust region, the orbital optimization will be restored to previous state and the step size of the orbital rotation needs to be reduced. scale_restoration controls how much to scale down the step size.

Saved results

e_totfloat

Total MCSCF energy (electronic energy plus nuclear repulsion)

e_casfloat

CAS space FCI energy

cindarray

CAS space FCI coefficients

mo_coeffndarray

Optimized CASSCF orbitals coefficients. When canonicalization is specified, the returned orbitals make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is enabled, the active segment of mo_coeff is transformed to natural orbitals.

mo_energyndarray

Diagonal elements of general Fock matrix (in mo_coeff representation).

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.kernel()[0]
-109.044401882238134
CISD = None#
DFMP2 = None#
Gradients(*args, **kwargs)#

Non-relativistic restricted open-shell Hartree-Fock gradients

MP2 = None#
TDA = None#
TDHF = None#
analyze(verbose=None, with_meta_lowdin=True, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis

canonicalize(mo_coeff, mo_occ, fock=None)#

Canonicalization diagonalizes the Fock matrix within occupied, open, virtual subspaces separatedly (without change occupancy).

check_sanity()#

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.

convert_from_(mf)#

Convert the input mean-field object to RHF/ROHF

dump_flags(verbose=None)[source]#
eig(fock, s)[source]#

Solver for generalized eigenvalue problem

\[HC = SCE\]
energy_elec(dm=None, h1e=None, vhf=None)#

Electronic part of Hartree-Fock energy, for given core hamiltonian and HF potential

… math:

E = \sum_{ij}h_{ij} \gamma_{ji}
  + \frac{1}{2}\sum_{ijkl} \gamma_{ji}\gamma_{lk} \langle ik||jl\rangle

Note this function has side effects which cause mf.scf_summary updated.

Args:

mf : an instance of SCF class

Kwargs:
dm2D ndarray

one-particle density matrix

h1e2D ndarray

Core hamiltonian

vhf2D ndarray

HF potential

Returns:

Hartree-Fock electronic energy and the Coulomb energy

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> dm = mf.make_rdm1()
>>> scf.hf.energy_elec(mf, dm)
(-1.5176090667746334, 0.60917167853723675)
>>> mf.energy_elec(dm)
(-1.5176090667746334, 0.60917167853723675)
gen_response(mo_coeff=None, mo_occ=None, with_j=True, hermi=0, max_memory=None)#

Generate a function to compute the product of UHF response function and UHF density matrices.

get_fock(h1e=None, s1e=None, vhf=None, dm=None, cycle=-1, diis=None, diis_start_cycle=None, level_shift_factor=None, damp_factor=None, fock_last=None)#

Build fock matrix based on Roothaan’s effective fock. See also get_roothaan_fock()

get_grad(mo_coeff, mo_occ, fock=None)[source]#

ROHF gradients is the off-diagonal block [co + cv + ov], where [ cc co cv ] [ oc oo ov ] [ vc vo vv ]

get_init_guess(mol=None, key='minao', **kwargs)#
get_occ(mo_energy=None, mo_coeff=None)#

Label the occupancies for each orbital. NOTE the occupancies are not assigned based on the orbital energy ordering. The first N orbitals are assigned to be occupied orbitals.

Examples:

>>> mol = gto.M(atom='H 0 0 0; O 0 0 1.1', spin=1)
>>> mf = scf.hf.SCF(mol)
>>> energy = numpy.array([-10., -1., 1, -2., 0, -3])
>>> mf.get_occ(energy)
array([2, 2, 2, 2, 1, 0])
get_veff(mol=None, dm=None, dm_last=0, vhf_last=0, hermi=1)[source]#

Unrestricted Hartree-Fock potential matrix of alpha and beta spins, for the given density matrix

\[\begin{split}V_{ij}^\alpha &= \sum_{kl} (ij|kl)(\gamma_{lk}^\alpha+\gamma_{lk}^\beta) - \sum_{kl} (il|kj)\gamma_{lk}^\alpha \\ V_{ij}^\beta &= \sum_{kl} (ij|kl)(\gamma_{lk}^\alpha+\gamma_{lk}^\beta) - \sum_{kl} (il|kj)\gamma_{lk}^\beta\end{split}\]
Args:

mol : an instance of Mole

dma list of ndarrays

A list of density matrices, stored as (alpha,alpha,…,beta,beta,…)

Kwargs:
dm_lastndarray or a list of ndarrays or 0

The density matrix baseline. When it is not 0, this function computes the increment of HF potential w.r.t. the reference HF potential matrix.

vhf_lastndarray or a list of ndarrays or 0

The reference HF potential matrix.

hermiint

Whether J, K matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

Returns:

\(V_{hf} = (V^\alpha, V^\beta)\). \(V^\alpha\) (and \(V^\beta\)) can be a list matrices, corresponding to the input density matrices.

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dmsa = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> dmsb = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> dms = numpy.vstack((dmsa,dmsb))
>>> dms.shape
(6, 2, 2)
>>> vhfa, vhfb = scf.uhf.get_veff(mol, dms, hermi=0)
>>> vhfa.shape
(3, 2, 2)
>>> vhfb.shape
(3, 2, 2)
init_guess_by_1e(mol=None)[source]#

Generate initial guess density matrix from core hamiltonian

Returns:

Density matrix, 2D ndarray

init_guess_by_atom(mol=None)[source]#

Generate initial guess density matrix from superposition of atomic HF density matrix. The atomic HF is occupancy averaged RHF

Returns:

Density matrix, 2D ndarray

init_guess_by_chkfile(chkfile=None, project=None)[source]#

Read the HF results from checkpoint file, then project it to the basis defined by mol

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

Returns:

Density matrix, 2D ndarray

init_guess_by_huckel(mol=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

Returns:

Density matrix, 2D ndarray

init_guess_by_minao(mol=None)[source]#

Generate initial guess density matrix based on ANO basis, then project the density matrix to the basis set defined by mol

Returns:

Density matrix, 2D ndarray

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> scf.hf.init_guess_by_minao(mol)
array([[ 0.94758917,  0.09227308],
       [ 0.09227308,  0.94758917]])
init_guess_by_mod_huckel(mol=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

In contrast to init_guess_by_huckel, this routine employs the updated GWH rule from doi:10.1021/ja00480a005 to form the guess.

Returns:

Density matrix, 2D ndarray

init_guess_by_sap(mol=None, **kwargs)[source]#

Generate initial guess density matrix from a superposition of atomic potentials (SAP), doi:10.1021/acs.jctc.8b01089. This is the Gaussian fit implementation, see doi:10.1063/5.0004046.

Args:
molMoleBase object

the molecule object for which the initial guess is evaluated

sap_basisdict

SAP basis in internal format (python dictionary)

Returns:
dm0ndarray

SAP initial guess density matrix

make_rdm1(mo_coeff=None, mo_occ=None, **kwargs)[source]#

One-particle density matrix. mo_occ is a 1D array, with occupancy 1 or 2.

property nelec#
property nelectron_alpha#
nuc_grad_method()[source]#

Hook to create object for analytical nuclear gradients.

spin_square(mo_coeff=None, s=None)[source]#

Spin square and multiplicity of RHF determinant

stability(internal=True, external=False, verbose=None, return_status=False, **kwargs)[source]#

ROHF/ROKS stability analysis.

See also pyscf.scf.stability.rohf_stability function.

Kwargs:
internalbool

Internal stability, within the RHF optimization space.

externalbool

External stability. It is not available in current version.

return_status: bool

Whether to return stable_i and stable_e

Returns:

If return_status is False (default), the return value includes two set of orbitals, which are more close to the stable condition. The first corresponds to the internal stability and the second corresponds to the external stability.

Else, another two boolean variables (indicating current status: stable or unstable) are returned. The first corresponds to the internal stability and the second corresponds to the external stability.

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.

to_ks(xc='HF')[source]#

Convert to ROKS object.

pyscf.scf.rohf.analyze(mf, verbose=5, with_meta_lowdin=True, origin=None, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis

pyscf.scf.rohf.canonicalize(mf, mo_coeff, mo_occ, fock=None)[source]#

Canonicalization diagonalizes the Fock matrix within occupied, open, virtual subspaces separatedly (without change occupancy).

pyscf.scf.rohf.energy_elec(mf, dm=None, h1e=None, vhf=None)[source]#
pyscf.scf.rohf.get_fock(mf, h1e=None, s1e=None, vhf=None, dm=None, cycle=-1, diis=None, diis_start_cycle=None, level_shift_factor=None, damp_factor=None, fock_last=None)[source]#

Build fock matrix based on Roothaan’s effective fock. See also get_roothaan_fock()

pyscf.scf.rohf.get_grad(mo_coeff, mo_occ, fock)[source]#

ROHF gradients is the off-diagonal block [co + cv + ov], where [ cc co cv ] [ oc oo ov ] [ vc vo vv ]

pyscf.scf.rohf.get_occ(mf, mo_energy=None, mo_coeff=None)[source]#

Label the occupancies for each orbital. NOTE the occupancies are not assigned based on the orbital energy ordering. The first N orbitals are assigned to be occupied orbitals.

Examples:

>>> mol = gto.M(atom='H 0 0 0; O 0 0 1.1', spin=1)
>>> mf = scf.hf.SCF(mol)
>>> energy = numpy.array([-10., -1., 1, -2., 0, -3])
>>> mf.get_occ(energy)
array([2, 2, 2, 2, 1, 0])
pyscf.scf.rohf.get_roothaan_fock(focka_fockb, dma_dmb, s)[source]#

Roothaan’s effective fock. Ref. http://www-theor.ch.cam.ac.uk/people/ross/thesis/node15.html

space

closed

open

virtual

closed

Fc

Fb

Fc

open

Fb

Fc

Fa

virtual

Fc

Fa

Fc

where Fc = (Fa + Fb) / 2

Returns:

Roothaan effective Fock matrix

pyscf.scf.rohf.init_guess_by_atom(mol)[source]#
pyscf.scf.rohf.init_guess_by_chkfile(mol, chkfile_name, project=None)[source]#

Read SCF chkfile and make the density matrix for ROHF initial guess.

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

pyscf.scf.rohf.init_guess_by_minao(mol)[source]#
pyscf.scf.rohf.init_guess_by_sap(mol, sap_basis, **kwargs)[source]#
pyscf.scf.rohf.make_rdm1(mo_coeff, mo_occ, **kwargs)[source]#

One-particle density matrix. mo_occ is a 1D array, with occupancy 1 or 2.

pyscf.scf.stability module#

Wave Function Stability Analysis

Ref. JCP 66, 3045 (1977); DOI:10.1063/1.434318 JCP 104, 9047 (1996); DOI:10.1063/1.471637

See also tddft/rhf.py and scf/newton_ah.py

pyscf.scf.stability.dhf_stability(mf, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#

Stability analysis for DHF/DKS method.

Args:

mf : DHF or DKS object

Kwargs:
return_status: bool

Whether to return stable_i and stable_e

nrootsint

Number of roots solved by Davidson solver

tolfloat

Convergence threshold for Davidson solver

Returns:

If return_status is False (default), the return value includes a new set of orbitals, which are more close to the stable condition.

Else, another one boolean variable (indicating current status: stable or unstable) is returned.

pyscf.scf.stability.dump_status(log, stable, method_class, stab_type)[source]#
pyscf.scf.stability.ghf_stability(mf, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#

Stability analysis for GHF/GKS method.

Args:

mf : GHF or GKS object

Kwargs:
return_status: bool

Whether to return stable_i and stable_e

nrootsint

Number of roots solved by Davidson solver

tolfloat

Convergence threshold for Davidson solver

Returns:

If return_status is False (default), the return value includes a new set of orbitals, which are more close to the stable condition.

Else, another one boolean variable (indicating current status: stable or unstable) is returned.

pyscf.scf.stability.rhf_external(mf, with_symmetry=True, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#
pyscf.scf.stability.rhf_internal(mf, with_symmetry=True, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#
pyscf.scf.stability.rhf_stability(mf, internal=True, external=False, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#

Stability analysis for RHF/RKS method.

Args:

mf : RHF or RKS object

Kwargs:
internalbool

Internal stability, within the RHF space.

externalbool

External stability. Including the RHF -> UHF and real -> complex stability analysis.

return_status: bool

Whether to return stable_i and stable_e

nrootsint

Number of roots solved by Davidson solver

tolfloat

Convergence threshold for Davidson solver

Returns:

If return_status is False (default), the return value includes two set of orbitals, which are more close to the stable condition. The first corresponds to the internal stability and the second corresponds to the external stability.

Else, another two boolean variables (indicating current status: stable or unstable) are returned. The first corresponds to the internal stability and the second corresponds to the external stability.

pyscf.scf.stability.rohf_external(mf, with_symmetry=True, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#
pyscf.scf.stability.rohf_internal(mf, with_symmetry=True, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#
pyscf.scf.stability.rohf_stability(mf, internal=True, external=False, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#

Stability analysis for ROHF/ROKS method.

Args:

mf : ROHF or ROKS object

Kwargs:
internalbool

Internal stability, within the RHF space.

externalbool

External stability. It is not available in current version.

return_status: bool

Whether to return stable_i and stable_e

nrootsint

Number of roots solved by Davidson solver

tolfloat

Convergence threshold for Davidson solver

Returns:

If return_status is False (default), the return value includes two set of orbitals, which are more close to the stable condition. The first corresponds to the internal stability and the second corresponds to the external stability.

Else, another two boolean variables (indicating current status: stable or unstable) are returned. The first corresponds to the internal stability and the second corresponds to the external stability.

pyscf.scf.stability.uhf_external(mf, with_symmetry=True, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#
pyscf.scf.stability.uhf_internal(mf, with_symmetry=True, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#
pyscf.scf.stability.uhf_stability(mf, internal=True, external=False, verbose=None, return_status=False, nroots=3, tol=0.0001)[source]#

Stability analysis for UHF/UKS method.

Args:

mf : UHF or UKS object

Kwargs:
internalbool

Internal stability, within the UHF space.

externalbool

External stability. Including the UHF -> GHF and real -> complex stability analysis.

return_status: bool

Whether to return stable_i and stable_e

nrootsint

Number of roots solved by Davidson solver

tolfloat

Convergence threshold for Davidson solver

Returns:

If return_status is False (default), the return value includes two set of orbitals, which are more close to the stable condition. The first corresponds to the internal stability and the second corresponds to the external stability.

Else, another two boolean variables (indicating current status: stable or unstable) are returned. The first corresponds to the internal stability and the second corresponds to the external stability.

pyscf.scf.stability_slow module#

Wave Function Stability Analysis

Ref. JCP, 66, 3045 (1977); DOI:10.1063/1.434318 JCP 104, 9047 (1996); DOI:10.1063/1.471637

pyscf.scf.stability_slow.rhf_external(mf, verbose=None)[source]#
pyscf.scf.stability_slow.rhf_internal(mf, verbose=None)[source]#
pyscf.scf.stability_slow.rhf_stability(mf, internal=True, external=False, verbose=None)[source]#
pyscf.scf.stability_slow.uhf_external(mf, verbose=None)[source]#
pyscf.scf.stability_slow.uhf_internal(mf, verbose=None)[source]#
pyscf.scf.stability_slow.uhf_stability(mf, internal=True, external=False, verbose=None)[source]#

pyscf.scf.ucphf module#

Unrestricted coupled perturbed Hartree-Fock solver

pyscf.scf.ucphf.kernel(fvind, mo_energy, mo_occ, h1, s1=None, max_cycle=50, tol=1e-09, hermi=False, verbose=2, level_shift=0)#
Args:
fvindfunction

Given density matrix, compute (ij|kl)D_{lk}*2 - (ij|kl)D_{jk}

pyscf.scf.ucphf.solve(fvind, mo_energy, mo_occ, h1, s1=None, max_cycle=50, tol=1e-09, hermi=False, verbose=2, level_shift=0)[source]#
Args:
fvindfunction

Given density matrix, compute (ij|kl)D_{lk}*2 - (ij|kl)D_{jk}

pyscf.scf.ucphf.solve_nos1(fvind, mo_energy, mo_occ, h1, max_cycle=20, tol=1e-09, hermi=False, verbose=2, level_shift=0)[source]#

For field independent basis. First order overlap matrix is zero

pyscf.scf.ucphf.solve_withs1(fvind, mo_energy, mo_occ, h1, s1, max_cycle=50, tol=1e-09, hermi=False, verbose=2, level_shift=0)[source]#

For field dependent basis. First order overlap matrix is non-zero. The first order orbitals are set to C^1_{ij} = -1/2 S1 e1 = h1 - s1*e0 + (e0_j-e0_i)*c1 + vhf[c1]

pyscf.scf.uhf module#

class pyscf.scf.uhf.HF1e(mol)[source]#

Bases: UHF

scf(*args)#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
spin_square(mo_coeff=None, s=None)[source]#

Spin square and multiplicity of UHF determinant

\[S^2 = \frac{1}{2}(S_+ S_- + S_- S_+) + S_z^2\]

where \(S_+ = \sum_i S_{i+}\) is effective for all beta occupied orbitals; \(S_- = \sum_i S_{i-}\) is effective for all alpha occupied orbitals.

  1. There are two possibilities for \(S_+ S_-\)
    1. same electron \(S_+ S_- = \sum_i s_{i+} s_{i-}\),

    \[\sum_i \langle UHF|s_{i+} s_{i-}|UHF\rangle = \sum_{pq}\langle p|s_+s_-|q\rangle \gamma_{qp} = n_\alpha\]

    2) different electrons \(S_+ S_- = \sum s_{i+} s_{j-}, (i\neq j)\). There are in total \(n(n-1)\) terms. As a two-particle operator,

    \[\langle S_+ S_- \rangle = \langle ij|s_+ s_-|ij\rangle - \langle ij|s_+ s_-|ji\rangle = -\langle i^\alpha|j^\beta\rangle \langle j^\beta|i^\alpha\rangle\]
  2. Similarly, for \(S_- S_+\)
    1. same electron

    \[\sum_i \langle s_{i-} s_{i+}\rangle = n_\beta\]
    1. different electrons

    \[\langle S_- S_+ \rangle = -\langle i^\beta|j^\alpha\rangle \langle j^\alpha|i^\beta\rangle\]
  3. For \(S_z^2\)
    1. same electron

    \[\langle s_z^2\rangle = \frac{1}{4}(n_\alpha + n_\beta)\]
    1. different electrons

    \[\begin{split}&\frac{1}{2}\sum_{ij}(\langle ij|2s_{z1}s_{z2}|ij\rangle -\langle ij|2s_{z1}s_{z2}|ji\rangle) \\ &=\frac{1}{4}(\langle i^\alpha|i^\alpha\rangle \langle j^\alpha|j^\alpha\rangle - \langle i^\alpha|i^\alpha\rangle \langle j^\beta|j^\beta\rangle - \langle i^\beta|i^\beta\rangle \langle j^\alpha|j^\alpha\rangle + \langle i^\beta|i^\beta\rangle \langle j^\beta|j^\beta\rangle) \\ &-\frac{1}{4}(\langle i^\alpha|j^\alpha\rangle \langle j^\alpha|i^\alpha\rangle + \langle i^\beta|j^\beta\rangle\langle j^\beta|i^\beta\rangle) \\ &=\frac{1}{4}(n_\alpha^2 - n_\alpha n_\beta - n_\beta n_\alpha + n_\beta^2) -\frac{1}{4}(n_\alpha + n_\beta) \\ &=\frac{1}{4}((n_\alpha-n_\beta)^2 - (n_\alpha+n_\beta))\end{split}\]

In total

\[\begin{split}\langle S^2\rangle &= \frac{1}{2} (n_\alpha-\sum_{ij}\langle i^\alpha|j^\beta\rangle \langle j^\beta|i^\alpha\rangle +n_\beta -\sum_{ij}\langle i^\beta|j^\alpha\rangle\langle j^\alpha|i^\beta\rangle) + \frac{1}{4}(n_\alpha-n_\beta)^2 \\\end{split}\]
Args:
moa list of 2 ndarrays

Occupied alpha and occupied beta orbitals

Kwargs:
sndarray

AO overlap

Returns:

A list of two floats. The first is the expectation value of S^2. The second is the corresponding 2S+1

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', charge=1, spin=1, verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.kernel()
-75.623975516256706
>>> mo = (mf.mo_coeff[0][:,mf.mo_occ[0]>0], mf.mo_coeff[1][:,mf.mo_occ[1]>0])
>>> print('S^2 = %.7f, 2S+1 = %.7f' % spin_square(mo, mol.intor('int1e_ovlp_sph')))
S^2 = 0.7570150, 2S+1 = 2.0070027
class pyscf.scf.uhf.UHF(mol)[source]#

Bases: SCF

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
Attributes for UHF:
nelec(int, int)

If given, freeze the number of (alpha,beta) electrons to the given value.

level_shiftnumber or two-element list

level shift (in Eh) for alpha and beta Fock if two-element list is given.

init_guess_breaksymint

This configuration controls the algorithm used to break the spin symmetry of the initial guess: - 0 to disable symmetry breaking in the initial guess. - 1 to use the default algorithm introduced in pyscf-1.7. - 2 to adjust the num. electrons for spin-up and spin-down density matrices (issue #1839).

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', charge=1, spin=1, verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.kernel()
-75.623975516256706
>>> print('S^2 = %.7f, 2S+1 = %.7f' % mf.spin_square())
S^2 = 0.7570150, 2S+1 = 2.0070027
CASCI = None#
CASSCF = None#
CISD(*args, **kwargs)#
DFMP2 = None#
Gradients(*args, **kwargs)#

Non-relativistic unrestricted Hartree-Fock gradients

MP2(*args, **kwargs)#
TDA(*args, **kwargs)#

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.

TDHF(*args, **kwargs)#
analyze(verbose=None, with_meta_lowdin=True, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis; Diople moment.

canonicalize(mo_coeff, mo_occ, fock=None)#

Canonicalization diagonalizes the UHF Fock matrix within occupied, virtual subspaces separatedly (without change occupancy).

convert_from_(mf)[source]#

Create UHF object based on the RHF/ROHF object

det_ovlp(mo1, mo2, occ1, occ2, ovlp=None)[source]#

Calculate the overlap between two different determinants. It is the product of single values of molecular orbital overlap matrix.

\[S_{12} = \langle \Psi_A | \Psi_B \rangle = (\mathrm{det}\mathbf{U}) (\mathrm{det}\mathbf{V^\dagger}) \prod\limits_{i=1}\limits^{2N} \lambda_{ii}\]

where \(\mathbf{U}, \mathbf{V}, \lambda\) are unitary matrices and single values generated by single value decomposition(SVD) of the overlap matrix \(\mathbf{O}\) which is the overlap matrix of two sets of molecular orbitals:

\[\mathbf{U}^\dagger \mathbf{O} \mathbf{V} = \mathbf{\Lambda}\]
Args:
mo1, mo22D ndarrays

Molecualr orbital coefficients

occ1, occ2: 2D ndarrays

occupation numbers

Return:
A list:

the product of single values: float (x_a, x_b): 1D ndarrays \(\mathbf{U} \mathbf{\Lambda}^{-1} \mathbf{V}^\dagger\) They are used to calculate asymmetric density matrix

dump_flags(verbose=None)[source]#
eig(fock, s)[source]#

Solver for generalized eigenvalue problem

\[HC = SCE\]
energy_elec(dm=None, h1e=None, vhf=None)#

Electronic energy of Unrestricted Hartree-Fock

Note this function has side effects which cause mf.scf_summary updated.

Returns:

Hartree-Fock electronic energy and the 2-electron part contribution

gen_response(mo_coeff=None, mo_occ=None, with_j=True, hermi=0, max_memory=None)#

Generate a function to compute the product of UHF response function and UHF density matrices.

get_fock(h1e=None, s1e=None, vhf=None, dm=None, cycle=-1, diis=None, diis_start_cycle=None, level_shift_factor=None, damp_factor=None, fock_last=None)#

F = h^{core} + V^{HF}

Special treatment (damping, DIIS, or level shift) will be applied to the Fock matrix if diis and cycle is specified (The two parameters are passed to get_fock function during the SCF iteration)

Kwargs:
h1e2D ndarray

Core hamiltonian

s1e2D ndarray

Overlap matrix, for DIIS

vhf2D ndarray

HF potential matrix

dm2D ndarray

Density matrix, for DIIS

cycleint

Then present SCF iteration step, for DIIS

diisan object of SCF.DIIS class

DIIS object to hold intermediate Fock and error vectors

diis_start_cycleint

The step to start DIIS. Default is 0.

level_shift_factorfloat or int

Level shift (in AU) for virtual space. Default is 0.

get_grad(mo_coeff, mo_occ, fock=None)[source]#

RHF orbital gradients

Args:
mo_coeff2D ndarray

Obital coefficients

mo_occ1D ndarray

Orbital occupancy

fock_ao2D ndarray

Fock matrix in AO representation

Returns:

Gradients in MO representation. It’s a num_occ*num_vir vector.

get_init_guess(mol=None, key='minao', **kwargs)[source]#
get_jk(mol=None, dm=None, hermi=1, with_j=True, with_k=True, omega=None)[source]#

Coulomb (J) and exchange (K)

Args:
dma list of 2D arrays or a list of 3D arrays

(alpha_dm, beta_dm) or (alpha_dms, beta_dms)

get_occ(mo_energy=None, mo_coeff=None)#

Label the occupancies for each orbital

Kwargs:
mo_energy1D ndarray

Obital energies

mo_coeff2D ndarray

Obital coefficients

Examples:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> energy = numpy.array([-10., -1., 1, -2., 0, -3])
>>> mf.get_occ(energy)
array([2, 2, 0, 2, 2, 2])
get_veff(mol=None, dm=None, dm_last=0, vhf_last=0, hermi=1)[source]#

Unrestricted Hartree-Fock potential matrix of alpha and beta spins, for the given density matrix

\[\begin{split}V_{ij}^\alpha &= \sum_{kl} (ij|kl)(\gamma_{lk}^\alpha+\gamma_{lk}^\beta) - \sum_{kl} (il|kj)\gamma_{lk}^\alpha \\ V_{ij}^\beta &= \sum_{kl} (ij|kl)(\gamma_{lk}^\alpha+\gamma_{lk}^\beta) - \sum_{kl} (il|kj)\gamma_{lk}^\beta\end{split}\]
Args:

mol : an instance of Mole

dma list of ndarrays

A list of density matrices, stored as (alpha,alpha,…,beta,beta,…)

Kwargs:
dm_lastndarray or a list of ndarrays or 0

The density matrix baseline. When it is not 0, this function computes the increment of HF potential w.r.t. the reference HF potential matrix.

vhf_lastndarray or a list of ndarrays or 0

The reference HF potential matrix.

hermiint

Whether J, K matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

Returns:

\(V_{hf} = (V^\alpha, V^\beta)\). \(V^\alpha\) (and \(V^\beta\)) can be a list matrices, corresponding to the input density matrices.

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dmsa = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> dmsb = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> dms = numpy.vstack((dmsa,dmsb))
>>> dms.shape
(6, 2, 2)
>>> vhfa, vhfb = scf.uhf.get_veff(mol, dms, hermi=0)
>>> vhfa.shape
(3, 2, 2)
>>> vhfb.shape
(3, 2, 2)
init_guess_breaksym = 1#
init_guess_by_1e(mol=None, breaksym=None)[source]#

Generate initial guess density matrix from core hamiltonian

Returns:

Density matrix, 2D ndarray

init_guess_by_atom(mol=None, breaksym=None)[source]#

Generate initial guess density matrix from superposition of atomic HF density matrix. The atomic HF is occupancy averaged RHF

Returns:

Density matrix, 2D ndarray

init_guess_by_chkfile(chkfile=None, project=None)[source]#

Read the HF results from checkpoint file, then project it to the basis defined by mol

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

Returns:

Density matrix, 2D ndarray

init_guess_by_huckel(mol=None, breaksym=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

Returns:

Density matrix, 2D ndarray

init_guess_by_minao(mol=None, breaksym=None)[source]#

Initial guess in terms of the overlap to minimal basis.

init_guess_by_mod_huckel(mol=None, breaksym=None)[source]#

Generate initial guess density matrix from a Huckel calculation based on occupancy averaged atomic RHF calculations, doi:10.1021/acs.jctc.8b01089

In contrast to init_guess_by_huckel, this routine employs the updated GWH rule from doi:10.1021/ja00480a005 to form the guess.

Returns:

Density matrix, 2D ndarray

init_guess_by_sap(mol=None, breaksym=None, **kwargs)[source]#

Generate initial guess density matrix from a superposition of atomic potentials (SAP), doi:10.1021/acs.jctc.8b01089. This is the Gaussian fit implementation, see doi:10.1063/5.0004046.

Args:
molMoleBase object

the molecule object for which the initial guess is evaluated

sap_basisdict

SAP basis in internal format (python dictionary)

Returns:
dm0ndarray

SAP initial guess density matrix

make_asym_dm(mo1, mo2, occ1, occ2, x)[source]#

One-particle asymmetric density matrix

Args:
mo1, mo22D ndarrays

Molecualr orbital coefficients

occ1, occ2: 2D ndarrays

Occupation numbers

x: 2D ndarrays

\(\mathbf{U} \mathbf{\Lambda}^{-1} \mathbf{V}^\dagger\). See also det_ovlp()

Return:

A list of 2D ndarrays for alpha and beta spin

Examples:

>>> mf1 = scf.UHF(gto.M(atom='H 0 0 0; F 0 0 1.3', basis='ccpvdz')).run()
>>> mf2 = scf.UHF(gto.M(atom='H 0 0 0; F 0 0 1.4', basis='ccpvdz')).run()
>>> s = gto.intor_cross('int1e_ovlp_sph', mf1.mol, mf2.mol)
>>> det, x = det_ovlp(mf1.mo_coeff, mf1.mo_occ, mf2.mo_coeff, mf2.mo_occ, s)
>>> adm = make_asym_dm(mf1.mo_coeff, mf1.mo_occ, mf2.mo_coeff, mf2.mo_occ, x)
>>> adm.shape
(2, 19, 19)
make_rdm1(mo_coeff=None, mo_occ=None, **kwargs)[source]#

One-particle density matrix in AO representation

Args:
mo_coefftuple of 2D ndarrays

Orbital coefficients for alpha and beta spins. Each column is one orbital.

mo_occtuple of 1D ndarrays

Occupancies for alpha and beta spins.

Returns:

A list of 2D ndarrays for alpha and beta spins

make_rdm2(mo_coeff=None, mo_occ=None, **kwargs)[source]#

Two-particle density matrix in AO representation

Args:
mo_coefftuple of 2D ndarrays

Orbital coefficients for alpha and beta spins. Each column is one orbital.

mo_occtuple of 1D ndarrays

Occupancies for alpha and beta spins.

Returns:

A tuple of three 4D ndarrays for alpha,alpha and alpha,beta and beta,beta spins

mulliken_meta(mol=None, dm=None, verbose=5, pre_orth_method='ANO', s=None)[source]#

Mulliken population analysis, based on meta-Lowdin AOs. In the meta-lowdin, the AOs are grouped in three sets: core, valence and Rydberg, the orthogonalization are carried out within each subsets.

Args:

mol : an instance of Mole

dmndarray or 2-item list of ndarray

Density matrix. ROHF dm is a 2-item list of 2D array

Kwargs:

verbose : int or instance of lib.logger.Logger

pre_orth_methodstr

Pre-orthogonalization, which localized GTOs for each atom. To obtain the occupied and unoccupied atomic shells, there are three methods

‘ano’ : Project GTOs to ANO basis
‘minao’ : Project GTOs to MINAO basis
‘scf’ : Symmetry-averaged fractional occupation atomic RHF
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

mulliken_meta_spin(mol=None, dm=None, verbose=5, pre_orth_method='ANO', s=None)[source]#
mulliken_pop(mol=None, dm=None, s=None, verbose=5)[source]#

Mulliken population analysis

\[M_{ij} = D_{ij} S_{ji}\]

Mulliken charges

\[\delta_i = \sum_j M_{ij}\]
Returns:

A list : pop, charges

popnparray

Mulliken population on each atomic orbitals

chargesnparray

Mulliken charges

mulliken_spin_pop(mol=None, dm=None, s=None, verbose=5)[source]#
property nelec#
property nelectron_alpha#
nuc_grad_method()[source]#

Hook to create object for analytical nuclear gradients.

spin_square(mo_coeff=None, s=None)[source]#

Spin square and multiplicity of UHF determinant

\[S^2 = \frac{1}{2}(S_+ S_- + S_- S_+) + S_z^2\]

where \(S_+ = \sum_i S_{i+}\) is effective for all beta occupied orbitals; \(S_- = \sum_i S_{i-}\) is effective for all alpha occupied orbitals.

  1. There are two possibilities for \(S_+ S_-\)
    1. same electron \(S_+ S_- = \sum_i s_{i+} s_{i-}\),

    \[\sum_i \langle UHF|s_{i+} s_{i-}|UHF\rangle = \sum_{pq}\langle p|s_+s_-|q\rangle \gamma_{qp} = n_\alpha\]

    2) different electrons \(S_+ S_- = \sum s_{i+} s_{j-}, (i\neq j)\). There are in total \(n(n-1)\) terms. As a two-particle operator,

    \[\langle S_+ S_- \rangle = \langle ij|s_+ s_-|ij\rangle - \langle ij|s_+ s_-|ji\rangle = -\langle i^\alpha|j^\beta\rangle \langle j^\beta|i^\alpha\rangle\]
  2. Similarly, for \(S_- S_+\)
    1. same electron

    \[\sum_i \langle s_{i-} s_{i+}\rangle = n_\beta\]
    1. different electrons

    \[\langle S_- S_+ \rangle = -\langle i^\beta|j^\alpha\rangle \langle j^\alpha|i^\beta\rangle\]
  3. For \(S_z^2\)
    1. same electron

    \[\langle s_z^2\rangle = \frac{1}{4}(n_\alpha + n_\beta)\]
    1. different electrons

    \[\begin{split}&\frac{1}{2}\sum_{ij}(\langle ij|2s_{z1}s_{z2}|ij\rangle -\langle ij|2s_{z1}s_{z2}|ji\rangle) \\ &=\frac{1}{4}(\langle i^\alpha|i^\alpha\rangle \langle j^\alpha|j^\alpha\rangle - \langle i^\alpha|i^\alpha\rangle \langle j^\beta|j^\beta\rangle - \langle i^\beta|i^\beta\rangle \langle j^\alpha|j^\alpha\rangle + \langle i^\beta|i^\beta\rangle \langle j^\beta|j^\beta\rangle) \\ &-\frac{1}{4}(\langle i^\alpha|j^\alpha\rangle \langle j^\alpha|i^\alpha\rangle + \langle i^\beta|j^\beta\rangle\langle j^\beta|i^\beta\rangle) \\ &=\frac{1}{4}(n_\alpha^2 - n_\alpha n_\beta - n_\beta n_\alpha + n_\beta^2) -\frac{1}{4}(n_\alpha + n_\beta) \\ &=\frac{1}{4}((n_\alpha-n_\beta)^2 - (n_\alpha+n_\beta))\end{split}\]

In total

\[\begin{split}\langle S^2\rangle &= \frac{1}{2} (n_\alpha-\sum_{ij}\langle i^\alpha|j^\beta\rangle \langle j^\beta|i^\alpha\rangle +n_\beta -\sum_{ij}\langle i^\beta|j^\alpha\rangle\langle j^\alpha|i^\beta\rangle) + \frac{1}{4}(n_\alpha-n_\beta)^2 \\\end{split}\]
Args:
moa list of 2 ndarrays

Occupied alpha and occupied beta orbitals

Kwargs:
sndarray

AO overlap

Returns:

A list of two floats. The first is the expectation value of S^2. The second is the corresponding 2S+1

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', charge=1, spin=1, verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.kernel()
-75.623975516256706
>>> mo = (mf.mo_coeff[0][:,mf.mo_occ[0]>0], mf.mo_coeff[1][:,mf.mo_occ[1]>0])
>>> print('S^2 = %.7f, 2S+1 = %.7f' % spin_square(mo, mol.intor('int1e_ovlp_sph')))
S^2 = 0.7570150, 2S+1 = 2.0070027
stability(internal=True, external=False, verbose=None, return_status=False, **kwargs)[source]#

Stability analysis for UHF/UKS method.

See also pyscf.scf.stability.uhf_stability function.

Args:

mf : UHF or UKS object

Kwargs:
internalbool

Internal stability, within the UHF space.

externalbool

External stability. Including the UHF -> GHF and real -> complex stability analysis.

return_status: bool

Whether to return stable_i and stable_e

Returns:

If return_status is False (default), the return value includes two set of orbitals, which are more close to the stable condition. The first corresponds to the internal stability and the second corresponds to the external stability.

Else, another two boolean variables (indicating current status: stable or unstable) are returned. The first corresponds to the internal stability and the second corresponds to the external stability.

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.

to_ks(xc='HF')[source]#

Convert to UKS object.

pyscf.scf.uhf.analyze(mf, verbose=5, with_meta_lowdin=True, **kwargs)[source]#

Analyze the given SCF object: print orbital energies, occupancies; print orbital coefficients; Mulliken population analysis; Dipole moment; Spin density for AOs and atoms;

pyscf.scf.uhf.canonicalize(mf, mo_coeff, mo_occ, fock=None)[source]#

Canonicalization diagonalizes the UHF Fock matrix within occupied, virtual subspaces separatedly (without change occupancy).

pyscf.scf.uhf.det_ovlp(mo1, mo2, occ1, occ2, ovlp)[source]#

Calculate the overlap between two different determinants. It is the product of single values of molecular orbital overlap matrix.

\[S_{12} = \langle \Psi_A | \Psi_B \rangle = (\mathrm{det}\mathbf{U}) (\mathrm{det}\mathbf{V^\dagger}) \prod\limits_{i=1}\limits^{2N} \lambda_{ii}\]

where \(\mathbf{U}, \mathbf{V}, \lambda\) are unitary matrices and single values generated by single value decomposition(SVD) of the overlap matrix \(\mathbf{O}\) which is the overlap matrix of two sets of molecular orbitals:

\[\mathbf{U}^\dagger \mathbf{O} \mathbf{V} = \mathbf{\Lambda}\]
Args:
mo1, mo22D ndarrays

Molecualr orbital coefficients

occ1, occ2: 2D ndarrays

occupation numbers

Return:
A list:

the product of single values: float (x_a, x_b): 1D ndarrays \(\mathbf{U} \mathbf{\Lambda}^{-1} \mathbf{V}^\dagger\) They are used to calculate asymmetric density matrix

pyscf.scf.uhf.energy_elec(mf, dm=None, h1e=None, vhf=None)[source]#

Electronic energy of Unrestricted Hartree-Fock

Note this function has side effects which cause mf.scf_summary updated.

Returns:

Hartree-Fock electronic energy and the 2-electron part contribution

pyscf.scf.uhf.get_fock(mf, h1e=None, s1e=None, vhf=None, dm=None, cycle=-1, diis=None, diis_start_cycle=None, level_shift_factor=None, damp_factor=None, fock_last=None)[source]#
pyscf.scf.uhf.get_grad(mo_coeff, mo_occ, fock_ao)[source]#

UHF Gradients

pyscf.scf.uhf.get_init_guess(mol, key='minao', **kwargs)[source]#
pyscf.scf.uhf.get_occ(mf, mo_energy=None, mo_coeff=None)[source]#
pyscf.scf.uhf.get_veff(mol, dm, dm_last=0, vhf_last=0, hermi=1, vhfopt=None)[source]#

Unrestricted Hartree-Fock potential matrix of alpha and beta spins, for the given density matrix

\[\begin{split}V_{ij}^\alpha &= \sum_{kl} (ij|kl)(\gamma_{lk}^\alpha+\gamma_{lk}^\beta) - \sum_{kl} (il|kj)\gamma_{lk}^\alpha \\ V_{ij}^\beta &= \sum_{kl} (ij|kl)(\gamma_{lk}^\alpha+\gamma_{lk}^\beta) - \sum_{kl} (il|kj)\gamma_{lk}^\beta\end{split}\]
Args:

mol : an instance of Mole

dma list of ndarrays

A list of density matrices, stored as (alpha,alpha,…,beta,beta,…)

Kwargs:
dm_lastndarray or a list of ndarrays or 0

The density matrix baseline. When it is not 0, this function computes the increment of HF potential w.r.t. the reference HF potential matrix.

vhf_lastndarray or a list of ndarrays or 0

The reference HF potential matrix.

hermiint

Whether J, K matrix is hermitian

0 : no hermitian or symmetric
1 : hermitian
2 : anti-hermitian
vhfopt :

A class which holds precomputed quantities to optimize the computation of J, K matrices

Returns:

\(V_{hf} = (V^\alpha, V^\beta)\). \(V^\alpha\) (and \(V^\beta\)) can be a list matrices, corresponding to the input density matrices.

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1')
>>> dmsa = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> dmsb = numpy.random.random((3,mol.nao_nr(),mol.nao_nr()))
>>> dms = numpy.vstack((dmsa,dmsb))
>>> dms.shape
(6, 2, 2)
>>> vhfa, vhfb = scf.uhf.get_veff(mol, dms, hermi=0)
>>> vhfa.shape
(3, 2, 2)
>>> vhfb.shape
(3, 2, 2)
pyscf.scf.uhf.init_guess_by_1e(mol, breaksym=None)[source]#
pyscf.scf.uhf.init_guess_by_atom(mol, breaksym=None)[source]#
pyscf.scf.uhf.init_guess_by_chkfile(mol, chkfile_name, project=None)[source]#

Read SCF chkfile and make the density matrix for UHF initial guess.

Kwargs:
projectNone or bool

Whether to project chkfile’s orbitals to the new basis. Note when the geometry of the chkfile and the given molecule are very different, this projection can produce very poor initial guess. In PES scanning, it is recommended to switch off project.

If project is set to None, the projection is only applied when the basis sets of the chkfile’s molecule are different to the basis sets of the given molecule (regardless whether the geometry of the two molecules are different). Note the basis sets are considered to be different if the two molecules are derived from the same molecule with different ordering of atoms.

pyscf.scf.uhf.init_guess_by_huckel(mol, breaksym=None)[source]#
pyscf.scf.uhf.init_guess_by_minao(mol, breaksym=None)[source]#

Generate initial guess density matrix based on ANO basis, then project the density matrix to the basis set defined by mol

Returns:

Density matrices, a list of 2D ndarrays for alpha and beta spins

pyscf.scf.uhf.init_guess_by_mod_huckel(mol, breaksym=None)[source]#
pyscf.scf.uhf.init_guess_by_sap(mol, sap_basis, breaksym=None, **kwargs)[source]#
pyscf.scf.uhf.make_asym_dm(mo1, mo2, occ1, occ2, x)[source]#

One-particle asymmetric density matrix

Args:
mo1, mo22D ndarrays

Molecualr orbital coefficients

occ1, occ2: 2D ndarrays

Occupation numbers

x: 2D ndarrays

\(\mathbf{U} \mathbf{\Lambda}^{-1} \mathbf{V}^\dagger\). See also det_ovlp()

Return:

A list of 2D ndarrays for alpha and beta spin

Examples:

>>> mf1 = scf.UHF(gto.M(atom='H 0 0 0; F 0 0 1.3', basis='ccpvdz')).run()
>>> mf2 = scf.UHF(gto.M(atom='H 0 0 0; F 0 0 1.4', basis='ccpvdz')).run()
>>> s = gto.intor_cross('int1e_ovlp_sph', mf1.mol, mf2.mol)
>>> det, x = det_ovlp(mf1.mo_coeff, mf1.mo_occ, mf2.mo_coeff, mf2.mo_occ, s)
>>> adm = make_asym_dm(mf1.mo_coeff, mf1.mo_occ, mf2.mo_coeff, mf2.mo_occ, x)
>>> adm.shape
(2, 19, 19)
pyscf.scf.uhf.make_rdm1(mo_coeff, mo_occ, **kwargs)[source]#

One-particle density matrix in AO representation

Args:
mo_coefftuple of 2D ndarrays

Orbital coefficients for alpha and beta spins. Each column is one orbital.

mo_occtuple of 1D ndarrays

Occupancies for alpha and beta spins.

Returns:

A list of 2D ndarrays for alpha and beta spins

pyscf.scf.uhf.make_rdm2(mo_coeff, mo_occ)[source]#

Two-particle density matrix in AO representation

Args:
mo_coefftuple of 2D ndarrays

Orbital coefficients for alpha and beta spins. Each column is one orbital.

mo_occtuple of 1D ndarrays

Occupancies for alpha and beta spins.

Returns:

A tuple of three 4D ndarrays for alpha,alpha and alpha,beta and beta,beta spins

pyscf.scf.uhf.mulliken_meta(mol, dm_ao, verbose=5, pre_orth_method='ANO', s=None)[source]#

Mulliken population analysis, based on meta-Lowdin AOs.

pyscf.scf.uhf.mulliken_meta_spin(mol, dm_ao, verbose=5, pre_orth_method='ANO', s=None)[source]#

Mulliken spin population analysis, based on meta-Lowdin AOs.

pyscf.scf.uhf.mulliken_pop(mol, dm, s=None, verbose=5)[source]#

Mulliken population analysis

pyscf.scf.uhf.mulliken_pop_meta_lowdin_ao(mol, dm_ao, verbose=5, pre_orth_method='ANO', s=None)#

Mulliken population analysis, based on meta-Lowdin AOs.

pyscf.scf.uhf.mulliken_spin_pop(mol, dm, s=None, verbose=5)[source]#

Mulliken spin density analysis

See Eq. 80 in https://arxiv.org/pdf/1206.2234.pdf and the surrounding text for more details.

\[M_{ij} = (D^a_{ij} - D^b_{ij}) S_{ji}\]

Mulliken charges

\[\delta_i = \sum_j M_{ij}\]
Returns:

A list : spin_pop, Ms

spin_popnparray

Mulliken spin density on each atomic orbitals

Msnparray

Mulliken spin density on each atom

pyscf.scf.uhf.mulliken_spin_pop_meta_lowdin_ao(mol, dm_ao, verbose=5, pre_orth_method='ANO', s=None)#

Mulliken spin population analysis, based on meta-Lowdin AOs.

pyscf.scf.uhf.spin_square(mo, s=1)[source]#

Spin square and multiplicity of UHF determinant

\[S^2 = \frac{1}{2}(S_+ S_- + S_- S_+) + S_z^2\]

where \(S_+ = \sum_i S_{i+}\) is effective for all beta occupied orbitals; \(S_- = \sum_i S_{i-}\) is effective for all alpha occupied orbitals.

  1. There are two possibilities for \(S_+ S_-\)
    1. same electron \(S_+ S_- = \sum_i s_{i+} s_{i-}\),

    \[\sum_i \langle UHF|s_{i+} s_{i-}|UHF\rangle = \sum_{pq}\langle p|s_+s_-|q\rangle \gamma_{qp} = n_\alpha\]

    2) different electrons \(S_+ S_- = \sum s_{i+} s_{j-}, (i\neq j)\). There are in total \(n(n-1)\) terms. As a two-particle operator,

    \[\langle S_+ S_- \rangle = \langle ij|s_+ s_-|ij\rangle - \langle ij|s_+ s_-|ji\rangle = -\langle i^\alpha|j^\beta\rangle \langle j^\beta|i^\alpha\rangle\]
  2. Similarly, for \(S_- S_+\)
    1. same electron

    \[\sum_i \langle s_{i-} s_{i+}\rangle = n_\beta\]
    1. different electrons

    \[\langle S_- S_+ \rangle = -\langle i^\beta|j^\alpha\rangle \langle j^\alpha|i^\beta\rangle\]
  3. For \(S_z^2\)
    1. same electron

    \[\langle s_z^2\rangle = \frac{1}{4}(n_\alpha + n_\beta)\]
    1. different electrons

    \[\begin{split}&\frac{1}{2}\sum_{ij}(\langle ij|2s_{z1}s_{z2}|ij\rangle -\langle ij|2s_{z1}s_{z2}|ji\rangle) \\ &=\frac{1}{4}(\langle i^\alpha|i^\alpha\rangle \langle j^\alpha|j^\alpha\rangle - \langle i^\alpha|i^\alpha\rangle \langle j^\beta|j^\beta\rangle - \langle i^\beta|i^\beta\rangle \langle j^\alpha|j^\alpha\rangle + \langle i^\beta|i^\beta\rangle \langle j^\beta|j^\beta\rangle) \\ &-\frac{1}{4}(\langle i^\alpha|j^\alpha\rangle \langle j^\alpha|i^\alpha\rangle + \langle i^\beta|j^\beta\rangle\langle j^\beta|i^\beta\rangle) \\ &=\frac{1}{4}(n_\alpha^2 - n_\alpha n_\beta - n_\beta n_\alpha + n_\beta^2) -\frac{1}{4}(n_\alpha + n_\beta) \\ &=\frac{1}{4}((n_\alpha-n_\beta)^2 - (n_\alpha+n_\beta))\end{split}\]

In total

\[\begin{split}\langle S^2\rangle &= \frac{1}{2} (n_\alpha-\sum_{ij}\langle i^\alpha|j^\beta\rangle \langle j^\beta|i^\alpha\rangle +n_\beta -\sum_{ij}\langle i^\beta|j^\alpha\rangle\langle j^\alpha|i^\beta\rangle) + \frac{1}{4}(n_\alpha-n_\beta)^2 \\\end{split}\]
Args:
moa list of 2 ndarrays

Occupied alpha and occupied beta orbitals

Kwargs:
sndarray

AO overlap

Returns:

A list of two floats. The first is the expectation value of S^2. The second is the corresponding 2S+1

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', charge=1, spin=1, verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.kernel()
-75.623975516256706
>>> mo = (mf.mo_coeff[0][:,mf.mo_occ[0]>0], mf.mo_coeff[1][:,mf.mo_occ[1]>0])
>>> print('S^2 = %.7f, 2S+1 = %.7f' % spin_square(mo, mol.intor('int1e_ovlp_sph')))
S^2 = 0.7570150, 2S+1 = 2.0070027

pyscf.scf.uhf_symm module#

Non-relativistic unrestricted Hartree-Fock with point group symmetry.

class pyscf.scf.uhf_symm.HF1e(mol)[source]#

Bases: SymAdaptedUHF

scf(*args)#

SCF main driver

Kwargs:
dm0ndarray

If given, it will be used as the initial guess density matrix

Examples:

>>> import numpy
>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1.1')
>>> mf = scf.hf.SCF(mol)
>>> dm_guess = numpy.eye(mol.nao_nr())
>>> mf.kernel(dm_guess)
converged SCF energy = -98.5521904482821
-98.552190448282104
class pyscf.scf.uhf_symm.SymAdaptedUHF(mol)[source]#

Bases: UHF

SCF base class. non-relativistic RHF.

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default equals to Mole.max_memory

chkfilestr

checkpoint file to save MOs, orbital energies etc. Writing to chkfile can be disabled if this attribute is set to None or False.

conv_tolfloat

converge threshold. Default is 1e-9

conv_tol_gradfloat

gradients converge threshold. Default is sqrt(conv_tol)

max_cycleint

max number of iterations. If max_cycle <= 0, SCF iteration will be skipped and the kernel function will compute only the total energy based on the initial guess. Default value is 50.

init_guessstr

initial guess method. It can be one of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘sap’, ‘chkfile’. Default is ‘minao’

sap_basisstr or dict

basis for SAP initial guess, either filename or path as str or internal format dictionary.

DIISDIIS class

The class to generate diis object. It can be one of diis.SCF_DIIS, diis.ADIIS, diis.EDIIS.

diisboolean or object of DIIS class defined in scf.diis.

Default is the object associated to the attribute self.DIIS. Set it to None/False to turn off DIIS. Note if this attribute is initialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS information (vector basis and error vector) will be held inside this object. When kernel function is called again, the old states (vector basis and error vector) will be reused.

diis_spaceint

DIIS space size. By default, 8 Fock matrices and errors vector are stored.

diis_dampfloat

DIIS damping factor. Default is 0.

diis_start_cycleint

The step to start DIIS. Default is 1.

diis_file: ‘str’

File to store DIIS vectors and error vectors.

level_shiftfloat or int

Level shift (in AU) for virtual space. Default is 0.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

Direct SCF cutoff threshold. Default is 1e-13.

callbackfunction(envs_dict) => None

callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current environment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

check_convergencefunction(envs) => bool

A hook for overloading convergence criteria in SCF iterations.

Saved results:

convergedbool

Whether the SCF iteration converged

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

cyclesint

The number of iteration cycles performed

Examples:

>>> mol = gto.M(atom='H 0 0 0; H 0 0 1.1', basis='cc-pvdz')
>>> mf = scf.hf.SCF(mol)
>>> mf.verbose = 0
>>> mf.level_shift = .4
>>> mf.scf()
-1.0811707843775884
Attributes for UHF:
nelec(int, int)

If given, freeze the number of (alpha,beta) electrons to the given value.

level_shiftnumber or two-element list

level shift (in Eh) for alpha and beta Fock if two-element list is given.

init_guess_breaksymint

This configuration controls the algorithm used to break the spin symmetry of the initial guess: - 0 to disable symmetry breaking in the initial guess. - 1 to use the default algorithm introduced in pyscf-1.7. - 2 to adjust the num. electrons for spin-up and spin-down density matrices (issue #1839).

Examples:

>>> mol = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='ccpvdz', charge=1, spin=1, verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.kernel()
-75.623975516256706
>>> print('S^2 = %.7f, 2S+1 = %.7f' % mf.spin_square())
S^2 = 0.7570150, 2S+1 = 2.0070027
Attributes for symmetry allowed UHF:
irrep_nelecdict

Specify the number of alpha/beta electrons for particular irrep {‘ir_name’:(int,int), …}. For the irreps not listed in these dicts, the program will choose the occupancy based on the orbital energies.

Examp