SDF

class hoomd.hpmc.compute.SDF(xmax, dx)

Bases: Compute

Compute the scale distribution function via volume perturbations.

Parameters:

• xmax (float) – Maximum x value at the right hand side of the rightmost bin $[\mathrm{length}]$.

• dx (float) – Bin width $[\mathrm{length}]$.

SDF computes the probability distributions $s_{\mathrm{comp}}(x)$ and $s_{\mathrm{exp}}(x)$ of particles overlapping as a function of separation for compressive and expansive perturbations, respectively. It estimates $s_{\mathrm{comp}}(x)$ and $s_{\mathrm{exp}}(x)$ numerically by computing histograms with $\lfloor x_\mathrm{max}/ \delta x \rfloor$ bins of width dx ($\delta x$).

See also

Anderson 2016 describes the theory relating SDF to the system pressure and its implementation in HOOMD-blue. Eppenga and Frenkel 1984 present a derivation relating the scale distribution function to the system pressure for hard, convex particles. Allen 2006 describes the theory for calculating the pressure in systems with discontinuous potential energy functions. The expansive perturbations are based on theory described in de Miguel and Jackson.

Implementation

SDF constructs two histograms, one for compressive volume perturbations and one for expansive volume perturbations. These histograms represent $s_{\mathrm{comp}}(x)$ and $s_{\mathrm{exp}}(x)$. The following discussion applies to compressive volume perturbations and the computation of $s_{\mathrm{comp}}(x)$; the computation of $s_{\mathrm{exp}}(x)$ proceeds similarly as noted throughout the description.

For each pair of particles $i$ and $j$ SDF scales the particle separation vector by the factor $(1 \pm x)$ ($+$ for expansive perturbations, $-$ for compressive perturbations) and finds the smallest positive value of $x$ leading to either an overlap of the particle shapes (a “hard overlap”) or a discontinuous change in the pair energy $U_{\mathrm{pair},ij}$ (a “soft overlap”). For compressive perturbations:

$$\begin{split} x_{ij}(\vec{A}) = \min \{ & x \in \mathbb{R}_{> 0} : \\ & \mathrm{overlap}\left( S_i(\mathbf{q}_i), S_j(\mathbf{q}_j, (1-x)(\vec{r}_j - (\vec{r}_i + \vec{A}))) \right) \ne \emptyset \\ &\lor \\ & U_{\mathrm{pair},ij}((1-x)(\vec{r}_j - (\vec{r}_i + \vec{A})), \mathbf{q}_i, \mathbf{q}_j) \ne U_{\mathrm{pair},ij}(\vec{r}_j - (\vec{r}_i + \vec{A}), \mathbf{q}_i, \mathbf{q}_j) \\ \} & \end{split}$$

where $\mathrm{overlap}$ is the shape overlap function defined in hoomd.hpmc.integrate, $S_i$ is the shape of particle $i$, and $\vec{A} = h\vec{a}_1 + k\vec{a}_2 + l\vec{a}_3$ is a vector that translates by periodic box images. For expansive perturbations,

$$\begin{split} x_{ij}(\vec{A}) = \max \{ & x \in \mathbb{R}_{< 0} : \\ & \mathrm{overlap}\left( S_i(\mathbf{q}_i), S_j(\mathbf{q}_j, (1+x)(\vec{r}_j - (\vec{r}_i + \vec{A}))) \right) \ne \emptyset \\ &\lor \\ & U_{\mathrm{pair},ij}((1+x)(\vec{r}_j - (\vec{r}_i + \vec{A})), \mathbf{q}_i, \mathbf{q}_j) \ne U_{\mathrm{pair},ij}(\vec{r}_j - (\vec{r}_i + \vec{A}), \mathbf{q}_i, \mathbf{q}_j) \\ \} & \end{split}$$

For particle $i$, SDF finds the the minimum (maximum for expansive perturbations) value $x_i$. For compressive perturbations:

$$x_i = \min \{ x_{ij} : \vec{A} \in B_\mathrm{images}, j \in [0,N_\mathrm{particles}) \}$$

where the set of box images includes all image vectors necessary to find overlaps between particles in the primary image with particles in periodic images. For expansive perturbations:

$$x_i = \max \{ x_{ij} : \vec{A} \in B_\mathrm{images}, j \in [0,N_\mathrm{particles}) \}$$

SDF adds a single count to each histogram for each particle $i$, weighted by a factor that is a function of the change in energy upon overlap. For compressive perturbations:

$$s_{\mathrm{comp}}(x + \delta x/2) = \frac{1}{N_\mathrm{particles} \cdot \delta x} \sum_{i=0}^{N_\mathrm{particles}-1} [x \le x_i < x + \delta x] \cdot (1 - \exp(-\beta \Delta U_{i}))$$

where $\Delta U_{i}$ is the change in energy associated with the first overlap involving particle $i$ ($\infty$ for hard overlaps), the square brackets denote the Iverson bracket, and $s_{\mathrm{comp}}(x + \delta x/2)$ is evaluated for $\{ x \in \mathbb{R}, 0 \le x < x_\mathrm{max}, x = k \cdot \delta x, k \in \mathbb{Z}^* \}$ for compressive perturbations. For expansive perturbations,

