1. Theory

Hybrid particle–field simulations switch out the ordinary particle–particle Lennard-Jones interactions with interactions between particles and a slowly varying density field. In this way, the most expensive part of normal molecular dynamics simulations is circumvented. Hybrid particle–field densities are defined as

\[\phi(\mathbf{r}) = \sum_{i=1}^NP(\mathbf{r}-\mathbf{r}_i),\]

where \(\mathbf{r}\) is a spatial coordinate, \(P\) is a window function used to distribute particle number densities onto a computational grid, and \(\mathbf{r}_i\) is the position of particle \(i\) (of total \(N\)). By default, HyMD uses a Cloud-In-Cell (CIC) window function. A Hamiltonian form for a system of \(N\) particles in \(M\) molecules is

\[\mathcal{H}(\{\mathbf{r}\})=\sum_{m=1}^MH_0(\{\mathbf{r},\mathbf{v}\}_m)+W[\phi(\mathbf{r})]\]

with \(H_0\) being a standard intramolecular Hamiltonian form (see Intramolecular bonds) including kinetic terms, while \(W\) is a density dependent interaction energy functional.

In the Hamiltonian hPF-MD formalism [Bore and Cascella, 2020], the density field is filtered using a grid-independent filtering function \(H\),

\[\tilde\phi(\mathbf{r})=\int\mathrm{d}\mathbf{x}\,\phi(\mathbf{x})H(\mathbf{r}-\mathbf{x}).\]

The filter smooths the density, ensuring that \(\tilde\phi\) and \(W[\tilde\phi([\phi])]\) both converge as the grid size is reduced.

1.1. External potential

The external potential acting on a particle is defined as the functional derivative of \(W\) with respect to \(\phi\). In the filtered formalism, the potential takes the form

\[\begin{split}V(\mathbf{r}) &= \int\mathrm{d}\mathbf{y}\,\frac{\delta w}{\delta\phi(\mathbf{r})} \\ &= \int\mathrm{d}\mathbf{y}\,\frac{\delta w}{\delta\tilde\phi(\mathbf{y})}\frac{\delta\tilde\phi(\mathbf{y})}{\delta\phi(\mathbf{r})},\end{split}\]

under the assumption of a local form of the interaction energy functional, \(W[\tilde\phi]=\int\mathrm{d}\mathbf{r}\,w[\tilde\phi(\mathbf{r})]\). Note that

\[\frac{\delta \tilde\phi(\mathbf{y})}{\delta \phi(\mathbf{r})} = H(\mathbf{y}-\mathbf{r}).\]

1.2. Force interpolation

The forces on particle \(i\) are obtained by differentiation of the external potential,

\[\mathbf{F}_i=-\int\mathrm{d}\mathbf{r}\,\nabla V(\mathbf{r})P(\mathbf{r}-\mathbf{r}_i).\]

1.3. Reciprocal space calculations

The field operations in HyMD are discretised and performed on a grid in reciprocal space using (discrete) fast Fourier transform algorithms. After interpolating the density \(\phi_{ijk}\) with CIC, we apply the filtering and obtain the discrete version of the external potential by

\[\tilde\phi_{ijk}=\mathrm{FFT}^{-1}\big[\mathrm{FFT}(\phi)\mathrm{FFT}(H)\big]\]

and

\[V_{ijk}=\mathrm{FFT}^{-1}\left[\mathrm{FFT}\left(\frac{\delta w(\tilde\phi)}{\delta \tilde\phi}\right)\mathrm{FFT}(H)\right].\]

The forces are obtained by differentiation of \(V\) in Fourier space as

\[\nabla V_{ijk} = \mathrm{FFT}^{-1}\left[i\mathbf{k}\mathrm{FFT}\left(\frac{\delta w(\tilde\phi)}{\delta \tilde\phi}\right)\mathrm{FFT}(H)\right].\]

1.4. Filter

By default, the filter used in HyMD is a simple Gaussian of the form

\[\begin{split}H(x) &= \frac{1}{\sqrt{2\pi}\sigma}\exp\left[\frac{-x^2}{2\sigma^2}\right] \\ \hat{H}(k) &= \exp\left[\frac{-\sigma^2k^2}{2}\right].\end{split}\]

For more details about the filtering, see Filtering.

1.5. Hamiltonian form

The default form of the interaction energy functional in HyMD is

\[W=\frac{1}{2\phi_0}\int\mathrm{d}\mathbf{r}\sum_{\text{i},\text{j}}\tilde\chi_{\text{i}-\text{j}}\tilde\phi_\text{i}(\mathbf{r})\tilde\phi_\text{j}(\mathbf{r}) + \frac{1}{2\kappa}\left(\sum_\text{k}\tilde\phi_\text{k}(\mathbf{r})-\phi_0\right)^2.\]

In the case of constant pressure simulations (NPT), the interaction energy functional becomes

\[W=\frac{1}{2\rho_0}\int\mathrm{d}\mathbf{r}\sum_{\text{i},\text{j}}\tilde\chi_{\text{i}-\text{j}}\tilde\phi_\text{i}(\mathbf{r})\tilde\phi_\text{j}(\mathbf{r}) + \frac{1}{2\kappa}\left(\sum_\text{k}\tilde\phi_\text{k}(\mathbf{r})-a\right)^2.\]

where, \(\rho_0= 1 / ν\) is an intrinsic parameter corresponding to the specific volume \((ν)\) of a coarse-grained particle, \(a\) is a calibrated parameter to obtain the correct average density at the target temperature and pressure. See Interaction energy functionals for details.