Orbitals and basis sets#

Electronic wave functions and molecular orbitals#

State vectors representing single-electron systems reside in a Hilbert space formed as a direct product

\[ \mathcal{V} = \mathcal{V}^\mathrm{o} \otimes \mathcal{V}^\mathrm{c} \otimes \mathcal{V}^\mathrm{s} \]

The underlying vector spaces are associated with orbital, charge, and spin degrees of freedom. The elements of \(\mathcal{V}\) can be expressed in terms of

\[\begin{split} \psi^{\mathrm{L}\alpha}(\mathbf{r},t) \otimes \begin{pmatrix} 1\\0 \end{pmatrix} \otimes \begin{pmatrix} 1\\0 \end{pmatrix}; \quad \psi^{\mathrm{L}\beta}(\mathbf{r},t) \otimes \begin{pmatrix} 1\\0 \end{pmatrix} \otimes \begin{pmatrix} 0\\1 \end{pmatrix} \end{split}\]
\[\begin{split} \psi^{\mathrm{S}\alpha}(\mathbf{r},t) \otimes \begin{pmatrix} 0\\1 \end{pmatrix} \otimes \begin{pmatrix} 1\\0 \end{pmatrix} ; \quad \psi^{\mathrm{S}\beta}(\mathbf{r},t) \otimes \begin{pmatrix} 0\\1 \end{pmatrix} \otimes \begin{pmatrix} 0\\1 \end{pmatrix} \end{split}\]

with a sum that equals a general wave function for a single-electron system, also known as a spinor,

\[\begin{split} \psi(\mathbf{r},t) = \begin{pmatrix} \psi^{\mathrm{L}\alpha} \\ \psi^{\mathrm{L}\beta} \\ \psi^{\mathrm{S}\alpha} \\ \psi^{\mathrm{S}\beta} \end{pmatrix} \end{split}\]

In this general case, the particle density becomes

\[ n(\mathbf{r},t) = \psi^\dagger \psi = \sum_\mathrm{c}^{\mathrm{L,S}} \sum_\sigma^{\alpha,\beta} \left| \psi^{\mathrm{c}\sigma}(\mathbf{r},t) \right|^2 \]

Spinor components \(\psi^{\mathrm{S}\alpha}\) and \(\psi^{\mathrm{S}\beta}\) are small and sometimes ignored, reducing the electronic wave function to two-component form

\[\begin{split} \psi(\mathbf{r},t) = \begin{pmatrix} \psi^{\mathrm{L}\alpha} \\ \psi^{\mathrm{L}\beta} \end{pmatrix} \end{split}\]

For a system described by a spin-independent Hamiltonian, the spatial and spin degrees of freedom are separable and we arrive at (suppressing the time variable) single-electron wave functions, or spin orbitals, of the form

\[\begin{split} \psi(\mathbf{r}) = \phi(\mathbf{r}) \sigma ; \quad \sigma = c_\alpha \begin{pmatrix} 1 \\ 0 \end{pmatrix} + c_\beta \begin{pmatrix} 0 \\ 1 \end{pmatrix} \end{split}\]

In the collinear approximation, all spin orbitals are eigenfunctions of the spin operator \(\hat{s}_z\) such that every spatial function (or molecular orbital), \(\phi(\mathrm{r})\), gives rise to two spin orbitals

\[\begin{split} \psi_{a}(\mathbf{r}) = \phi_a(\mathbf{r}) \begin{pmatrix} 1 \\ 0 \end{pmatrix} ; \quad \psi_{\bar{a}}(\mathbf{r}) = \phi_a(\mathbf{r}) \begin{pmatrix} 0 \\ 1 \end{pmatrix} \end{split}\]

referred to as \(\alpha\)- and \(\beta\)-spin orbitals, respectively.

Linear combination of atomic orbitals#

In quantum chemistry, molecular orbitals (MOs) are normally expanded in a set of atom-centered basis functions, or localized atomic orbitals (AOs)

