Backflow

Backflow displacement is the sum of homogeneous, isotropic electron–electron terms \(\eta(r_{ij})\), isotropic electron–nucleus terms \(\mu(r_{iI})\) centered on the nuclei, isotropic electron–electron–nucleus terms \(\Phi(r_{iI}, r_{jI}, r_{ij})\), \(\Theta(r_{iI}, r_{jI}, r_{ij})\), also centered on the nuclei:

\[\mathbf{\xi} = \sum_{i \neq j}^{N_e} \eta(r_{ij})\mathbf{r}_{ij} + \sum_{I=1}^{N_I} \mu(r_{iI})\mathbf{r}_{iI} + \sum_{I=1}^{N_I}\sum_{i \neq j}^{N_e} \Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij} + \Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI}\]

When all-electron atoms are present, the electron-nucleus Kato cusp conditions cannot be fulfilled unless the backflow displacement at the nuclear position is zero. This can be obtained by applying smooth cutoffs around such atoms. An artificial multiplicative cutoff function \(g(r_{iI})\) is applied to all contributions to the backflow displacement of particle \(i\) that do not depend on the distance \(r_{iI}\) to the all-electron atom \(I\). This includes the homogeneous backflow term displacement like \(\eta(r_{ij})\mathbf{r}_{ij}\) and the inhomogeneous contributions centered on other atom \(J \neq I\) like \(\Phi(r_{iJ}, r_{jJ}, r_{ij})\mathbf{r}_{ij}\) or \(\Theta(r_{iJ}, r_{jJ}, r_{ij})\mathbf{r}_{ij}\). The \(g(r_{iI})\) function used is of the form:

\[g(r_{iI}) = \left(\frac{r_{iI}}{L_{gI}}\right)^2 \left[6 - 8 \left(\frac{r_{iI}}{L_{gI}}\right) + 3 \left(\frac{r_{iI}}{L_{gI}}\right)^2 \right]\]

Backflow part of wavefunction is represented by the casino.Backflow class.

It must be initialized from the configuration files:

from casino.readers import CasinoConfig
from casino.backflow import Backflow

config_path = <path to a directory containing input file>
config = CasinoConfig(config_path)
config.read()
backflow = Backflow(config)

To prevent code duplication, we need to prepare the necessary intermediate data:

import numpy as np

neu, ned = config.input.neu, config.input.ned
ne = neu + ned
r_e = np.random.uniform(-1, 1, ne * 3).reshape((ne, 3))
atom_positions = config.wfn.atom_positions
e_vectors = np.expand_dims(r_e, 1) - np.expand_dims(r_e, 0)
n_vectors = np.expand_dims(r_e, 0) - np.expand_dims(atom_positions, 1)
e_powers = backflow.ee_powers(e_vectors)
n_powers = backflow.en_powers(n_vectors)
r_ij = np.linalg.norm(e_vectors, axis=2)
r_iI = np.linalg.norm(n_vectors, axis=2)
ee_spin_mask = np.ones(shape=(ne, ne), dtype=int)
ee_spin_mask[:neu, :neu] = 0
ee_spin_mask[neu:, neu:] = 2
en_spin_mask = np.zeros(shape=(ne,), dtype=int)
en_spin_mask[neu:] = 1
C = backflow.trunc

Summary of Methods

Backflow class has a following methods:

Method	Output	Shape
eta_term	\(\eta(r_{ij})\mathbf{r}_{ij}\)	\((3N_e,)\)
mu_term	\(\mu(r_{ij})\mathbf{r}_{ij}\)	\((3N_e,)\)
phi_term	\(\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij} + \Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI}\)	\((3N_e,)\)
eta_term_gradient	\(\nabla \eta(r_{ij})\mathbf{r}_{ij}\)	\((3N_e, 3N_e)\)
mu_term_gradient	\(\nabla \mu(r_{ij})\mathbf{r}_{ij}\)	\((3N_e, 3N_e)\)
phi_term_gradient	\(\nabla (\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij} + \Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI})\)	\((3N_e, 3N_e)\)
eta_term_laplacian	\(\Delta \eta(r_{ij})\mathbf{r}_{ij}\)	\((3N_e,)\)
mu_term_laplacian	\(\Delta \mu(r_{ij})\mathbf{r}_{ij}\)	\((3N_e,)\)
phi_term_laplacian	\(\Delta (\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij} + \Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI})\)	\((3N_e,)\)

