10.25. scf — Self-consistent field methods

The scf module implements restricted and unrestricted, closed shell and open shell Hartree-Fock methods.

10.25.2. Addons

Special treatments may be required to the SCF methods in some situations. These special treatments cannot be universally applied for all SCF models. They were defined in the scf.addons module. For example, in an UHF calculation, we may want the \(S_z\) value to be changed (the numbers of alpha and beta electrons not conserved) during SCF iteration while conserving the total number of electrons. scf.addons.dynamic_sz_() can provide this functionality:

from pyscf import gto, scf
mol = gto.M(atom='O 0 0 0; O 0 0 1')
mf = scf.UHF(mol)
mf.verbose=4
mf = scf.addons.dynamic_sz_(mf)
mf.kernel()
print('S^2 = %s, 2S+1 = %s' % mf.spin_square())

This function automatically converges the ground sate of oxygen molecule to triplet state although we didn’t specify spin state in the mol object.

Note

Function scf.addons.dynamic_sz_() has side effects. It changes the underlying mean-field object.

The addons mechanism increases the flexibility of PySCf program. You can define various addons to customize the default behaviour of pyscf program. For example, if you’d like to track the changes of the density (the diagonal term of density matrix) of certain basis during the SCF iteration, you can write the following addon to output the required density:

def output_density(mf, basis_label):
    ao_labels = mf.mol.ao_labels()
    old_make_rdm1 = mf.make_rdm1
    def make_rdm1(mo_coeff, mo_occ):
        dm = old_make_rdm1(mo_coeff, mo_occ)
        print('AO         alpha             beta')
        for i,s in enumerate(ao_labels):
            if basis_label in s:
                print(s, dm[0][i,i], dm[1][i,i])
        return dm
    mf.make_rdm1 = make_rdm1
    return mf
from pyscf import gto, scf
mol = gto.M(atom='O 0 0 0; O 0 0 1')
mf = scf.UHF(mol)
mf.verbose=4
mf = scf.addons.dynamic_sz_(mf)
mf = output_density(mf, 'O 2p')
mf.kernel()

10.25.3. Caching two-electron integrals

When memory is enough (specified by the max_memory of SCF object), the SCF object generates all two-electron integrals in memory and cache them in _eri array. The default max_memory (defined in lib.parameters.MAX_MEMORY, see max_mem) is 4 GB. It roughly corresponds to two-electron real integrals for 250 orbitals. For small systems, the cached integrals usually provide the best performance. If you have enough main memory in your computer, you can increase the max_memory of SCF object to cache the integrals in memory.

The cached integrals _eri are treated as a dense tensor. When system becomes larger and the two-electron integral tensor becomes sparse, caching integrals may lose performance advantage. This is mainly due to the fact that the implementation of J/K build for the cached integrals did not utilize the sparsity of the integral tensor. Also, the data locality was not considered in the implementation which sometimes leads to bad OpenMP multi-threading speed up. For large system, the AO-driven direct SCF method is more favorable.

10.25.4. Customizing Hamiltonian

This integral object _eri is not merely used by the mean-field calculation. Along with the get_hcore() method, this two-electron integral object will be treated as the Hamiltonian in the post-SCF code whenever possible. This mechanism provides a way to model arbitrary fermion system in PySCF. You can customize a system by changing the 1-electron Hamiltonian and the mean-field _eri attribute. For example, the following code solves a model system:

import numpy
from pyscf import gto, scf, ao2mo, ccsd
mol = gto.M()
n = 10
mol.nelectron = 10
mf = scf.RHF(mol)

t = numpy.zeros((n,n))
for i in range(n-1):
    t[i,i+1] = t[i+1,i] = -1.0
t[n-1,0] = t[0,n-1] = 1.0  # anti-PBC
eri = numpy.zeros((n,n,n,n))
for i in range(n):
    eri[i,i,i,i] = 4.0

