23. lo — Orbital localization

23.1. Foster-Boys, Edmiston-Ruedenberg, Pipek-Mezey localization

Foster-Boys localization

Edmiston-Ruedenberg localization

Pipek-Mezey localization

pyscf.lo.pipek.atomic_pops(mol, mo_coeff, method='meta_lowdin')[source]
Kwargs:
method
: string
one of mulliken, lowdin, meta_lowdin
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.

23.2. Meta-Lowdin

23.3. Natural atomic orbitals

Natural atomic orbitals Ref:

  1. Weinhold et al., J. Chem. Phys. 83(1985), 735-746
pyscf.lo.nao.set_atom_conf(element, description)[source]

Change the default atomic core and valence configuration to the one given by “description”. See data/elements.py for the default configuration.

Args:
element
: str or int
Element symbol or nuclear charge
description
: str 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

23.4. Intrinsic Atomic Orbitals

Intrinsic Atomic Orbitals ref. JCTC, 9, 4834

pyscf.lo.iao.fast_iao_mullikan_pop(mol, dm, iaos, verbose=5)[source]
Args:

mol : the molecule or cell object

iaos
: 2D array
(orthogonal or non-orthogonal) IAO orbitals
Returns:
mullikan population analysis in the basis IAO
pyscf.lo.iao.iao(mol, orbocc, minao='minao')[source]

Intrinsic Atomic Orbitals. [Ref. JCTC, 9, 4834]

Args:

mol : the molecule or cell object

orbocc
: 2D 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')[source]

Create a molecule which uses reference minimal basis