eta-term

\(\eta(r_{ij})\mathbf{r}_{ij}\) consists of a complete power expansion in electron-electron distances \(r_{ij}\):

\[\eta(r_{ij}) = (1 - r_{ij}/L_\eta)^C\Theta(L_\eta - r_{ij}) \sum_{k=0}^{N_\eta}c_kr^k_{ij}\]

where \(\Theta\) is the Heaviside function. Electron-electron Kato cusp conditions at \(r_{ij} = 0\) satisfied by constraint for spin-like electrons only:

\[L_\eta c_1 = C c_0\]

For certain electron coordinates, \(\eta\) term can be obtained with casino.Backflow.eta_term() method:

backflow.eta_term(e_powers, e_vectors)[1]

this is equivalent to (continues from):

from numpy.polynomial.polynomial import polyval
poly = polyval(r_ij, backflow.eta_parameters.T)
cutoff = np.maximum(1 - r_ij / backflow.eta_cutoff[0], 0) ** C
np.fill_diagonal(cutoff, 0)
eta = np.expand_dims(cutoff * np.choose(ee_spin_mask, poly, mode='wrap'), -1)
np.sum(-e_vectors * eta, axis=0)

mu-term

\(\mu(r_{iI})\mathbf{r}_{iI}\) term consists of a complete power expansion in electron-nucleus distances \(r_{iI}\):

\[\mu(r_{iI}) = (1 - r_{iI}/L_\mu)^C\Theta(L_\mu - r_{iI}) \sum_{k=0}^{N_\mu}d_kr^k_{iI}\]

where \(\Theta\) is the Heaviside function. The electron-nucleus Kato cusp conditions at \(r_{iI} = 0\) satisfied if

\[L_{\mu I} d_{1 I} = C d_{0 I}\]

for all atoms, and in addition,

\[d_{0 I} = 0\]

for all-electron atoms.

For certain electron coordinates, \(\mu\) term can be obtained with casino.Backflow.mu_term() method:

backflow.mu_term(n_powers, n_vectors)[1]

this is equivalent to (continues from):

from numpy.polynomial.polynomial import polyval
poly = polyval(r_iI, backflow.mu_parameters[0].T)
cutoff = np.maximum(1 - r_iI / backflow.mu_cutoff, 0) ** C
n_vectors * np.expand_dims(cutoff[0] * np.choose(en_spin_mask, poly, mode='wrap'), -1)

phi-term

\(\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij}\) and \(\Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI}\) terms describe two-electron backflow displacements in terms of \(r_{ij}\) , \(r_{iI}\) , and \(r_{jI}\) and vectors \(\mathbf{r}_{ij}\) , \(\mathbf{r}_{iI}\):

\[\Phi(r_{iI}, r_{jI}, r_{ij}) = \sum_{k=0}^{N_{\Phi I}^{eN}}\sum_{l=0}^{N_{\Phi I}^{eN}}\sum_{m=0}^{N_{\Phi I}^{ee}} \Phi_{klmI}\]

\[\Theta(r_{iI}, r_{jI}, r_{ij}) = \sum_{k=0}^{N_{\Phi I}^{eN}}\sum_{l=0}^{N_{\Phi I}^{eN}}\sum_{m=0}^{N_{\Phi I}^{ee}} \Theta_{klmI}\]

where

\[\Phi_{klmI} = (1 - r_{iI}/L_{\Phi I})^C(1 - r_{jI}/L_{\Phi I})^C\Theta(L_{\Phi I} - r_{iI})\Theta(L_{\Phi I} - r_{jI})\phi_{klmI}r_{iI}^kr_{jI}^lr_{ij}^m\]

\[\Theta_{klmI} = (1 - r_{iI}/L_{\Phi I})^C(1 - r_{jI}/L_{\Phi I})^C\Theta(L_{\Phi I} - r_{iI})\Theta(L_{\Phi I} - r_{jI})\theta_{klmI}r_{iI}^kr_{jI}^lr_{ij}^m\]