mf.get_hcore = lambda *args: t
mf.get_ovlp = lambda *args: numpy.eye(n)
# ao2mo.restore(8, eri, n) to get 8-fold symmetry of the integrals
# ._eri only supports the 2-electron integrals in 4-fold or 8-fold symmetry.
mf._eri = ao2mo.restore(8, eri, n)

mf.kernel()

mycc = ccsd.RCCSD(mf).run()
e,v = mycc.ipccsd(nroots=3)
print('IP = ', e)
e,v = mycc.eaccsd(nroots=3)
print('EA = ', e)

Some post-SCF methods require the 4-index MO integrals. Depending the available memory (affected by the value of max_memory in each class), these methods may not use the “AO integrals” cached in _eri. To ensure the post mean-field methods to use the _eri integrals no matter whether the actual memory usage is over the max_memory limite, you can set the flag incore_anyway in Mole class to True before calling the kernel() function of the post-SCF methods. In the following example, without setting incore_anyway=True, the CCSD calculations may crash:

import numpy
from pyscf import gto, scf, ao2mo, ccsd
mol = gto.M()
n = 10
mol.nelectron = n
mol.max_memory = 0
mf = scf.RHF(mol)

t = numpy.zeros((n,n))
for i in range(n-1):
    t[i,i+1] = t[i+1,i] = 1.0
t[n-1,0] = t[0,n-1] = -1.0
eri = numpy.zeros((n,n,n,n))
for i in range(n):
    eri[i,i,i,i] = 4.0

mf.get_hcore = lambda *args: t
mf.get_ovlp = lambda *args: numpy.eye(n)
mf._eri = ao2mo.restore(8, eri, n)
mf.kernel()

mol.incore_anyway = True
mycc = ccsd.RCCSD(mf).run()
e,v = mycc.ipccsd(nroots=3)
print('IP = ', e)
e,v = mycc.eaccsd(nroots=3)
print('EA = ', e)

Holding the entire two-particle interactions matrix elements in memory often leads to high memory usage. In the SCF calculation, the memory usage can be optimized if _eri is sparse. The SCF iterations requires only the Fock matrix which in turn calls the J/K build function SCF.get_jk() to compute the Coulomb and HF-exchange matrix. Overwriting the SCF.get_jk() function can reduce the memory footprint of the SCF part in the above example:

import numpy
from pyscf import gto, scf, ao2mo, ccsd
mol = gto.M()
n = 10
mol.nelectron = n
mol.max_memory = 0
mf = scf.RHF(mol)

t = numpy.zeros((n,n))
for i in range(n-1):
    t[i,i+1] = t[i+1,i] = 1.0
t[n-1,0] = t[0,n-1] = -1.0

mf.get_hcore = lambda *args: t
mf.get_ovlp = lambda *args: numpy.eye(n)
def get_jk(mol, dm, *args):
    j = numpy.diag(dm.diagonal()) * 4.
    k = numpy.diag(dm.diagonal()) * 4.
    return j, k
mf.get_jk = get_jk
mf.kernel()

Another way to handle the two-particle interactions of large model system is to use the density fitting/Cholesky decomposed integrals. See also Saving/Loading DF integral tensor.

10.25.5. Program reference

10.25.5.1. Hartree-Fock

Simple usage:

>>> from pyscf import gto, scf
>>> mol = gto.M(atom='H 0 0 0; H 0 0 1')
>>> mf = scf.RHF(mol).run()

scf.RHF() returns an instance of SCF class. There are some parameters to control the SCF method.

verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

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

chkfilestr

checkpoint file to save MOs, orbital energies etc.

conv_tolfloat

converge threshold. Default is 1e-10

max_cycleint

max number of iterations. Default is 50

init_guessstr

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

DIISclass listed in scf.diis

Default is diis.SCF_DIIS. Set it to None/False to turn off DIIS.

diisbool

whether to do DIIS. Default is True.

diis_spaceint

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

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.

