2. Quick Start

PYSCF-ITA is a pyscf extension for interpreting the outputs of molecular quantum chemistry calculations with Information-Theoretic Approach (ITA).

To use PYSCF-ITA, you should first perform an electronic structure calculation with pyscf. Here we give a quick introduction of basic pyscf calculation. For more information, please check the offical website of is Pyscf.

This quick start section provides a brief introduction to the use of PySCF-ITA in common quantum chemical simulations. These make reference to specific examples within the dedicated examples directory. For brevity, and so as to not repeat a number of function calls, please note that the cells below often share objects in-between one another. The selection below is far from exhaustive: additional details on the various modules are presented in the accompanying user guide and within the examples directory.

2.1. Molecular build

The fist step of pyscf calculation is to build the gto.Mole object, it include attributes of molecular geometry, basis, integrals and so on. Besides the following required parameters, there are many other optional parameters that you can visit pyscf offical website for more detail.

Required parameters:

• mol = gto.Mole() Initialize an instance of gto.Mole class to variable mol.

• mol.atom = ... Assign the geometry of molecule.

• mol.basis = "6-31G" Define the basis set of molecule.

• mol.charge = 0 Charge of molecule, by default 0.

• mol.multiplicity = 1 Multiplicity of molecule, by default 1.

• mol.unit = 'A' Length unit of molecular geometry; either ‘A’, ‘Angstrom’ for Angstrom; ‘B’, ‘Bohr’ for Bohr; by default ‘Angstrom’.

• mol.build() Call the method to build the molecule.

from pyscf import gto

mol = gto.Mole()
mol.atom = """
8  0.0  0.0  0.0
O  0.0  0.0  1.5
H  1.0  0.0  0.0
H  0.0  0.7  1.0
"""
mol.basis = "6-31G"
mol.charge = 0
mol.multiplicity = 1
mol.unit = 'A'
mol.build()

<pyscf.gto.mole.Mole at 0x1109e86d0>

2.2. Pyscf Calculation

2.2.1. Hartree Fock

The Hartree–Fock method is an uncorrelated mean-field theory that offers a qualitative description of chemical systems. Although Hartree–Fock theory is only qualitatively correct, it forms the basis for more accurate models and becomes the cornerstone of ab initio quantum chemistry.

Minimal examples for HF calculations are shown below.

from pyscf import scf

hf_mf = scf.HF(mol) # likewise for RHF, UHF, ROHF and GHF
hf_mf.run()

converged SCF energy = -150.585033780838

RHF object of <class 'pyscf.scf.hf.RHF'>

2.2.2. Kohn-Sham Density Functional Theory

Running a KS-DFT calculation is as straightforward as the above HF counterpart. In Libxc you can find various types of functionals: LDA, GGA, and meta-GGA (mGGA) functionals. LDAs, GGAs, and meta-GGAs depend on local information, in the sense that the value of the density functional part of the energy density at a given point depends only on the values of the density, the gradient of the density, and the kinetic energy density and/or the density laplacian, respectively, at the given point.

Minimal examples for each type of calculations are shown below.

from pyscf import dft

dft_mf = dft.KS(mol) # likewise for RKS, UKS, ROKS and GKS
dft_mf.xc = 'm062x'
dft_mf.run()

converged SCF energy = -151.305840677392

RKS object of <class 'pyscf.dft.rks.RKS'>

2.2.3. Single reference Post-HF Method

In computational chemistry, post–Hartree–Fock (post-HF) methods are the set of methods developed to improve on the Hartree–Fock (HF), or self-consistent field (SCF) method. They add electron correlation which is a more accurate way of including the repulsions between electrons than in the Hartree–Fock method where repulsions are only averaged.

Minimal examples for each type of calculations are shown below

from pyscf import mp, ci, cc

mymp = mp.MP2(hf_mf).run()
myci = ci.CISD(hf_mf).run()
mycc = cc.CCSD(hf_mf).run()

E(MP2) = -150.854045557244  E_corr = -0.269011776406098
E(SCS-MP2) = -150.850347103892  E_corr = -0.26531332305376
E(RCISD) = -150.8389090064467  E_corr = -0.2538752256083305
E(CCSD) = -150.8591226999828  E_corr = -0.2740889191444854

2.2.4. Multi reference Post-HF Method

Multiconfigurational self-consistent field (MCSCF) methods go beyond the single-determinant Hartree-Fock (HF) method by allowing the wave function to become a linear combination of multiple determinants. While the configurations i.e. determinants can in principle be chosen in an arbitrary number of ways, PySCF focuses on the complete active space (CAS) family of methods, where the set of electron configurations is defined in terms of a given set of active orbitals, also known as the “active space”. The CAS method generates all possible electron configurations that can be formed from the set of the active orbitals, and is therefore equivalent to an FCI procedure on a subset of the molecular orbitals; please see Configuration interaction (CISD and FCI) for a discussion on the FCI method. The use of MCSCF methods is crucial for reliable modeling of systems that exhibit nearly degenerate orbitals i.e. static correlation, such as transition metal complexes. For a general discussion of MCSCF methods, we direct the reader to References and, and for specific details about PySCF’s implementation of CASSCF.

