Data stored in TREXIO

We consider here basis functions centered on nuclei. Hence, it is possibile to define dummy atoms to place basis functions in arbitrary positions.

The atomic basis set is defined as a list of shells. Each shell \(s\) is centered on a center \(A\), possesses a given angular momentum \(l\) and a radial function \(R_s\). The radial function is a linear combination of \(N_{\text{prim}}\) primitive functions that can be of type Slater (\(p=1\)) or Gaussian (\(p=2\)), parameterized by exponents \(\gamma_{ks}\) and coefficients \(a_{ks}\): \[ R_s(\mathbf{r}) = \mathcal{N}_s \vert\mathbf{r}-\mathbf{R}_A\vert^{n_s} \sum_{k=1}^{N_{\text{prim}}} a_{ks}\, f_{ks}(\gamma_{ks},p)\, \exp \left( - \gamma_{ks} \vert \mathbf{r}-\mathbf{R}_A \vert ^p \right). \]

In the case of Gaussian functions, \(n_s\) is always zero.

Different codes normalize functions at different levels. Computing normalization factors requires the ability to compute overlap integrals, so the normalization factors should be written in the file to ensure that the file is self-contained and does not need the client program to have the ability to compute such integrals.

Some codes assume that the contraction coefficients are for a linear combination of normalized primitives. This implies that a normalization constant for the primitive \(ks\) needs to be computed and stored. If this normalization factor is not required, \(f_{ks}=1\).

Some codes assume that the basis function are normalized. This implies the computation of an extra normalization factor, \(\mathcal{N}_s\). If the the basis function is not considered normalized, \(\mathcal{N}_s=1\).

All the basis set parameters are stored in one-dimensional arrays.

A plane wave is defined as

\[ \chi_j(\mathbf{r}) = \exp \left( -i \mathbf{G}_j \cdot \mathbf{r} \right) \]

The basis set is defined as the array of $k$-points in the reciprocal space \(\mathbf{G}_j\), defined in the pbc group. The kinetic energy cutoff e_cut is the only input data relevant to plane waves.

For example, consider H₂ with the following basis set (in GAMESS format), where both the AOs and primitives are considered normalized:

HYDROGEN
S   5
1         3.387000E+01           6.068000E-03
2         5.095000E+00           4.530800E-02
3         1.159000E+00           2.028220E-01
4         3.258000E-01           5.039030E-01
5         1.027000E-01           3.834210E-01
S   1
1         3.258000E-01           1.000000E+00
S   1
1         1.027000E-01           1.000000E+00
P   1
1         1.407000E+00           1.000000E+00
P   1
1         3.880000E-01           1.000000E+00
D   1
1         1.057000E+00           1.000000E+00

In TREXIO representaion we have:

type  = "Gaussian"
prim_num   = 20
shell_num   = 12

# 6 shells per H atom
nucleus_index =
[ 0, 0, 0, 0, 0, 0,
  1, 1, 1, 1, 1, 1 ]

# 3 shells in S (l=0), 2 in P (l=1), 1 in D (l=2)
shell_ang_mom =
[ 0, 0, 0, 1, 1, 2,
  0, 0, 0, 1, 1, 2 ]

# no need to renormalize shells
shell_factor =
[ 1., 1., 1., 1., 1., 1.,
  1., 1., 1., 1., 1., 1. ]

# 5 primitives for the first S shell and then 1 primitive per remaining shells in each H atom
shell_index =
[ 0, 0, 0, 0, 0, 1, 2, 3, 4, 5,
  6, 6, 6, 6, 6, 7, 8, 9, 10, 11 ]

# parameters of the primitives (10 per H atom)
exponent =
[ 33.87, 5.095, 1.159, 0.3258, 0.1027, 0.3258, 0.1027, 1.407, 0.388, 1.057,
  33.87, 5.095, 1.159, 0.3258, 0.1027, 0.3258, 0.1027, 1.407, 0.388, 1.057 ]

coefficient =
[ 0.006068, 0.045308, 0.202822, 0.503903, 0.383421, 1.0, 1.0, 1.0, 1.0, 1.0,
  0.006068, 0.045308, 0.202822, 0.503903, 0.383421, 1.0, 1.0, 1.0, 1.0, 1.0 ]

prim_factor =
[ 1.0006253235944540e+01, 2.4169531573445120e+00, 7.9610924849766440e-01
  3.0734305383061117e-01, 1.2929684417481876e-01, 3.0734305383061117e-01,
  1.2929684417481876e-01, 2.1842769845268308e+00, 4.3649547399719840e-01,
  1.8135965626177861e+00, 1.0006253235944540e+01, 2.4169531573445120e+00,
  7.9610924849766440e-01, 3.0734305383061117e-01, 1.2929684417481876e-01,
  3.0734305383061117e-01, 1.2929684417481876e-01, 2.1842769845268308e+00,
  4.3649547399719840e-01, 1.8135965626177861e+00 ]

For example, consider H₂ molecule with the following effective core potential (in GAMESS input format for the H atom):