direct_scfbool

Direct SCF is used by default.

direct_scf_tolfloat

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

callbackfunction

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 envrionment.

conv_checkbool

An extra cycle to check convergence after SCF iterations.

nelec(int,int), for UHF/ROHF class

freeze the number of (alpha,beta) electrons.

irrep_nelecdict, for symmetry- RHF/ROHF/UHF class only

to indicate the number of electrons for each irreps. In RHF, give {‘ir_name’:int, …} ; In ROHF/UHF, give {‘ir_name’:(int,int), …} . It is effective when Mole.symmetry is set True.

auxbasisstr, for density fitting SCF only

Auxiliary basis for density fitting.

>>> mf = scf.density_fit(scf.UHF(mol))
>>> mf.scf()

Density fitting can be applied to all non-relativistic HF class.

with_ssssbool, for Dirac-Hartree-Fock only

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

with_gauntbool, for Dirac-Hartree-Fock only

If False, ignore Gaunt interaction. Default is False.

Saved results

convergedbool

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

pyscf.scf.DHF(mol, *args)[source]

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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, 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
pyscf.scf.GHF(mol, *args)[source]

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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.

pyscf.scf.HF(mol, *args)[source]

A wrap function to create SCF class (RHF or 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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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
pyscf.scf.RHF(mol, *args)[source]

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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
pyscf.scf.ROHF(mol, *args)[source]

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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
pyscf.scf.UHF(mol, *args)[source]

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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_breaksymlogical

If given, overwrite BREAKSYM.

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
pyscf.scf.X2C(mol, *args)[source]

X2C UHF (in testing)

10.25.5.2. Non-relativistic Hartree-Fock

Hartree-Fock

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

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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)
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 eletron numbers.

natorbbool

Whether to transform natural orbital in active space. Be cautious of this parameter when CASCI/CASSCF are combined with DMRG solver or selected CI solver because DMRG and selected CI are not invariant to the rotation in active space. False by default.

canonicalizationbool

Whether to canonicalize orbitals. Note that canonicalization does not change the orbitals in active space by default. It only diagonalizes core and external space of the general Fock matirx. To get the natural orbitals in active space, attribute natorb need 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 orbitls.

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
CASSCF(*args, **kwargs)

CASCI

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 eletron numbers.

natorbbool

Whether to transform natural orbital in active space. Be cautious of this parameter when CASCI/CASSCF are combined with DMRG solver or selected CI solver because DMRG and selected CI are not invariant to the rotation in active space. False by default.

canonicalizationbool

Whether to canonicalize orbitals. Note that canonicalization does not change the orbitals in active space by default. It only diagonalizes core and external space of the general Fock matirx. To get the natural orbitals in active space, attribute natorb need 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 orbitls.

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
CASSCF

Extra attributes for CASSCF:

conv_tolfloat

Converge threshold. Default is 1e-7

conv_tol_gradfloat

Converge threshold for CI gradients and orbital rotation gradients. Default is 1e-4

max_stepsizefloat

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

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 3.

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 0.2.

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 2.

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 envrionment.

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 specified, the active segment of the mo_coeff is natural orbitls.

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
Gradients(*args, **kwargs)

Non-relativistic restricted Hartree-Fock gradients

Hessian(*args, **kwargs)

Non-relativistic restricted Hartree-Fock hessian

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 stablize 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()
>>> # 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_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-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.

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_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-seperated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-seperated 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)[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.

Returns:

New orbitals that are more close to the stable condition. The return value includes two set of orbitals. The first corresponds to the internal stability and the second corresponds to the external stability.

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

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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 stablize 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()
>>> # freeze 2 core orbitals and 3 high lying unoccupied orbitals
>>> mycc.set(frozen = [0,1,16,17,18]).run()

Saved results

convergedbool

CCSD converged or not

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)

DIIS