Minimal examples for each type of calculations are shown below

from pyscf import mcscf

ncas, nelecas = (6,8)
mycasscf = mcscf.CASSCF(hf_mf, ncas, nelecas).run()
mycasci = mcscf.CASCI(hf_mf, ncas, nelecas).run()

CASSCF energy = -150.679332235502
CASCI E = -150.679332235502  E(CI) = -21.6884571730814  S^2 = 0.0000000
CASCI E = -150.626755311336  E(CI) = -18.5296721415734  S^2 = 0.0000000

2.3. Grids

For the integral of DFT, since the analytical form is not as easy to derive, numerical integration is used. Numerical integration is the approximate computation of an integral using numerical techniques. The numerical computation of an integral is sometimes called quadrature. to the desired precision The basic problem in numerical integration is to compute an approximate solution to a definite integral of a function \(f(q)\) as a weighted sum of function values at specified points within the domain of integration:

\[\int_a^b f(q) dq = \sum_i^n f(q_i)\omega(i) \tag{1}\]

These grids are constructed starting from a “parent” grid \((N_r, N_{\Omega})\) consisting of Nr radial spheres with \(N_{\Omega}\) angular (Lebedev) grid points on each, then systematically pruning the number of angular points in regions where sophisticated angular quadrature is not necessary, such as near the nuclei where the charge density is nearly spherically symmetric and at long distance from the nuclei where it varies slowly. A large number of points are retained in the valence region where angular accuracy is critical.

Table 1. Deafult radial and angular fineness of atomic grids.

Elements	Radial Part	Angular Part
H, He	50	302
Li - Ne	75	302
Na - Ar	80	434
K - Kr	90	434
Rb - Xe	95	434
Cs - Rn	100	434

Required parameters:

• grids = dft.Grids(mol) Initialize an instance of dft.Grids class to variable grids.

• grids.atom_grid = (75, 302) Radial and angular fineness of atomic grids; by default using the fineness in Table 2 .

• grids.becke_schme Atom to molecule Reweight scheme; by default original_becke scheme.

• grids.prune = True If appy grids prune scheme; by default True.

• grids.build() Call the method to build the grids.

from pyscf import dft

grids = dft.Grids(mol)
grids.atom_grid = (75, 302)

grids.becke_scheme
grids.prune = True
grids.build()

(90600,)

2.4. Compute ITA

The information-theoretic quantities developed by Shubin Liu et al. is becoming more and more popular in predicting and understanding many chemical relevant problems, such as reactivity, regioselectivity, aromaticity, pKa and so on. See Acta Phys. -Chim. Sin., 32, 98 (2016) and WIREs Comput Mol Sci., e1461 (2019) for reviews.

Table 2. ITA Code table.

Name	Symbol	Dependence	Code	Minimal Rung
Shannon entropy	\(S_S\)	\(n(\mathbf{r})\)	11	[1,1]
Fisher information	\(I_F\)	\(n(\mathbf{r}), \mathbf{\nabla} n(\mathbf{r})\)	12	[2,1]
Alternative Fisher information	\(I^{\prime}_F\)	\(n(\mathbf{r}), \mathbf{\nabla} n(\mathbf{r}),\nabla^2 n(\mathbf{r})\)	13	[3,1]
Renyi entropy	\(R_n\)	\(n(\mathbf{r})\)	14	[1,1]
Tsallis entropy	\(T_n\)	\(n(\mathbf{r})\)	15	[1,1]
Onicescu information	\(E_n\)	\(n(\mathbf{r})\)	16	[1,1]
Ghosh–Berkowitz–Parr entropy	\(S_\text{GBP}\)	\(n(\mathbf{r}), \mathbf{\nabla} n(\mathbf{r}),\nabla^2 n(\mathbf{r}), \tau(\mathbf{r})\)	17	[3,1]
Relative Shannon entropy	\(S^r_S\)	\(n(\mathbf{r})\), \(n_0(\mathbf{r})\)	21	[1,1]
Relative Fisher Information	\(I^r_F\)	\(n(\mathbf{r})\), \(\nabla n(\mathbf{r})\), \(n_0(\mathbf{r})\), \(\nabla n_0(\mathbf{r})\)	22	[2,2]
Relative Alternative Fisher Information	\(I^{\prime r}_F\)	\(n(\mathbf{r})\), \(\nabla^2 n(\mathbf{r})\), \(n_0(\mathbf{r})\)	23	[3,1]
Relative Renyi entropy	\(R^r_n\)	\(n(\mathbf{r})\), \(n_0(\mathbf{r})\)	24	[1,1]
Relative Tsallis entropy	\(T^r_n\)	\(n(\mathbf{r})\), \(n_0(\mathbf{r})\)	25	[1,1]
Relative Onicescu information	\(E_n^r\)	\(n(\mathbf{r})\), \(n_0(\mathbf{r})\)	26	[1,1]
\(G_1\)	\(G_1\)	\(n_A(\mathbf{r})\), \(\nabla^2_A n(\mathbf{r})\), \(n^0_A(\mathbf{r})\)	28	[3,1]
\(G_2\)	\(G_2\)	\(n_A(\mathbf{r})\), \(\nabla^2_A n(\mathbf{r})\), \(n^0_A(\mathbf{r})\), \(\nabla^2 n^0_A(\mathbf{r})\)	29	[3,3]
\(G_3\)	\(G_3\)	\(n_A(\mathbf{r})\), \(\nabla_A n(\mathbf{r})\), \(n^0_A(\mathbf{r})\), \(\nabla n^0_A(\mathbf{r})\)	30	[2,2]