H-ccECP GEN 0 1
3
1.00000000000000    1 21.24359508259891
21.24359508259891   3 21.24359508259891
-10.85192405303825  2 21.77696655044365
1
0.00000000000000    2 1.000000000000000

In TREXIO representation this would be:

num = 8

# lmax+1 per atom
max_ang_mom_plus_1 = [ 1, 1 ]

# number of core electrons to remove per atom
zcore = [ 0, 0 ]

# first 4 ECP elements correspond to the first H atom ; the remaining 4 elements are for the second H atom
nucleus_index = [
  0, 0, 0, 0,
  1, 1, 1, 1
  ]

# 3 first ECP elements correspond to potential of the P orbital (l=1), then 1 element for the S orbital (l=0) ; similar for the second H atom
ang_mom = [
  1, 1, 1, 0,
  1, 1, 1, 0
  ]

# ECP quantities that can be attributed to atoms and/or angular momenta based on the aforementioned ecp_nucleus and ecp_ang_mom arrays
coefficient = [
  1.00000000000000, 21.24359508259891, -10.85192405303825, 0.00000000000000,
  1.00000000000000, 21.24359508259891, -10.85192405303825, 0.00000000000000
  ]

exponent = [
  21.24359508259891, 21.24359508259891, 21.77696655044365, 1.000000000000000,
  21.24359508259891, 21.24359508259891, 21.77696655044365, 1.000000000000000
  ]

power = [
  -1, 1, 0, 0,
  -1, 1, 0, 0
  ]

• \[ \hat{V}_{\text{ne}} = \sum_{A=1}^{N_\text{nucl}} \sum_{i=1}^{N_\text{elec}} \frac{-Z_A }{\vert \mathbf{R}_A - \mathbf{r}_i \vert} \] : electron-nucleus attractive potential,
• \[ \hat{T}_{\text{e}} = \sum_{i=1}^{N_\text{elec}} -\frac{1}{2}\hat{\Delta}_i \] : electronic kinetic energy
• \(\hat{h} = \hat{T}_\text{e} + \hat{V}_\text{ne} + \hat{V}_\text{ECP}\) : core electronic Hamiltonian

The one-electron integrals for a one-electron operator \(\hat{O}\) are \[ \langle p \vert \hat{O} \vert q \rangle \], returned as a matrix over atomic orbitals.

The two-electron integrals for a two-electron operator \(\hat{O}\) are \[ \langle p q \vert \hat{O} \vert r s \rangle \] in physicists notation, where \(p,q,r,s\) are indices over atomic orbitals.

• \[ \hat{W}_{\text{ee}} = \sum_{i=2}^{N_\text{elec}} \sum_{j=1}^{i-1} \frac{1}{\vert \mathbf{r}_i - \mathbf{r}_j \vert} \] : electron-electron repulsive potential operator.
• \[ \hat{W}^{lr}_{\text{ee}} = \sum_{i=2}^{N_\text{elec}} \sum_{j=1}^{i-1} \frac{\text{erf}(\mu\, \vert \mathbf{r}_i - \mathbf{r}_j \vert)}{\vert \mathbf{r}_i - \mathbf{r}_j \vert} \] : electron-electron long range potential

The Cholesky decomposition of the integrals can also be stored:

\[ \langle ij | kl \rangle = \sum_{\alpha} G_{ik\alpha} G_{jl\alpha} \]

The operators as the same as those defined in the AO one-electron integrals section. Here, the integrals are given in the basis of molecular orbitals.

The operators are the same as those defined in the AO two-electron integrals section. Here, the integrals are given in the basis of molecular orbitals.

The first form of Jastrow factor is the one used in the CHAMP program:

\[ J(\mathbf{r},\mathbf{R}) = J_{\text{eN}}(\mathbf{r},\mathbf{R}) + J_{\text{ee}}(\mathbf{r}) + J_{\text{eeN}}(\mathbf{r},\mathbf{R}) \]

\(J_{\text{eN}}\) contains electron-nucleus terms:

\[ J_{\text{eN}}(\mathbf{r},\mathbf{R}) = \sum_{i=1}^{N_\text{elec}} \sum_{\alpha=1}^{N_\text{nucl}}\left[ \frac{a_{1,\alpha}\, f_\alpha(R_{i\alpha})}{1+a_{2,\alpha}\, f_\alpha(R_{i\alpha})} + \sum_{p=2}^{N_\text{ord}^a} a_{p+1,\alpha}\, [f_\alpha(R_{i\alpha})]^p - J_{\text{eN}}^\infty \right] \]

\(J_{\text{ee}}\) contains electron-electron terms:

\[ J_{\text{ee}}(\mathbf{r}) = \sum_{i=1}^{N_\text{elec}} \sum_{j=1}^{i-1} \left[ \frac{\frac{1}{2}\big(1 + \delta^{\uparrow\downarrow}_{ij}\big)\,b_1\, f_{\text{ee}}(r_{ij})}{1+b_2\, f_{\text{ee}}(r_{ij})} + \sum_{p=2}^{N_\text{ord}^b} b_{p+1}\, [f_{\text{ee}}(r_{ij})]^p - J_{\text{ee},ij}^\infty \right] \]