alias of pyscf.scf.diis.CDIIS

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
canonicalize(mo_coeff, mo_occ, fock=None)

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

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

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-partical 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_tot(dm=None, h1e=None, vhf=None)

Total Hartree-Fock energy, electronic part plus nuclear repulstion 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

Returns:

Density matrix, 2D ndarray

from_fcidump(filename, molpro_orbsym=False)

Update the SCF object with the quantities defined in FCIDUMP file

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)

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_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-seperated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-seperated 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_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_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

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]])
kernel(*args, **kwargs)

An alias to method scf Function Signature: kernel(self, 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
make_rdm1(mo_coeff=None, mo_occ=None, **kwargs)[source]

One-particle density matrix in AO representation

Args:
mo_coeff2D ndarray

Orbital coefficients. Each column is one orbital.

mo_occ1D ndarray

Occupancy

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 carreid 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 carreid 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

nuc_grad_method()[source]

Hook to create object for analytical nuclear gradients.

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 carreid 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

reset(mol=None)[source]

Reset mol and relevant attributes associated to the old mol object

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
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_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_.

pyscf.scf.hf.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; 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.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.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.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-partical 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 repulstion 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)[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')[source]

Generate density matrix for initial guess

Kwargs:
keystr

One of ‘minao’, ‘atom’, ‘huckel’, ‘hcore’, ‘1e’, ‘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-seperated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-seperated 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

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.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 envrionment.

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

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 carreid 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 carreid 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.uniq_var_indices(mo_occ)[source]

Indicies 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.

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]
scf(*args)[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.hf_symm.RHF

alias of pyscf.scf.hf_symm.SymAdaptedRHF

pyscf.scf.hf_symm.ROHF

alias of pyscf.scf.hf_symm.SymAdaptedROHF

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

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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}
CASSCF(*args, **kwargs)

CASCI

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 eletron numbers.

natorbbool

Whether to transform natural orbital in active space. Be cautious of this parameter when CASCI/CASSCF are combined with DMRG solver or selected CI solver because DMRG and selected CI are not invariant to the rotation in active space. False by default.

canonicalizationbool

Whether to canonicalize orbitals. Note that canonicalization does not change the orbitals in active space by default. It only diagonalizes core and external space of the general Fock matirx. To get the natural orbitals in active space, attribute natorb need 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 orbitls.

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
CASSCF

Extra attributes for CASSCF:

conv_tolfloat

Converge threshold. Default is 1e-7

conv_tol_gradfloat

Converge threshold for CI gradients and orbital rotation gradients. Default is 1e-4

max_stepsizefloat

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

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 3.

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 0.2.

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 2.

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 envrionment.

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 specified, the active segment of the mo_coeff is natural orbitls.

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

canonicalize(mo_coeff, mo_occ, fock=None)

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

eig(h, s)

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.

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

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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)}
CASSCF(*args, **kwargs)

CASCI

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 eletron numbers.

natorbbool

Whether to transform natural orbital in active space. Be cautious of this parameter when CASCI/CASSCF are combined with DMRG solver or selected CI solver because DMRG and selected CI are not invariant to the rotation in active space. False by default.

canonicalizationbool

Whether to canonicalize orbitals. Note that canonicalization does not change the orbitals in active space by default. It only diagonalizes core and external space of the general Fock matirx. To get the natural orbitals in active space, attribute natorb need 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 orbitls.

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
CASSCF

Extra attributes for CASSCF:

conv_tolfloat

Converge threshold. Default is 1e-7

conv_tol_gradfloat

Converge threshold for CI gradients and orbital rotation gradients. Default is 1e-4

max_stepsizefloat

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

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 3.

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 0.2.

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 2.

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 envrionment.

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 specified, the active segment of the mo_coeff is natural orbitls.

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; Mulliken population analysis

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

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

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_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])
make_rdm1(mo_coeff=None, mo_occ=None, **kwargs)[source]