\[ \phi_a(\mathbf{r}) = \sum_\alpha c_{\alpha a} \chi_\alpha(\mathbf{r} - \mathbf{R}_\alpha) \]

where \(\mathbf{R}_\alpha\) denotes the atomic position center of basis function \(\chi_\alpha\), and the expansion coefficients \(c_{\alpha a}\) are known as molecular orbital (MO) coefficients. Using the Dirac notation, this linear combination of atomic orbitals (LCAO) expansion takes the form

\[ | \phi_a \rangle = \sum_\alpha | \chi_\alpha \rangle c_{\alpha a} \]

It can be convenient to collect the set of MOs (and AOs) as a row vector and introduce a notation

\[ | \overline{\phi} \rangle = (\ldots, | \phi_a \rangle, \ldots) \]

The LCAO expansion for the entire set of MOs can then be compactly written as

\[ | \overline{\phi} \rangle = | \overline{\chi} \rangle C \]

where \(C\) is the MO coefficient matrix.

import veloxchem as vlx

mol_str = """
C        0.00000000    0.00000000    0.00000000
O        0.00000000    0.00000000    1.43
"""
molecule = vlx.Molecule.read_str(mol_str, units='angstrom')
basis = vlx.MolecularBasis.read(molecule, "sto-3g")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.compute(molecule, basis)

Ordering of atomic orbitals#

The VeloxChem program orders AO basis functions based on

  1. Angular momentum quantum number, \(l\)

  2. Angular momentum quantum number, \(m_l\)

  3. User-defined order of atoms

In our present example of carbon monoxide using a mimal basis set, the order of basis functions become

\(\chi_{1s}^{C}\), \(\chi_{2s}^{C}\), \(\chi_{1s}^{O}\), \(\chi_{2s}^{O}\) \(\chi_{2py}^{C}\), \(\chi_{2py}^{O}\), \(\chi_{2pz}^{C}\), \(\chi_{2pz}^{O}\), \(\chi_{2px}^{C}\), \(\chi_{2px}^{O}\)

import numpy as np
np.set_printoptions(precision = 3, suppress = True, linewidth = 170)

C = scf_drv.scf_tensors['C']
print('MO coefficient matrix:\n', C)
MO coefficient matrix:
 [[ 0.001  0.994  0.091  0.228  0.109  0.     0.     0.     0.     0.105]
 [-0.004  0.024 -0.26  -0.836 -0.47  -0.    -0.    -0.    -0.    -0.536]
 [ 0.995  0.     0.231 -0.114  0.036 -0.    -0.     0.    -0.    -0.089]
 [ 0.023 -0.005 -0.856  0.526 -0.196  0.     0.    -0.     0.     0.547]
 [ 0.     0.     0.    -0.     0.    -0.095 -0.42  -0.034  0.912 -0.   ]
 [-0.    -0.    -0.     0.     0.    -0.187 -0.826  0.02  -0.548  0.   ]
 [-0.004  0.004 -0.141  0.003  0.621  0.     0.     0.    -0.    -0.942]
 [-0.003  0.002  0.112  0.38  -0.578 -0.    -0.     0.    -0.    -0.852]
 [-0.     0.    -0.    -0.     0.    -0.42   0.095 -0.912 -0.034 -0.   ]
 [-0.    -0.    -0.    -0.     0.    -0.826  0.187  0.548  0.02   0.   ]]

Overlap matrix#

The overlap between two basis functions is

\[ S_{\alpha\beta} = \langle \chi_\alpha | \chi_\beta \rangle \]

and the overlap matrix can be written as

\[ S = \langle \overline{\chi} | \overline{\chi} \rangle \]

where

