Short Range Action

The exact pair density matrix for any interaction potential in infinite volume can be calculated by the matrix squaring technique by Storer (1968). This method is applied to the short-range part generated by the Ewald break-up of the Coulomb potential. Since the potential is short range the periodicity is irrelevant. First, one factorizes the density matrix into a center-of-mass term and a term depending on the relative coordinates. The latter term is equivalent to the density matrix for a particle with the reduced mass $\mu^{-1} = m_1^{-1} + m_2^{-1}$ in an external potential. The one expands the pair density matrix in partial waves. In $D=3$ dimensions, it reads,

$$\rho({\bf r},{\bf r}';\tau) = \frac{1}{4 \pi r r'} \; \sum_... ... l + 1) \; \rho_l(r,r';\tau) \; {\rm P}_l( \cos \theta) \quad,$$ (41)

where $\theta$ is the angle between ${\bf r}$ and ${\bf r}'$ and P$_l$ denotes the $l$th Legendre polynomial. For spherically symmetric interaction potentials, different partial wave components $\rho_l$ do not mix and can be derived from independent matrix squaring procedures. The six dimensional pair density matrix is reduced to a sum of two dimensional objects. Each component $\rho_l$ satisfies a 1 dimensional Bloch equation with an additional term,

$$- \frac{\partial \rho_l(r,r';\beta)}{\partial \beta} = \lef... ... \frac{ \lambda}{r^2} \: l(l+1) \: \right] \rho_l(r,r';\beta)$$ (42)

and also fulfills the convolution equation,

$$\rho_l(r,r';\beta) = \int_0^\infty \!\!\!\! dr'' \; \rho_l(r,r'';\beta/2) \; \rho_l(r'',r';\beta/2) \quad.$$ (43)

This is a one dimensional integral for a given pair of $r$ and $r'$, which can be calculated numerically. In order to derive the pair density matrix for a time step $\tau^{-1}=10^6\,\rm {K}$, one typically performs of the order of $m=12$ matrix squarings starting at the inverse temperature $\tau_0^{-1} = 2^m\times10^6\,\rm {K}$. The partial waves are initialized using semi-classical expression analogous to Eq. 2.29, for details see Ceperley (1995) and Magro (1996). The resulting pair density matrix can be verified by using the Feynman-Kac formula Eq. 2.27 in a separate MC simulation (Pollock and Ceperley, 1984).

The pair density matrix is between two particles at initial position $({\bf r}_i,{\bf r}_j)$ and final position $({\bf r}'_i,{\bf r}'_j)$ needs to be evaluated very frequently in a PIMC simulation. Using the fact that initial and final position cannot be too far apart, one can expand the action in a power series. It is convenient to use the three distance

$$q = \frac{1}{2}( \vert{\bf r}\vert+\vert{\bf r}'\vert)$$ (44)

$$s = \vert{\bf r}-{\bf r}'\vert$$ (45)

$$z = \vert{\bf r}\vert - \vert{\bf r}'\vert,$$ (46)

where ${\bf r}={\bf r}_i-{\bf r}_j$ and ${\bf r}'={\bf r}'_i-{\bf r}'_j$. The variables $s$ and $z$ are of the order of $\sqrt{\lambda \tau}$. The action can then be expanded as,

$$u({\bf r},{\bf r}';\tau) = \frac{1}{2}\left[ u_0(r;\tau) + u... ...} \sum_{j=0}^{k} u_{kj}(q;\tau) \, z^{2j} \, s^{2(k-j)} \quad.$$ (47)

where $n_{\rm A}$ denotes the order of the expansion. In zeroth order, only the first term, the end-point action, is considered. The following terms are off-diagonal contributions, which are important because they allow to reduce the number of time slices in a PIMC simulation. The same expansion formula is used for the contributions to the energy given by $\beta$ derivative of the action. $n_{\rm E}$ denotes the order in this expansion.