One-particle densit matrix. mo_occ is a 1D array, with occupancy 1 or 2.

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.eig(mf, h, 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.

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.so2ao_mo_coeff(so, irrep_mo_coeff)[source]

Transfer the basis of MO coefficients, from spin-adapted basis to AO basis

Restricted Open-shell Hartree-Fock

class pyscf.scf.rohf.HF1e(mol)[source]
scf(*args)[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
class pyscf.scf.rohf.ROHF(mol)[source]

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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)
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 eletron numbers.

natorbbool

Whether to transform natural orbital in active space. Be cautious of this parameter when CASCI/CASSCF are combined with DMRG solver or selected CI solver because DMRG and selected CI are not invariant to the rotation in active space. False by default.

canonicalizationbool

Whether to canonicalize orbitals. Note that canonicalization does not change the orbitals in active space by default. It only diagonalizes core and external space of the general Fock matirx. To get the natural orbitals in active space, attribute natorb need 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 orbitls.

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
CASSCF(*args, **kwargs)

CASCI

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 eletron numbers.

natorbbool

Whether to transform natural orbital in active space. Be cautious of this parameter when CASCI/CASSCF are combined with DMRG solver or selected CI solver because DMRG and selected CI are not invariant to the rotation in active space. False by default.

canonicalizationbool

Whether to canonicalize orbitals. Note that canonicalization does not change the orbitals in active space by default. It only diagonalizes core and external space of the general Fock matirx. To get the natural orbitals in active space, attribute natorb need 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 orbitls.

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
CASSCF

Extra attributes for CASSCF:

conv_tolfloat

Converge threshold. Default is 1e-7

conv_tol_gradfloat

Converge threshold for CI gradients and orbital rotation gradients. Default is 1e-4

max_stepsizefloat

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

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 3.

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 0.2.

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 2.

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 envrionment.

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 specified, the active segment of the mo_coeff is natural orbitls.

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
Gradients(*args, **kwargs)

Non-relativistic ROHF gradients

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.

eig(fock, s)[source]

Solver for generalized eigenvalue problem

\[HC = SCE\]
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)

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_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

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]])
make_rdm1(mo_coeff=None, mo_occ=None, **kwargs)[source]

One-particle densit matrix. mo_occ is a 1D array, with occupancy 1 or 2.

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)[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.

Returns:

The return value includes two set of orbitals which are more close to the required stable condition.

pyscf.scf.rohf.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

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.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)[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.make_rdm1(mo_coeff, mo_occ, **kwargs)[source]

One-particle densit matrix. mo_occ is a 1D array, with occupancy 1 or 2.

class pyscf.scf.uhf.HF1e(mol)[source]
scf(*args)[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
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]

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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_breaksymlogical

If given, overwrite BREAKSYM.

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
Gradients(*args, **kwargs)

Non-relativistic unrestricted Hartree-Fock gradients

Hessian(*args, **kwargs)

Non-relativistic UHF hessian

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

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_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_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_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, breaksym=True)[source]

Generate initial guess density matrix from core hamiltonian

Returns:

Density matrix, 2D ndarray

init_guess_by_atom(mol=None, breaksym=True)[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

Returns:

Density matrix, 2D ndarray

init_guess_by_huckel(mol=None, breaksym=True)[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=True)[source]

Initial guess in terms of the overlap to minimal basis.

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

Returns:

A list of 2D ndarrays for alpha and 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 carreid 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 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)[source]

Stability analysis for RHF/RKS 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.

Returns:

New orbitals that are more close to the stable condition. The return value includes two set of orbitals. The first corresponds to the internal stability and the second corresponds to the external stability.

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_grad(mo_coeff, mo_occ, fock_ao)[source]