\(\delta^{\uparrow\downarrow}_{ij}\) is zero when the electrons \(i\) and \(j\) have the same spin, and one otherwise.

\(J_{\text{eeN}}\) contains electron-electron-Nucleus terms:

\[ J_{\text{eeN}}(\mathbf{r},\mathbf{R}) = \sum_{\alpha=1}^{N_{\text{nucl}}} \sum_{i=1}^{N_{\text{elec}}} \sum_{j=1}^{i-1} \sum_{p=2}^{N_{\text{ord}}} \sum_{k=0}^{p-1} \sum_{l=0}^{p-k-2\delta_{k,0}} c_{lkp\alpha} \left[ g_{\text{ee}}({r}_{ij}) \right]^k \nonumber \\ \left[ \left[ g_\alpha({R}_{i\alpha}) \right]^l + \left[ g_\alpha({R}_{j\alpha}) \right]^l \right] \left[ g_\alpha({R}_{i\,\alpha}) \, g_\alpha({R}_{j\alpha}) \right]^{(p-k-l)/2} \] \(c_{lkp\alpha}\) are non-zero only when \(p-k-l\) is even.

The terms \(J_{\text{ee},ij}^\infty\) and \(J_{\text{eN}}^\infty\) are shifts to ensure that \(J_{\text{eN}}\) and \(J_{\text{ee}}\) have an asymptotic value of zero:

\[ J_{\text{eN}}^{\infty} = \frac{a_{1,\alpha}\, \kappa_\alpha^{-1}}{1+a_{2,\alpha}\, \kappa_\alpha^{-1}} + \sum_{p=2}^{N_\text{ord}^a} a_{p+1,\alpha}\, \kappa_\alpha^{-p} \] \[ J_{\text{ee},ij}^{\infty} = \frac{\frac{1}{2}\big(1 + \delta^{\uparrow\downarrow}_{ij}\big)\,b_1\, \kappa_{\text{ee}}^{-1}}{1+b_2\, \kappa_{\text{ee}}^{-1}} + \sum_{p=2}^{N_\text{ord}^b} b_{p+1}\, \kappa_{\text{ee}}^{-p} \]

\(f\) and \(g\) are scaling function defined as

\[ f_\alpha(r) = \frac{1-e^{-\kappa_\alpha\, r}}{\kappa_\alpha} \text{ and } g_\alpha(r) = e^{-\kappa_\alpha\, r}, \]

Mu-Jastrow is based on a one-parameter correlation factor that has been introduced in the context of transcorrelated methods. This correlation factor imposes the electron-electron cusp, and it is built such that the leading order in \(1/r_{12}\) of the effective two-electron potential reproduces the long-range interaction of the range-separated density functional theory. Its analytical expression reads

\[ J(\mathbf{r}, \mathbf{R}) = J_{\text{eeN}}(\mathbf{r}, \mathbf{R}) + J_{\text{eN}}(\mathbf{r}, \mathbf{R}) \].

The electron-electron cusp is incorporated in the three-body term

\[ J_\text{eeN} (\mathbf{r}, \mathbf{R}) = \sum_{i=1}^{N_\text{elec}} \sum_{j=1}^{i-1} \, u\left(\mu, r_{ij}\right) \, \Pi_{\alpha=1}^{N_{\text{nucl}}} \, E_\alpha({R}_{i\alpha}) \, E_\alpha({R}_{j\alpha}), \]

where ww\(u\) is an electron-electron function

\[ u\left(\mu, r\right) = \frac{r}{2} \, \left[ 1 - \text{erf}(\mu\, r) \right] - \frac{1}{2 \, \mu \, \sqrt{\pi}} \exp \left[ -(\mu \, r)^2 \right]. \]

This electron-electron term is tuned by the parameter \(\mu\) which controls the depth and the range of the Coulomb hole between electrons.

An envelope function has been introduced to cancel out the Jastrow effects between two-electrons when at least one is close to a nucleus (to perform a frozen-core calculation). The envelope function is given by

\[ E_\alpha(R) = 1 - \exp\left( - \gamma_{\alpha} \, R^2 \right). \]

In particular, if the parameters \(\gamma_\alpha\) tend to zero, the Mu-Jastrow factor becomes a two-body Jastrow factor:

\[ J_{\text{ee}}(\mathbf{r}) = \sum_{i=1}^{N_\text{elec}} \sum_{j=1}^{i-1} \, u\left(\mu, r_{ij}\right) \]

and for large \(\gamma_\alpha\) it becomes zero.

To increase the flexibility of the Jastrow and improve the electron density the following electron-nucleus term is added

\[ J_{\text{eN}}(\mathbf{r},\mathbf{R}) = \sum_{i=1}^{N_\text{elec}} \sum_{\alpha=1}^{N_\text{nucl}} \, \left[ \exp\left( a_{\alpha} R_{i \alpha}^2 \right) - 1\right]. \]

The parameter \(\mu\) is stored in the ee array, the parameters \(\gamma_\alpha\) are stored in the een array, and the parameters \(a_\alpha\) are stored in the en array.