and \(\Theta\) is the Heaviside function. To ensure electron–electron Kato cusp conditions folowing \(3(N_{\Phi I}^{ee} + N_{\Phi I}^{en} + 1)\) constraints is applied:

\[\sum_{l,m}^{l+m=\alpha}(C\phi_{0lmI} - L_{\phi I}\phi_{1lmI}) = \sum_{k,m}^{k+m=\alpha}(C\phi_{k0mI} - L_{\phi I}\phi_{k1mI}) = \sum_{k,m}^{k+m=\alpha}(C\theta_{k0mI} - L_{\phi I}\theta_{k1mI}) = 0\]

another \(2N_{\Phi I}^{en} + 1\) constraints from the electron-electron Kato cusp conditions:

\[\sum_{k,l}^{k+l=\alpha}\theta_{kl1I} = 0\]

and extra \(2N_{\Phi I}^{en} + 1\) constraints for spin-like electrons:

\[\sum_{k,l}^{k+l=\alpha}\phi_{kl1I} = 0\]

for all-electron atoms there are \(4(N_{\Phi I}^{ee} + N_{\Phi I}^{en})+2\) constraints on \(\phi_{klm}\)

\[\sum_{l,m}^{l+m=\alpha}\phi_{0lmI} = \sum_{l,m}^{l+m=\alpha}m\phi_{0lmI} = \sum_{k,m}^{k+m=\alpha}\phi_{k0mI} = \sum_{k,m}^{k+m=\alpha}m\phi_{k0mI} = 0\]

for all-electron atoms there are \(3(N_{\Phi I}^{ee} + N_{\Phi I}^{en})+2\) constraints on \(\theta_{klm}\)

\[\sum_{l,m}^{l+m=\alpha}\theta_{0lmI} = \sum_{l,m}^{l+m=\alpha}m\theta_{0lmI} = \sum_{k,m}^{k+m=\alpha}m\theta_{k0mI} = 0\]

For certain electron coordinates, \(\Phi\) and \(\Theta\) terms can be obtained with casino.Backflow.phi_term() method:

backflow.phi_term(e_powers, n_powers, e_vectors, n_vectors)[1]

this is equivalent to (continues from):

from numpy.polynomial.polynomial import polyval3d
r_ijI = np.tile(r_iI[0], (ne, 1))
cutoff = np.maximum(1 - r_iI/backflow.phi_cutoff, 0) ** C
phi_poly = polyval3d(r_ijI, r_ijI.T, r_ij, backflow.phi_parameters[0].T)
theta_poly = polyval3d(r_ijI, r_ijI.T, r_ij, backflow.theta_parameters[0].T)
phi = np.outer(cutoff[0], cutoff[0]) * np.choose(ee_spin_mask, phi_poly, mode='wrap')
theta = np.outer(cutoff[0], cutoff[0]) * np.choose(ee_spin_mask, theta_poly, mode='wrap')
np.fill_diagonal(theta, 0)
np.sum(-e_vectors * np.expand_dims(phi, -1) + n_vectors * np.expand_dims(theta, -1), axis=0)

eta-term gradient

Considering that vector gradient of spherically symmetric vector function (in 3-D space) is:

\[\nabla (f(r)\mathbf{r}) = f'(r) \mathbf{\hat r} \otimes \mathbf{r} + f \cdot \mathbf{I}\]

There is only two non-zero terms of \(\eta(r_{ij})\) gradient, i.e. by \(i\)-th or \(j\)-th electron coordinates:

\[\nabla_{e_i} (\eta(r_{ij})\mathbf{r}_{ij}) = (1 - r_{ij}/L_\eta)^C\Theta(L_\eta - r_{ij}) \sum_{k=0}^{N_\eta} \left[\left(\frac{k}{r_{ij}} - \frac{C}{L_\eta - r_{ij}}\right) \mathbf{\hat r}_{ij} \otimes \mathbf{r}_{ij} + \mathbf{I} \right] c_kr^k_{ij}\]

\[\nabla_{e_j} (\eta(r_{ij})\mathbf{r}_{ij}) = - \nabla_{e_i} (\eta(r_{ij})\mathbf{r}_{ij})\]

where \(\mathbf{\hat r}_{ij}\) is the unit vector in the direction of the \(\mathbf{r}_{ij}\)

For certain electron coordinates, \(\eta\) term gradient can be obtained with casino.Backflow.eta_term_gradient() method:

backflow.eta_term_gradient(e_powers, e_vectors)[1]

this is equivalent to (continues from):

from numpy.polynomial.polynomial import polyval
L = backflow.eta_cutoff
k = np.arange(backflow.eta_parameters.shape[1])
cutoff = np.maximum(1 - r_ij / backflow.eta_cutoff[0], 0) ** C
np.fill_diagonal(cutoff, 0)
poly = polyval(r_ij, backflow.eta_parameters.T)
poly_k = polyval(r_ij, (k * backflow.eta_parameters).T)
unit_e_vectors = np.nan_to_num(e_vectors/np.expand_dims(r_ij, -1))
t1 = cutoff * np.choose(ee_spin_mask, poly, mode='wrap')
t2 = cutoff * np.choose(ee_spin_mask, poly * C / (r_ij - L) + poly_k / r_ij, mode='wrap')
tt1 = np.einsum('ij,ab->aibj', np.eye(3), np.diag(np.sum(t1, axis=0)) - t1)
np.fill_diagonal(t2, 0)
tt2 = np.einsum('abi,abj,ab->aibj', e_vectors, unit_e_vectors, -t2)
(tt1 + tt2).reshape(ne*3, ne*3)

mu-term gradient

Considering that vector gradient of spherically symmetric vector function (in 3-D space) is:

\[\nabla (f(r)\mathbf{r}) = f'(r) \mathbf{\hat r} \otimes \mathbf{r} + f \cdot \mathbf{I}\]

There is only one non-zero term of \(\mu(r_{iI})\) gradient, i.e. by \(i\)-th electron coordinates:

\[\nabla_{e_i} (\mu(r_{iI})\mathbf{r}_{iI}) = (1 - r_{iI}/L_\mu)^C\Theta(L_\mu - r_{iI}) \sum_{k=0}^{N_\mu} \left[\left(\frac{k}{r_{iI}} - \frac{C}{L_\mu - r_{iI}}\right) \mathbf{\hat r}_{iI} \otimes \mathbf{r}_{iI} + \mathbf{I}\right] d_kr^k_{ij}\]

where \(\mathbf{\hat r}_{iI}\) is the unit vector in the direction of the \(\mathbf{r}_{iI}\)

For certain electron coordinates, \(\mu\) term gradient can be obtained with casino.Backflow.mu_term_gradient() method:

backflow.mu_term_gradient(n_powers, n_vectors)[1]

this is equivalent to (continues from):

from numpy.polynomial.polynomial import polyval
L = backflow.mu_cutoff
k = np.arange(backflow.mu_parameters[0].shape[1])
cutoff = np.maximum(1 - r_iI / L, 0) ** C
poly = polyval(r_iI, backflow.mu_parameters[0].T)
poly_k = polyval(r_iI, (k * backflow.mu_parameters[0]).T)
unit_n_vectors = n_vectors/np.expand_dims(r_iI, -1)
t1 = cutoff[0] * np.choose(en_spin_mask, poly, mode='wrap')
t2 = cutoff[0] * np.choose(en_spin_mask, poly * C / (r_iI - L) + poly_k / r_iI, mode='wrap')
tt1 = np.einsum('ij,ab,Ia->aibj', np.eye(3), np.eye(ne), t1)
tt2 = np.einsum('Iai,Iaj,ab,Ia->aibj', n_vectors, unit_n_vectors, np.eye(ne), t2)
(tt1 + tt2).reshape(ne*3, ne*3)

phi-term gradient

Considering that vector gradient of spherically symmetric vector function (in 3-D space) is:

\[\nabla (f(r)\mathbf{r}) = f'(r) \mathbf{\hat r} \otimes \mathbf{r} + f \cdot \mathbf{I}\]

There is only two non-zero terms of \(\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij}\) gradient, i.e. by \(i\)-th:

\[\nabla_{e_i} (\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij}) = \sum_{k=0}^{N_{\Phi I}^{eN}} \sum_{l=0}^{N_{\Phi I}^{eN}} \sum_{m=0}^{N_{\Phi I}^{ee}} \left[\left(\frac{k}{r_{iI}} - \frac{C}{L_{\Phi I} - r_{iI}} \right) \mathbf{\hat r}_{iI} \otimes \mathbf{r}_{ij} + \left(\frac{m}{r_{ij}} \right) \mathbf{\hat r}_{ij} \otimes \mathbf{r}_{ij} + \mathbf{I} \right] \Phi_{klmI}\]

or \(j\)-th electron coordinates:

\[\nabla_{e_j} (\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij}) = \sum_{k=0}^{N_{\Phi I}^{eN}} \sum_{l=0}^{N_{\Phi I}^{eN}} \sum_{m=0}^{N_{\Phi I}^{ee}} \left[\left(\frac{l}{r_{jI}} - \frac{C}{L_{\Phi I} - r_{jI}} \right) \mathbf{\hat r}_{jI} \otimes \mathbf{r}_{ij} - \left(\frac{m}{r_{ij}} \right) \mathbf{\hat r}_{ij} \otimes \mathbf{r}_{ij} - \mathbf{I} \right] \Phi_{klmI}\]