UHF Gradients

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_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 swith 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_minao(mol, breaksym=True)[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.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

Returns:

A list of 2D ndarrays for alpha and 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

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

class pyscf.scf.uhf_symm.HF1e(mol)[source]
scf(*args)[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
class pyscf.scf.uhf_symm.SymAdaptedUHF(mol)[source]

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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_breaksymlogical

If given, overwrite BREAKSYM.

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.

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.623975516256692
>>> 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.429189192031131
>>> mf.get_irrep_nelec()
{'A1': (3, 3), 'A2': (0, 0), 'B1': (1, 0), 'B2': (1, 1)}
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 in occupied, virtual subspaces separatedly (without change occupancy).

eig(h, 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_irrep_nelec(mol=None, mo_coeff=None, mo_occ=None, s=None)[source]

Alpha/beta electron numbers for each irreducible representation.

Args:
molan instance of Mole

To provide irrep_id, and spin-adapted basis

mo_occa list of 1D ndarray

Regular occupancy, without grouping for irreps

mo_coeffa list of 2D ndarray

Regular orbital coefficients, without grouping for irreps

Returns:
irrep_nelecdict

The number of alpha/beta electrons for each irrep {‘ir_name’:(int,int), …}.

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.UHF(mol)
>>> mf.scf()
-75.623975516256721
>>> scf.uhf_symm.get_irrep_nelec(mol, mf.mo_coeff, mf.mo_occ)
{'A1': (3, 3), 'A2': (0, 0), 'B1': (1, 1), 'B2': (1, 0)}
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.

pyscf.scf.uhf_symm.UHF

alias of pyscf.scf.uhf_symm.SymAdaptedUHF

pyscf.scf.uhf_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.uhf_symm.get_irrep_nelec(mol, mo_coeff, mo_occ, s=None)[source]

Alpha/beta electron numbers for each irreducible representation.

Args:
molan instance of Mole

To provide irrep_id, and spin-adapted basis

mo_occa list of 1D ndarray

Regular occupancy, without grouping for irreps

mo_coeffa list of 2D ndarray

Regular orbital coefficients, without grouping for irreps

Returns:
irrep_nelecdict

The number of alpha/beta electrons for each irrep {‘ir_name’:(int,int), …}.

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.UHF(mol)
>>> mf.scf()
-75.623975516256721
>>> scf.uhf_symm.get_irrep_nelec(mol, mf.mo_coeff, mf.mo_occ)
{'A1': (3, 3), 'A2': (0, 0), 'B1': (1, 1), 'B2': (1, 0)}
class pyscf.scf.atom_hf.AtomHF1e(mol)[source]
eig(f, s)

Solver for generalized eigenvalue problem

\[HC = SCE\]
class pyscf.scf.atom_hf.AtomSphericAverageRHF(mol)[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

Non-relativistic generalized Hartree-Fock

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

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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.

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', 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

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_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-seperated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-seperated 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
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

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]])
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 carreid 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
pyscf.scf.ghf.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

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.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-seperated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-seperated 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.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 swith 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

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

class pyscf.scf.ghf_symm.GHF(mol)[source]

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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.

canonicalize(mo_coeff, mo_occ, fock=None)

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

eig(h, 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_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.

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).

Restricted coupled pertubed Hartree-Fock solver

pyscf.scf.cphf.kernel(fvind, mo_energy, mo_occ, h1, s1=None, max_cycle=20, tol=1e-09, hermi=False, verbose=2)
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.

pyscf.scf.cphf.solve(fvind, mo_energy, mo_occ, h1, s1=None, max_cycle=20, tol=1e-09, hermi=False, verbose=2)[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.

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

For field independent basis. First order overlap matrix is zero

pyscf.scf.cphf.solve_withs1(fvind, mo_energy, mo_occ, h1, s1, max_cycle=20, tol=1e-09, hermi=False, verbose=2)[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:
hermiboolean

Whether the matrix defined by fvind is Hermitian or not.

Returns:

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

Unrestricted coupled pertubed Hartree-Fock solver

pyscf.scf.ucphf.kernel(fvind, mo_energy, mo_occ, h1, s1=None, max_cycle=20, tol=1e-09, hermi=False, verbose=2)
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=20, tol=1e-09, hermi=False, verbose=2)[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)[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=20, tol=1e-09, hermi=False, verbose=2)[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]

10.25.5.3. Relativistic Hartree-Fock

Dirac Hartree-Fock

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

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 skiped and the kernel function will compute only the total energy based on the intial guess. Default value is 50.

init_guessstr

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

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 inialized as a DIIS object, the SCF driver will use this object in the iteration. The DIIS informations (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_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 envrionment.

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

SCF converged or not

e_totfloat

Total HF energy (electronic energy plus nuclear repulsion)

mo_energy :

Orbital energies

mo_occ

Orbital occupancy

mo_coeff

Orbital coefficients

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, 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

analyze(verbose=None)[source]

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

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

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_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-seperated Coulomb operator: erf( omega * r12 ) / r12. If specified, integration are evaluated based on the long-range part of the range-seperated 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_veff(mol=None, dm=None, dm_last=0, vhf_last=0, hermi=1)[source]

Dirac-Coulomb

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

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.

make_rdm1(mo_coeff=None, mo_occ=None, **kwargs)[source]

One-particle density matrix in AO representation

Args:
mo_coeff2D ndarray

Orbital coefficients. Each column is one orbital.

mo_occ1D ndarray

Occupancy

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
class pyscf.scf.dhf.HF1e(mol)[source]
scf(*args)[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
class pyscf.scf.dhf.RDHF(mol)[source]

Kramers restricted Dirac-Hartree-Fock

pyscf.scf.dhf.RHF

alias of pyscf.scf.dhf.RDHF

pyscf.scf.dhf.UDHF

alias of pyscf.scf.dhf.DHF

pyscf.scf.dhf.UHF

alias of pyscf.scf.dhf.DHF

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.get_grad(mo_coeff, mo_occ, fock_ao)[source]

DHF Gradients

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

Generate density matrix for initial guess

Kwargs:
keystr

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

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 swith 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.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)]^*

10.25.5.4. Addons

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_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.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. It can be applied to unrestricted HF/KS and restricted open-shell HF/KS.

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

Use maximum overlap method to determine occupation number for each orbital in every iteration. It can be applied to unrestricted HF/KS and restricted open-shell HF/KS.

pyscf.scf.addons.partial_cholesky_orth_(S, cholthr=1e-09, canthr=1e-07)[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)
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)[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.chkfile.dump_scf(mol, chkfile, e_tot, mo_energy, mo_coeff, mo_occ, overwrite_mol=True)[source]

save temporary results

DIIS

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

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

update(s, d, f, mf, h1e, vhf)[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)[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 pyscf.scf.diis.CDIIS

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

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

update(s, d, f, mf, h1e, vhf)[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 pyscf.scf.diis.CDIIS

pyscf.scf.diis.SCF_DIIS

alias of pyscf.scf.diis.CDIIS

pyscf.scf.diis.get_err_vec(s, d, f)[source]

error vector = SDF - FDS

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))

10.25.5.5. Stability analysis

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.rhf_stability(mf, internal=True, external=False, verbose=None)[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.

Returns:

New orbitals that are more close to the stable condition. The return value includes two set of orbitals. The first corresponds to the internal stability and the second corresponds to the external stability.

pyscf.scf.stability.rohf_stability(mf, internal=True, external=False, verbose=None)[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.

Returns:

The return value includes two set of orbitals which are more close to the required stable condition.

pyscf.scf.stability.uhf_stability(mf, internal=True, external=False, verbose=None)[source]

Stability analysis for RHF/RKS 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.

Returns:

New orbitals that are more close to the stable condition. The return value includes two set of orbitals. The first corresponds to the internal stability and the second corresponds to the external stability.

Wave Function Stability Analysis

Ref. JCP, 66, 3045 (1977); DOI:10.1063/1.434318 JCP 104, 9047 (1996); DOI:10.1063/1.471637