\[ \langle \overline{\chi} | = (\ldots, \langle \chi_\alpha |, \ldots)^T \]
S = scf_drv.scf_tensors['S']
print('Overlap matrix:\n', S)
Overlap matrix:
 [[ 1.     0.248  0.     0.017  0.     0.     0.    -0.031  0.     0.   ]
 [ 0.248  1.     0.022  0.255  0.     0.     0.    -0.251  0.     0.   ]
 [ 0.     0.022  1.     0.237  0.     0.     0.036  0.     0.     0.   ]
 [ 0.017  0.255  0.237  1.     0.     0.     0.337  0.     0.     0.   ]
 [ 0.     0.     0.     0.     1.     0.132  0.     0.     0.     0.   ]
 [ 0.     0.     0.     0.     0.132  1.     0.     0.     0.     0.   ]
 [ 0.     0.     0.036  0.337  0.     0.     1.    -0.285  0.     0.   ]
 [-0.031 -0.251  0.     0.     0.     0.    -0.285  1.     0.     0.   ]
 [ 0.     0.     0.     0.     0.     0.     0.     0.     1.     0.132]
 [ 0.     0.     0.     0.     0.     0.     0.     0.     0.132  1.   ]]

Basis sets#

Our software comes distibuted with several of the standard basis sets and families of basis sets. For a given element, the available basis sets in the VeloxChem library can be provided from the MolecularBasis class object. Should a basis set not be included in the software distribution, it can be readily supplied by the user in a separate file in the working directory.

print('Available basis sets for carbon:\n', vlx.MolecularBasis.get_avail_basis('C'))
Available basis sets for carbon:
 ['6-31++G', '6-31+G', '6-31+G**', '6-311++G', '6-311++G(2D,2P)', '6-311++G**', '6-311+G', '6-311+G(2D,P)', '6-311G', '6-311G(2DF,2PD)', '6-31G', '6-31G(2DF,P)', 'AUG-CC-PCVDZ', 'AUG-CC-PCVQZ', 'AUG-CC-PCVTZ', 'AUG-CC-PVDZ', 'AUG-CC-PVTZ', 'CC-PCVDZ', 'CC-PVDZ', 'CC-PVDZ-RI', 'CC-PVTZ', 'CC-PVTZ-RI', 'D-AUG-CC-PVDZ', 'D-AUG-CC-PVQZ', 'D-AUG-CC-PVTZ', 'DEF2-QZVP', 'DEF2-QZVPD', 'DEF2-QZVPP', 'DEF2-QZVPPD', 'DEF2-RI-J', 'DEF2-SV(P)', 'DEF2-SVP', 'DEF2-SVPD', 'DEF2-TZVP', 'DEF2-TZVPP', 'DEF2-TZVPPD', 'MIN-CC-PVDZ', 'SADLEJ-PVTZ', 'STO-3G', 'STO-6G']

Basis set assignment#

A basis set can be assigned to a molecule with the read method.

basis = vlx.MolecularBasis.read(molecule, 'cc-pVDZ')
print(basis.get_string('cc-pVDZ', molecule))
Molecular Basis (cc-pVDZ)
===========================

Basis: CC-PVDZ                                        

  Atom Contracted GTOs          Primitive GTOs           

  O   (3S,2P,1D)               (17S,4P,1D)              
  C   (3S,2P,1D)               (17S,4P,1D)              

Contracted Basis Functions : 28                       
Primitive Basis Functions  : 68                       

The moleular basis set label is available with use of the get_label method.

vlx.MolecularBasis.get_label(basis)
'CC-PVDZ'

Number of basis functions#

The number of contracted and primitive basis functions are available from the get_dimensions_of_basis and get_dimensions_of_primitive_basis methods, respectively.

nbas = vlx.MolecularBasis.get_dimensions_of_basis(basis, molecule)
nprim = vlx.MolecularBasis.get_dimensions_of_primitive_basis(basis, molecule)

print('Number of contracted basis functions:', nbas)
print('Number of primitive basis functions:', nprim)
Number of contracted basis functions: 28
Number of primitive basis functions: 68