Notes

• In Ghosh–Berkowitz–Parr entropy equation, k is the Boltzmann constant, by default 1.0 for convenience.

• In Renyi entropy equation, \(\log_{10}\) is used instead of \(\ln\) for convenience.

• The G1, G2 and G3 ITA works only for atoms-in-molecules representation.

Required parameters:

• ita = ITA() Initialize an instance of ITA class to variable ita.

• ita.method = ... Assign the pyscf method to ITA instance.

• ita.grids = grids Assign the pyscf grids to ITA instance.

• ita.spin = 'ab' Type of density to compute; either ‘ab’ for total, ‘a’ for alpha-spin, ‘b’ for beta-spin, or ‘m’ for magnetization; by default ‘ab’.

• ita.promolecule = promol Assign instance of ProMolecule to ITA instance, by default None.

• ita.build() Call the method to build the ITA class.

Minimal examples for each type of calculations are shown below.

from pyscf.ita.ita import ITA

ita = ITA(hf_mf, grids)
ita.build()

ita.shannon_entropy()
ita.fisher_information()
ita.alternative_fisher_information()

917.4670485618817

ita.renyi_entropy(n=2)
ita.tsallis_entropy(n=2)
ita.onicescu_information(n=2)
ita.GBP_entropy()

117.73791583391437

Required parameters:

• promol = ProMolecule(hf_mf) Initialize an instance of ProMolecule class to variable promol.

• promol.pro_charge = {'H':0, 'O':0} Charge of proatoms in promolecule; if None, set all the charge to 0.

• promol.pro_multiplicity =  {'H':2, 'O':3} Multiplicity of protaoms in promolecule; if None, set all the multiplicity to corresponding 0 charge atom.

• promol.build() Call the method to build the ProMolecule class.

Notes

• If the attributes of ita instance is updated, one must call the build() method again.

from pyscf.ita.promolecule import ProMolecule

promol = ProMolecule(hf_mf)
promol.pro_charge = {'H':0, 'O':0}
promol.pro_multiplicity = {'H':2, 'O':3}
promol.build()

# Call the build method again to update ita instance.
ita.promolecule = promol
ita.build()

ita.relative_shannon_entropy()
ita.relative_fisher_information()
ita.relative_alternative_fisher_information()
ita.relative_renyi_entropy(n=2)
ita.relative_tsallis_entropy(n=2)
ita.relative_onicescu_information(n=2)

18.60935839551751

Required parameters:

• ita_code = ... List of ITA code to calculate. Check the ITA code table above for detail.

• representation = 'electron density' Representation; either ‘electron density’, ‘shape function’ or ‘atoms in molecules’; by default ‘electron density’.

• partition = 'Hirshfeld' Atom in molecules partition scheme; either ‘Becke’ or ‘Hirshfeld’; by default ‘Hirshfeld’.

• filename = './ita.log' File path and name of output, by default ‘ita.log’.

Notes

• If you wanna use Becke partition in following calcuation, before build the instance of ITA class, one must update the grids attributes coords and weights by: grids.coords, grids.weights = grids.gen_partition(mol)

• TODO: Implement the Bader’s quantum theory of atoms in molecules (QTAIM) partition.

ita.batch_compute(ita_code=[11,12,13,14,15,16,17], representation='electron density', partition = 'hirshfeld', filename='./ita_ed.log')

ita_ed.log

ita.batch_compute(ita_code=[11,12,13,14,15,16,17], representation='shape function', partition = 'hirshfeld', filename='./ita_sf.log')

ita_sf.log

ita.batch_compute(ita_code=[21,22,23,24,25,26,28,29,30], representation='atoms in molecules', partition = 'hirshfeld', filename='./ita_aim.log')

ita_aim.log

Notes

• The G1, G2 and G3 ITA works only for atoms-in-molecules representation.