There is only two non-zero terms of \(\Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI}\) gradient, i.e. by \(i\)-th:

\[\nabla_{e_i} (\Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI}) = \sum_{k=0}^{N_{\Phi I}^{eN}} \sum_{l=0}^{N_{\Phi I}^{eN}} \sum_{m=0}^{N_{\Phi I}^{ee}} \left[\left(\frac{k}{r_{iI}} -\frac{C}{L_{\Phi I} - r_{iI}}\right) \mathbf{\hat r}_{iI} \otimes \mathbf{r}_{iI} + \left(\frac{m}{r_{ij}}\right) \mathbf{\hat r}_{ij} \otimes \mathbf{r}_{iI} + \mathbf{I} \right] \Theta_{klmI}\]

or \(j\)-th electron coordinates:

\[\nabla_{e_j} (\Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI}) = \sum_{k=0}^{N_{\Phi I}^{eN}} \sum_{l=0}^{N_{\Phi I}^{eN}} \sum_{m=0}^{N_{\Phi I}^{ee}} \left[\left(\frac{l}{r_{jI}} - \frac{C}{L_{\Phi I} - r_{jI}}\right) \mathbf{\hat r}_{jI} \otimes \mathbf{r}_{iI} - \left(\frac{m}{r_{ij}}\right) \mathbf{\hat r}_{ij} \otimes \mathbf{r}_{iI} \right] \Theta_{klmI}\]

where

\[\Phi_{klmI} = (1 - r_{iI}/L_{\Phi I})^C(1 - r_{jI}/L_{\Phi I})^C\Theta(L_{\Phi I} - r_{iI})\Theta(L_{\Phi I} - r_{jI})\phi_{klmI}r_{iI}^kr_{jI}^lr_{ij}^m\]

\[\Theta_{klmI} = (1 - r_{iI}/L_{\Phi I})^C(1 - r_{jI}/L_{\Phi I})^C\Theta(L_{\Phi I} - r_{iI})\Theta(L_{\Phi I} - r_{jI})\theta_{klmI}r_{iI}^kr_{jI}^lr_{ij}^m\]

where \(\mathbf{\hat r}_{ij}\) is the unit vector in the direction of the \(\mathbf{r}_{ij}\) and \(\mathbf{\hat r}_{iI}\) is the unit vector in the direction of the \(\mathbf{r}_{iI}\)

For certain electron coordinates, \(\phi\) term gradient can be obtained with casino.Backflow.phi_term_gradient() method:

backflow.phi_term_gradient(e_powers, n_powers, e_vectors, n_vectors)[1]

this is equivalent to (continues from):

from numpy.polynomial.polynomial import polyval3d

eta-term laplacian

Considering that vector laplacian of spherically symmetric vector function (in 3-D space) is:

\[\Delta (f(r)\mathbf{r}) = \left(f''(r) + \frac{4}{r} f'(r)\right) \mathbf{r}\]

There is only two non-zero terms of \(\eta(r_{ij})\mathbf{r}_{ij}\) laplacian, i.e. by \(i\)-th or \(j\)-th electron coordinates:

\[\Delta_{e_i} (\eta(r_{ij})\mathbf{r}_{ij}) = (1 - r_{ij}/L_\eta)^C\Theta(L_\eta - r_{ij}) \mathbf{r}_{ij}\sum_{k=0}^{N_\eta} \left[\frac{C(C-1)}{(L_\eta - r_{ij})^2} - \frac{2C(k+2)}{r_{ij}(L_\eta - r_{ij})} + \frac{k(k+3)}{r_{ij}^2} \right] c_kr^k_{ij}\]

\[\Delta_{e_j} (\eta(r_{ij})\mathbf{r}_{ij}) = - \Delta_{e_i} (\eta(r_{ij})\mathbf{r}_{ij})\]

For certain electron coordinates, \(\eta\) laplacian term can be obtained with casino.Backflow.eta_term_laplacian() method:

backflow.eta_term_laplacian(e_powers, e_vectors)[1]

this is equivalent to (continues from):

from numpy.polynomial.polynomial import polyval

mu-term laplacian

Considering that vector laplacian of spherically symmetric vector function (in 3-D space) is:

\[\Delta (f(r)\mathbf{r}) = \left(f''(r) + \frac{4}{r} f'(r)\right) \mathbf{r}\]

There is only one non-zero term of \(\mu(r_{iI})\mathbf{r}_{iI}\) laplacian, i.e. by \(i\)-th electron coordinates:

\[\Delta_{e_i} (\mu(r_{iI})\mathbf{r}_{iI}) = (1 - r_{iI}/L_\mu)^C\Theta(L_\mu - r_{iI}) \mathbf{r}_{iI}\sum_{k=0}^{N_\mu} \left[\frac{C(C-1)}{(L_\mu - r_{iI})^2} - \frac{2C(k+2)}{r_{iI}(L_\mu - r_{iI})} + \frac{k(k+3)}{r_{iI}^2} \right]d_kr^k_{iI}\]

For certain electron coordinates, \(\mu\) term laplacian can be obtained with casino.Backflow.mu_term_laplacian() method:

backflow.mu_term_laplacian(n_powers, n_vectors)[1]

this is equivalent to (continues from):

from numpy.polynomial.polynomial import polyval

phi-term laplacian

Considering that gradient of spherically symmetric function (in 3-D space) is:

\[\nabla f(r) = f'(r) \mathbf{\hat r}\]

and laplacian of spherically symmetric vector function (in 3-D space) is:

\[\Delta f(r) = f''(r) + \frac{2}{r} f'(r)\]

and \(\Phi\) term addent is a product of \(f=r_{ij}^m\), \(g=(1−r_{iI}/L_{\Phi I})^С r_{iI}^k\), \(h=(1−r_{jI}/L_{\Phi I})^С r_{jI}^l\), \(\phi_{klmI}\) and \(\mathbf{r}_{ij}\) so using:

\[\frac{\Delta_{e_i} \Phi}{\Phi} = \frac{\Delta_{e_i} (fg\mathbf{r}_{ij})}{fg} = \left(\frac{\Delta_{e_i} f}{f} + 2\frac{\nabla_{e_i} f}{f} \cdot \frac{\nabla_{e_i} g}{g} + \frac{\Delta_{e_i} g}{g}\right)\mathbf{r}_{ij} + 2\left(\frac{\nabla_{e_i} f}{f} + \frac{\nabla_{e_i} g}{g}\right) \nabla_{e_i} \mathbf{r}_{ij}\]

\[\frac{\Delta_{e_j} \Phi}{\Phi} = \frac{\Delta_{e_j} (fh\mathbf{r}_{ij})}{fh} = \left(\frac{\Delta_{e_j} f}{f} + 2\frac{\nabla_{e_j} f}{f} \cdot \frac{\nabla_{e_j} h}{h} + \frac{\Delta_{e_j} h}{h}\right)\mathbf{r}_{ij} + 2\left(\frac{\nabla_{e_j} f}{f} + \frac{\nabla_{e_j} h}{h}\right) \nabla_{e_j} \mathbf{r}_{ij}\]

There is only two non-zero terms of \(\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij}\) laplacian, i.e. by \(i\)-th:

\[l_{kmI} = \frac{\Delta_{e_i} f}{f} + \frac{\Delta_{e_i} g}{g} = \left( \frac{m(m+1)}{r_{ij}^2} + \frac{k(k+1)}{r_{iI}^2} + \frac{C(C+1)}{(L_{\Phi I} - r_{iI})^2} - \frac{2C(k+1)}{r_{iI}(L_{\Phi I} - r_{iI})} \right)\]

\[d_{kmI} = \frac{\nabla_{e_i} f}{f} \cdot \frac{\nabla_{e_i} g}{g} = \frac{m}{r_{ij}} \left( \frac{k}{r_{iI}} - \frac{C}{L_{\Phi I} - r_{iI}} \right) (\mathbf{\hat r}_{ij} \cdot \mathbf{\hat r}_{iI})\]

\[g_{kmI} = \frac{\nabla_{e_i} f}{f} + \frac{\nabla_{e_i} g}{g} = \frac{m}{r_{ij}} \mathbf{\hat r}_{ij} + \left( \frac{k}{r_{iI}} - \frac{C}{L_{\Phi I} - r_{iI}} \right)\mathbf{\hat r}_{iI}\]

\[\Delta_{e_i} (\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij}) = \sum_{k=0}^{N_{fI}^{eN}}\sum_{l=0}^{N_{fI}^{eN}}\sum_{m=0}^{N_{fI}^{ee}} (l_{kmI}\mathbf{r}_{ij} + 2d_{kmI}\mathbf{r}_{ij} + 2g_{kmI}) \Phi_{klmI}\]

or \(j\)-th electron coordinates:

\[l_{lmI} = \frac{\Delta_{e_j} f}{f} + \frac{\Delta_{e_j} g}{g} = \left( \frac{m(m+1)}{r_{ij}^2} + \frac{l(l+1)}{r_{jI}^2} + \frac{C(C+1)}{(L_{\Phi I} - r_{jI})^2} - \frac{2C(l+1)}{r_{iI}(L_{\Phi I} - r_{jI})} \right)\]

\[d_{lmI} = \frac{\nabla_{e_j} f}{f} \cdot \frac{\nabla_{e_j} h}{h} = - \frac{m}{r_{ij}} \left( \frac{l}{r_{jI}} - \frac{C}{L_{\Phi I} - r_{jI}} \right) (\mathbf{\hat r}_{ij} \cdot \mathbf{\hat r}_{jI})\]

\[g_{lmI} = \frac{\nabla_{e_j} f}{f} + \frac{\nabla_{e_j} h}{h} = \frac{m}{r_{ij}} \mathbf{\hat r}_{ij} - \left( \frac{l}{r_{jI}} - \frac{C}{L_{\Phi I} - r_{jI}} \right)\mathbf{\hat r}_{jI}\]

\[\Delta_{e_j} (\Phi(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{ij}) = \sum_{k=0}^{N_{fI}^{eN}}\sum_{l=0}^{N_{fI}^{eN}}\sum_{m=0}^{N_{fI}^{ee}} (l_{lmI}\mathbf{r}_{ij} + 2d_{lmI}\mathbf{r}_{ij} + 2g_{lmI}) \Phi_{klmI}\]

\(\Theta\) term addent is a product of \(f=r_{ij}^m\), \(g=(1−r_{iI}/L_{\Phi I})^С r_{iI}^k\), \(h=(1−r_{jI}/L_{\Phi I})^С r_{jI}^l\), \(\theta_{klmI}\) and \(\mathbf{r}_{iI}\) so using:

\[\frac{\Delta_{e_i} \Theta}{\Theta} = \frac{\Delta_{e_i} (fg\mathbf{r}_{ij})}{fg} = \left(\frac{\Delta_{e_i} f}{f} + 2\frac{\nabla_{e_i} f}{f} \cdot \frac{\nabla_{e_i} g}{g} + \frac{\Delta_{e_i} g}{g}\right)\mathbf{r}_{iI} + 2\left(\frac{\nabla_{e_i} f}{f} + \frac{\nabla_{e_i} g}{g}\right) \nabla_{e_i} \mathbf{r}_{iI}\]

\[\frac{\Delta_{e_j} \Theta}{\Theta} = \frac{\Delta_{e_j} (fh\mathbf{r}_{ij})}{fh} = \left(\frac{\Delta_{e_j} f}{f} + 2\frac{\nabla_{e_j} f}{f} \cdot \frac{\nabla_{e_j} h}{h} + \frac{\Delta_{e_j} h}{h}\right)\mathbf{r}_{iI} + 2\left(\frac{\nabla_{e_j} f}{f} + \frac{\nabla_{e_j} h}{h}\right) \nabla_{e_j} \mathbf{r}_{iI}\]

There is only two non-zero terms of \(\Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI}\) laplacian, i.e. by \(i\)-th:

\[l_{kmI} = \left( \frac{m(m+1)}{r_{ij}^2} + \frac{k(k+1)}{r_{iI}^2} + \frac{C(C+1)}{(L_{\Phi I} - r_{iI})^2} - \frac{2C(k+1)}{r_{iI}(L_{\Phi I} - r_{jI})} \right)\]

\[d_{kmI} = \frac{m}{r_{ij}} \left( \frac{k}{r_{iI}} - \frac{C}{L_{\Phi I} - r_{iI}} \right) (\mathbf{\hat r}_{ij} \cdot \mathbf{\hat r}_{iI})\]

\[g_{kmI} = \frac{m}{r_{ij}} \mathbf{\hat r}_{ij} + \left( \frac{k}{r_{iI}} - \frac{C}{L_{\Phi I} - r_{iI}} \right)\mathbf{\hat r}_{iI}\]

\[\Delta_{e_i} (\Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI}) = \sum_{k=0}^{N_{fI}^{eN}}\sum_{l=0}^{N_{fI}^{eN}}\sum_{m=0}^{N_{fI}^{ee}} (l_{kmI}\mathbf{r}_{iI} + 2d_{kmI}\mathbf{r}_{iI} + 2g_{kmI}) \Theta_{klmI}\]

or \(j\)-th electron coordinates:

\[l_{lmI} = \left( \frac{m(m+1)}{r_{ij}^2} + \frac{l(l+1)}{r_{jI}^2} + \frac{C(C+1)}{(L_{\Phi I} - r_{jI})^2} - \frac{2C(k+1)}{r_{iI}(L_{\Phi I} - r_{jI})} \right)\]

\[d_{lmI} = - \frac{m}{r_{ij}} \left( \frac{l}{r_{jI}} - \frac{C}{L_{\Phi I} - r_{jI}} \right) (\mathbf{\hat r}_{ij} \cdot \mathbf{\hat r}_{jI})\]

\[\Delta_{e_j} (\Theta(r_{iI}, r_{jI}, r_{ij})\mathbf{r}_{iI}) = \sum_{k=0}^{N_{fI}^{eN}}\sum_{l=0}^{N_{fI}^{eN}}\sum_{m=0}^{N_{fI}^{ee}} (l_{lmI}\mathbf{r}_{iI} + 2d_{lmI}\mathbf{r}_{iI}) \Theta_{klmI}\]

where

\[\Phi_{klmI} = (1 - r_{iI}/L_{\Phi I})^C(1 - r_{jI}/L_{\Phi I})^C\Theta(L_{\Phi I} - r_{iI})\Theta(L_{\Phi I} - r_{jI})\phi_{klmI}r_{iI}^kr_{jI}^lr_{ij}^m\]

\[\Theta_{klmI} = (1 - r_{iI}/L_{\Phi I})^C(1 - r_{jI}/L_{\Phi I})^C\Theta(L_{\Phi I} - r_{iI})\Theta(L_{\Phi I} - r_{jI})\theta_{klmI}r_{iI}^kr_{jI}^lr_{ij}^m\]

For certain electron coordinates, \(\phi\) term laplacian can be obtained with casino.Backflow.phi_term_laplacian() method:

backflow.phi_term_laplacian(e_powers, n_powers, e_vectors, n_vectors)

this is equivalent to (continues from):

from numpy.polynomial.polynomial import polyval3d