$$s_{\mathrm{exp}}(x - \delta x/2) = \frac{1}{N_\mathrm{particles} \cdot \delta x} \sum_{i=0}^{N_\mathrm{particles}-1} [x - \delta x \le x_i < x] \cdot (1 - \exp(-\beta \Delta U_{i}))$$

where $s_{\mathrm{exp}}(x - \delta x/2)$ is evaluated for $\{ x \in \mathbb{R}, -|x_\mathrm{max}| < x \le 0, x = (k - \lfloor x_\mathrm{max} / \delta x \rfloor + 1) \cdot \delta x, k \in \mathbb{Z}^* \}$.

Pressure

The pressure $P$ is related to the one-sided limits $s_{\mathrm{comp}}(0+)$ and $s_{\mathrm{exp}}(0-)$, computed by fitting and extrapolating $s_{\mathrm{comp}}$ and $s_{\mathrm{exp}}$ to $x = 0$.

$$\beta P = \rho \left( 1 + \frac{s_{\mathrm{comp}}(0+) - s_{\mathrm{exp}}(0-)}{2d} \right)$$

where $d$ is the dimensionality of the system, $\rho$ is the number density, and $\beta = \frac{1}{kT}$. This measurement of the pressure is inherently noisy due to the nature of the sampling. Average P over many timesteps to obtain accurate results.

Assuming particle diameters are ~1, these parameter values typically achieve good results:

• xmax = 0.02

• dx = 1e-4

In systems near densest packings, dx=1e-5 may be needed along with smaller xmax. Check that $\sum_k s_\mathrm{comp}(x_k) \cdot dx \approx 0.5$.

Important

SDF samples pair configurations at discrete separations. Therefore, the computed pressure is correct only for potentials with constant values and step discontinuities.

Note

SDF always runs on the CPU.

Mixed precision

SDF uses reduced precision floating point arithmetic when checking for particle overlaps in the local particle reference frame.

Box images

SDF does not apply the minimum image convention. It supports small boxes where particles may overlap with non-primary images of other particles, including self overlap.

---

Members inherited from AutotunedObject:

property kernel_parameters

Kernel parameters. Read more...

property is_tuning_complete

Check if kernel parameter tuning is complete. Read more...

tune_kernel_parameters()

Start tuning kernel parameters. Read more...

---

Members defined in SDF:

xmax

Maximum x value at the right hand side of the rightmost bin $[\mathrm{length}]$.

Type:

float

dx

Bin width $[\mathrm{length}]$.

Type:

float

property P

Pressure in NVT simulations $\left[ \mathrm{energy} \ \mathrm{length}^{-d} \right]$.

$$P = \frac{1}{\beta} \beta P$$

where $\beta P$ is given by betaP.

Attention

In MPI parallel execution, P is available on rank 0 only. P is None on ranks >= 1.

(Loggable: category=”scalar”)

Type:

float

property betaP

Beta times pressure in NVT simulations $\left[ \mathrm{length}^{-d} \right]$.

Uses a polynomial curve fit of degree 5 to estimate $s_\mathrm{comp}(0+)$ (and $s_\mathrm{exp}(0-)$ if required) and computes the pressure via:

$$\beta P = \rho \left(1 + \frac{s_\mathrm{comp}(0+)}{2d} + \frac{s_\mathrm{exp}(0-)}{2d} \right)$$

where $d$ is the dimensionality of the system, $\rho$ is the number density, and $\beta = \frac{1}{kT}$.

Attention

In MPI parallel execution, betaP is available on rank 0 only. betaP is None on ranks >= 1.

(Loggable: category=”scalar”)

Type:

float

property sdf_compression

$s_\mathrm{comp}[k]$ - The scale distribution function for compression moves $[\mathrm{probability\ density}]$.

See also

x_compression defines the bin center locations.

Attention

In MPI parallel execution, the array is available on rank 0 only. sdf_compression is None on ranks >= 1.

(Loggable: category=”sequence”)

Type:

(N_bins,) numpy.ndarray of float)

property sdf_expansion

$s_\mathrm{exp}[k]$ - The scale distribution function for the expansion moves $[\mathrm{probability\ density}]$.

See also

x_expansion defines the bin center locations..

Attention

In MPI parallel execution, the array is available on rank 0 only. sdf_expansion is None on ranks >= 1.

(Loggable: category=”sequence”)

Type:

(N_bins,) numpy.ndarray of float)

property x_compression

The x values at the center of each bin corresponding to the scale distribution function for the compressive perturbations $[\mathrm{length}]$.

(Loggable: category=”sequence”)

Type:

(N_bins,) numpy.ndarray of float)

property x_expansion

The x values at the center of each bin corresponding to the scale distribution function for the expansion moves $[\mathrm{length}]$.

(Loggable: category=”sequence”)

Type:

(N_bins,) numpy.ndarray of float)