pyscf.lo package¶
Submodules¶
pyscf.lo.boys module¶
Foster-Boys localization
-
pyscf.lo.boys.
BF
¶ alias of
pyscf.lo.boys.Boys
-
class
pyscf.lo.boys.
Boys
(mol, mo_coeff=None)¶ Bases:
pyscf.soscf.ciah.CIAHOptimizer
The Foster-Boys localization optimizer that maximizes the orbital dipole
- Args:
mol : Mole object
- Kwargs:
- mo_coeffsize (N,N) np.array
The orbital space to localize for Boys localization. When initializing the localization optimizer
bopt = Boys(mo_coeff)
,Note these orbitals
mo_coeff
may or may not be used as initial guess, depending on the attribute.init_guess
. If.init_guess
is set to None, themo_coeff
will be used as initial guess. If.init_guess
is ‘atomic’, a few atomic orbitals will be constructed inside the space of the input orbitals and the atomic orbitals will be used as initial guess.Note when calling .kernel(orb) method with a set of orbitals as argument, the orbitals will be used as initial guess regardless of the value of the attributes .mo_coeff and .init_guess.
- Attributes for Boys class:
- 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 1e-6
- conv_tol_gradfloat
Converge threshold for orbital rotation gradients. Default 1e-3
- max_cycleint
The max. number of macro iterations. Default 100
- max_itersint
The max. number of iterations in each macro iteration. Default 20
- max_stepsizefloat
The step size for orbital rotation. Small step (0.005 - 0.05) is prefered. Default 0.03.
- init_guessstr or None
Initial guess for optimization. If set to None, orbitals defined by the attribute .mo_coeff will be used as initial guess. If set to ‘atomic’, atomic orbitals will be used as initial guess. Default is ‘atomic’
Saved results
- mo_coeffndarray
Localized orbitals
-
ah_max_cycle
= 40¶
-
ah_start_tol
= 1000000000.0¶
-
ah_trust_region
= 3¶
-
conv_tol
= 1e-06¶
-
conv_tol_grad
= None¶
-
cost_function
(u=None)¶
-
dump_flags
(verbose=None)¶
-
gen_g_hop
(u)¶
-
get_grad
(u=None)¶
-
get_init_guess
(key='atomic')¶ Generate initial guess for localization.
- Kwargs:
- keystr or bool
If key is ‘atomic’, initial guess is based on the projected atomic orbitals. False
-
init_guess
= 'atomic'¶
-
kernel
(mo_coeff=None, callback=None, verbose=None)¶ Kernel function is the main driver of a method. Every method should define the kernel function as the entry of the calculation. Note the return value of kernel function is not strictly defined. It can be anything related to the method (such as the energy, the wave-function, the DFT mesh grids etc.).
-
max_cycle
= 100¶
-
max_iters
= 20¶
-
max_stepsize
= 0.05¶
-
pyscf.lo.boys.
FB
¶ alias of
pyscf.lo.boys.Boys
-
pyscf.lo.boys.
atomic_init_guess
(mol, mo_coeff)¶
-
pyscf.lo.boys.
dipole_integral
(mol, mo_coeff)¶
-
pyscf.lo.boys.
kernel
(localizer, mo_coeff=None, callback=None, verbose=None)¶
pyscf.lo.edmiston module¶
Edmiston-Ruedenberg localization
-
pyscf.lo.edmiston.
ER
¶ alias of
pyscf.lo.edmiston.EdmistonRuedenberg
-
pyscf.lo.edmiston.
Edmiston
¶ alias of
pyscf.lo.edmiston.EdmistonRuedenberg
pyscf.lo.iao module¶
Intrinsic Atomic Orbitals ref. JCTC, 9, 4834
-
pyscf.lo.iao.
fast_iao_mullikan_pop
(mol, dm, iaos, verbose=5)¶ - Args:
mol : the molecule or cell object
- iaos2D array
(orthogonal or non-orthogonal) IAO orbitals
- Returns:
mullikan population analysis in the basis IAO
-
pyscf.lo.iao.
iao
(mol, orbocc, minao='minao', kpts=None)¶ Intrinsic Atomic Orbitals. [Ref. JCTC, 9, 4834]
- Args:
mol : the molecule or cell object
- orbocc2D array
occupied orbitals
- Returns:
non-orthogonal IAO orbitals. Orthogonalize them as C (C^T S C)^{-1/2}, eg using
orth.lowdin()
>>> orbocc = mf.mo_coeff[:,mf.mo_occ>0] >>> c = iao(mol, orbocc) >>> numpy.dot(c, orth.lowdin(reduce(numpy.dot, (c.T,s,c))))
-
pyscf.lo.iao.
reference_mol
(mol, minao='minao')¶ Create a molecule which uses reference minimal basis
pyscf.lo.ibo module¶
Intrinsic Bonding Orbitals ref. JCTC, 9, 4834
Below here is work done by Paul Robinson. much of the code below is adapted from code published freely on the website of Gerald Knizia Ref: JCTC, 2013, 9, 4834-4843
-
pyscf.lo.ibo.
MakeAtomIbOffsets
(Atoms)¶ calcualte offset of first orbital of individual atoms in the valence minimal basis (IB)
-
pyscf.lo.ibo.
MakeAtomInfos
()¶
-
pyscf.lo.ibo.
PM
(mol, orbocc, iaos, s, exponent)¶ Note this localization is slightly different to Knizia’s implementation. The localization here reserves orthogonormality during optimization. Orbitals are projected to IAO basis first and the Mulliken pop is calculated based on IAO basis (in function atomic_pops). A series of unitary matrices are generated and applied on the input orbitals. The intemdiate orbitals in the optimization and the finally localized orbitals are all orthogonormal.
Examples:
>>> from pyscf import gto, scf >>> from pyscf.lo import ibo >>> mol = gto.M(atom='H 0 0 0; F 0 0 1', >>> basis='unc-sto3g') >>> mf = scf.RHF(mol).run() >>> pm = ibo.PM(mol, mf.mo_coeff[:,mf.mo_occ>0]) >>> loc_orb = pm.kernel()
-
pyscf.lo.ibo.
Pipek
(mol, orbocc, iaos, s, exponent)¶ Note this localization is slightly different to Knizia’s implementation. The localization here reserves orthogonormality during optimization. Orbitals are projected to IAO basis first and the Mulliken pop is calculated based on IAO basis (in function atomic_pops). A series of unitary matrices are generated and applied on the input orbitals. The intemdiate orbitals in the optimization and the finally localized orbitals are all orthogonormal.
Examples:
>>> from pyscf import gto, scf >>> from pyscf.lo import ibo >>> mol = gto.M(atom='H 0 0 0; F 0 0 1', >>> basis='unc-sto3g') >>> mf = scf.RHF(mol).run() >>> pm = ibo.PM(mol, mf.mo_coeff[:,mf.mo_occ>0]) >>> loc_orb = pm.kernel()
-
pyscf.lo.ibo.
PipekMezey
(mol, orbocc, iaos, s, exponent)¶ Note this localization is slightly different to Knizia’s implementation. The localization here reserves orthogonormality during optimization. Orbitals are projected to IAO basis first and the Mulliken pop is calculated based on IAO basis (in function atomic_pops). A series of unitary matrices are generated and applied on the input orbitals. The intemdiate orbitals in the optimization and the finally localized orbitals are all orthogonormal.
Examples:
>>> from pyscf import gto, scf >>> from pyscf.lo import ibo >>> mol = gto.M(atom='H 0 0 0; F 0 0 1', >>> basis='unc-sto3g') >>> mf = scf.RHF(mol).run() >>> pm = ibo.PM(mol, mf.mo_coeff[:,mf.mo_occ>0]) >>> loc_orb = pm.kernel()
-
pyscf.lo.ibo.
ibo
(mol, orbocc, locmethod='IBO', iaos=None, s=None, exponent=4, grad_tol=1e-08, max_iter=200, verbose=3)¶ Intrinsic Bonding Orbitals
This function serves as a wrapper to the underlying localization functions ibo_loc and PipekMezey to create IBOs.
- Args:
mol : the molecule or cell object
orbocc : occupied molecular orbital coefficients
- Kwargs:
- locmethodstring
the localization method ‘PM’ for Pipek Mezey localization or ‘IBO’ for the IBO localization
- iaos2D array
the array of IAOs
- s2D array
the overlap array in the ao basis
- Returns:
IBOs in the basis defined in mol object.
-
pyscf.lo.ibo.
ibo_loc
(mol, orbocc, iaos, s, exponent, grad_tol, max_iter, verbose=3)¶ Intrinsic Bonding Orbitals. [Ref. JCTC, 9, 4834]
This implementation follows Knizia’s implementation execept that the resultant IBOs are symmetrically orthogonalized. Note the IBOs of this implementation do not strictly maximize the IAO Mulliken charges.
IBOs can also be generated by another implementation (see function pyscf.lo.ibo.PM). In that function, PySCF builtin Pipek-Mezey localization module was used to maximize the IAO Mulliken charges.
- Args:
mol : the molecule or cell object
- orbocc2D array or a list of 2D array
occupied molecular orbitals or crystal orbitals for each k-point
- Kwargs:
- iaos2D array
the array of IAOs
- exponentinteger
Localization power in PM scheme
- grad_tolfloat
convergence tolerance for norm of gradients
- Returns:
IBOs in the big basis (the basis defined in mol object).
-
pyscf.lo.ibo.
shell_str
(l, n_cor, n_val)¶ Help function to define core and valence shells for shell with different l
pyscf.lo.nao module¶
Natural atomic orbitals Ref:
Weinhold et al., J. Chem. Phys. 83(1985), 735-746
-
pyscf.lo.nao.
nao
(mol, mf, s=None, restore=True)¶
-
pyscf.lo.nao.
prenao
(mol, dm)¶
-
pyscf.lo.nao.
set_atom_conf
(element, description)¶ Change the default atomic core and valence configuration to the one given by “description”. See data/elements.py for the default configuration.
- Args:
- elementstr or int
Element symbol or nuclear charge
- descriptionstr or a list of str
- “double p” : double p shell“double d” : double d shell“double f” : double f shell“polarize” : add one polarized shell“1s1d” : keep core unchanged and set 1 s 1 d shells for valence(“3s2p”,”1d”) : 3 s, 2 p shells for core and 1 d shells for valence
pyscf.lo.orth module¶
-
pyscf.lo.orth.
orth_ao
(mf_or_mol, method='meta_lowdin', pre_orth_ao=None, scf_method=None, s=None)¶ Orthogonalize AOs
- Kwargs:
- methodstr
One of | lowdin : Symmetric orthogonalization | meta-lowdin : Lowdin orth within core, valence, virtual space separately (JCTC, 10, 3784) | NAO
- pre_orth_ao: numpy.ndarray
Coefficients to restore AO characters for arbitrary basis. If not given, the coefficients are generated based on ANO basis. You can skip this by setting pre_orth_ao to identity matrix, e.g. np.eye(mol.nao).
-
pyscf.lo.orth.
pre_orth_ao
(mol, method='ANO')¶ Restore AO characters. Possible methods include the ANO/MINAO projection or fraction-averaged atomic RHF calculation
-
pyscf.lo.orth.
pre_orth_ao_atm_scf
(mol)¶
-
pyscf.lo.orth.
restore_ao_character
(mol, method='ANO')¶ Restore AO characters. Possible methods include the ANO/MINAO projection or fraction-averaged atomic RHF calculation
-
pyscf.lo.orth.
schmidt
(s)¶
-
pyscf.lo.orth.
vec_lowdin
(c, s=1)¶ lowdin orth for the metric c.T*s*c and get x, then c*x
-
pyscf.lo.orth.
vec_schmidt
(c, s=1)¶ schmidt orth for the metric c.T*s*c and get x, then c*x
pyscf.lo.pipek module¶
Pipek-Mezey localization
ref. JCTC, 10, 642 (2014); DOI:10.1021/ct401016x
-
pyscf.lo.pipek.
PM
¶ alias of
pyscf.lo.pipek.PipekMezey
-
pyscf.lo.pipek.
Pipek
¶ alias of
pyscf.lo.pipek.PipekMezey
-
class
pyscf.lo.pipek.
PipekMezey
(mol, mo_coeff=None, mf=None)¶ Bases:
pyscf.lo.boys.Boys
The Pipek-Mezey localization optimizer that maximizes the orbital population
- Args:
mol : Mole object
- Kwargs:
- mo_coeffsize (N,N) np.array
The orbital space to localize for PM localization. When initializing the localization optimizer
bopt = PM(mo_coeff)
,Note these orbitals
mo_coeff
may or may not be used as initial guess, depending on the attribute.init_guess
. If.init_guess
is set to None, themo_coeff
will be used as initial guess. If.init_guess
is ‘atomic’, a few atomic orbitals will be constructed inside the space of the input orbitals and the atomic orbitals will be used as initial guess.Note when calling .kernel(orb) method with a set of orbitals as argument, the orbitals will be used as initial guess regardless of the value of the attributes .mo_coeff and .init_guess.
- Attributes for PM class:
- 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 1e-6
- conv_tol_gradfloat
Converge threshold for orbital rotation gradients. Default 1e-3
- max_cycleint
The max. number of macro iterations. Default 100
- max_itersint
The max. number of iterations in each macro iteration. Default 20
- max_stepsizefloat
The step size for orbital rotation. Small step (0.005 - 0.05) is prefered. Default 0.03.
- init_guessstr or None
Initial guess for optimization. If set to None, orbitals defined by the attribute .mo_coeff will be used as initial guess. If set to ‘atomic’, atomic orbitals will be used as initial guess. Default ‘atomic’
- pop_methodstr
How the orbital population is calculated. By default, meta-lowdin population (JCTC, 10, 3784) is used. It can be set to ‘mulliken’, or ‘lowdin’ for other population definition
- exponentint
The power to define norm. It can be 2 or 4. Default 2.
Saved results
- mo_coeffndarray
Localized orbitals
-
atomic_pops
(mol, mo_coeff, method=None)¶ - Kwargs:
- methodstring
The atomic population projection scheme. It can be mulliken, lowdin, meta_lowdin, iao, or becke
- Returns:
A 3-index tensor [A,i,j] indicates the population of any orbital-pair density |i><j| for each species (atom in this case). This tensor is used to construct the population and gradients etc.
You can customize the PM localization wrt other population metric, such as the charge of a site, the charge of a fragment (a group of atoms) by overwriting this tensor. See also the example pyscf/examples/loc_orb/40-hubbard_model_PM_localization.py for the PM localization of site-based population for hubbard model.
-
conv_tol
= 1e-06¶
-
cost_function
(u=None)¶
-
dump_flags
(verbose=None)¶
-
exponent
= 2¶
-
gen_g_hop
(u)¶
-
get_grad
(u=None)¶
-
pop_method
= 'meta_lowdin'¶
-
pyscf.lo.pipek.
atomic_pops
(mol, mo_coeff, method='meta_lowdin', mf=None)¶ - Kwargs:
- methodstring
The atomic population projection scheme. It can be mulliken, lowdin, meta_lowdin, iao, or becke
- Returns:
A 3-index tensor [A,i,j] indicates the population of any orbital-pair density |i><j| for each species (atom in this case). This tensor is used to construct the population and gradients etc.
You can customize the PM localization wrt other population metric, such as the charge of a site, the charge of a fragment (a group of atoms) by overwriting this tensor. See also the example pyscf/examples/loc_orb/40-hubbard_model_PM_localization.py for the PM localization of site-based population for hubbard model.
pyscf.lo.vvo module¶
Valence Virtual Orbitals ref. 10.1021/acs.jctc.7b00493
-
pyscf.lo.vvo.
livvo
(mol, orbocc, orbvirt, locmethod='IBO', iaos=None, s=None, exponent=4, grad_tol=1e-08, max_iter=200, verbose=None)¶ Localized Intrinsic Valence Virtual Orbitals ref. 10.1021/acs.jctc.7b00493
Localized Intrinsic valence virtual orbitals are formed when the valence virtual orbitals are localized using an IBO-type of localization. Here the VVOs are created in the IAO basis then the IBO localization functions are called to localize the VVOs.
- Args:
mol : the molecule or cell object
orbocc : occupied molecular orbital coefficients
orbvirt : virtual molecular orbital coefficients
- Kwargs:
- locmethodstring
the localization method ‘PM’ for Pipek Mezey localization or ‘IBO’ for the IBO localization
- iaos2D array
the array of IAOs
- s2D array
the overlap array in the ao basis
- Returns:
LIVVOs in the basis defined in mol object.
-
pyscf.lo.vvo.
vvo
(mol, orbocc, orbvirt, iaos=None, s=None, verbose=None)¶ Valence Virtual Orbitals ref. 10.1021/acs.jctc.7b00493
Valence virtual orbitals can be formed from the singular value decomposition of the overlap between the canonical molecular orbitals and an accurate underlying atomic basis set. This implementation uses the intrinsic atomic orbital as this underlying set. VVOs can also be formed from the null space of the overlap of the canonical molecular orbitals and the underlying atomic basis sets (IAOs). This is not implemented here.
- Args:
mol : the molecule or cell object
orbocc : occupied molecular orbital coefficients
orbvirt : virtual molecular orbital coefficients
- Kwargs:
- iaos2D array
the array of IAOs
- s2D array
the overlap array in the ao basis
- Returns:
VVOs in the basis defined in mol object.