Density Module

Overview

freud.density.CorrelationFunction

Computes the complex pairwise correlation function.

freud.density.GaussianDensity

Computes the density of a system on a grid.

freud.density.LocalDensity

Computes the local density around a particle.

freud.density.RDF

Computes the RDF \(g \left( r \right)\) for supplied data.

Details

The freud.density module contains various classes relating to the density of the system. These functions allow evaluation of particle distributions with respect to other particles.

class freud.density.CorrelationFunction

Bases: freud.locality._SpatialHistogram1D

Computes the complex pairwise correlation function.

The correlation function is given by \(C(r) = \left\langle s^*_1(0) \cdot s_2(r) \right\rangle\) between two sets of points \(p_1\) (points) and \(p_2\) (query_points) with associated values \(s_1\) (values) and \(s_2\) (query_values). Computing the correlation function results in an array of the expected (average) product of all values at a given radial distance \(r\).

The values of \(r\) where the correlation function is computed are controlled by the r_max and dr parameters to the constructor. r_max determines the maximum distance at which to compute the correlation function and dr is the step size for each bin.

Note

Self-correlation: It is often the case that we wish to compute the correlation function of a set of points with itself. If query_points is the same as points, not provided, or None, we omit accumulating the self-correlation value in the first bin.

Parameters
  • bins (unsigned int) – The number of bins in the RDF.

  • r_max (float) – Maximum pointwise distance to include in the calculation.

bin_centers

The centers of each bin in the histogram.

Type

\((N_{bins}, )\) numpy.ndarray

property bin_counts

The bin counts in the histogram.

Type

numpy.ndarray

bin_edges

The edges of each bin in the histogram. Is one element larger becauseeach bin has a lower and upper bound.

Type

\((N_{bins}+1, )\) numpy.ndarray

bounds

A tuple indicating upper and lower bounds of the histogram.

Type

tuple

property box

The box object used in the last computation.

Type

freud.box.Box

compute

Calculates the correlation function and adds to the current histogram.

Parameters
  • system – Any object that is a valid argument to freud.locality.NeighborQuery.from_system.

  • values ((\(N_{points}\)) numpy.ndarray) – Values associated with the system points used to calculate the correlation function.

  • query_points ((\(N_{query\_points}\), 3) numpy.ndarray, optional) – Query points used to calculate the correlation function. Uses the system’s points if None (Default value = None).

  • query_values ((\(N_{query\_points}\)) numpy.ndarray, optional) – Query values used to calculate the correlation function. Uses values if None. (Default value = None).

  • neighbors (freud.locality.NeighborList or dict, optional) – Either a NeighborList of neighbor pairs to use in the calculation, or a dictionary of query arguments (Default value: None).

  • reset (bool) – Whether to erase the previously computed values before adding the new computation; if False, will accumulate data (Default value: True).

property correlation

Expected (average) product of all values at a given radial distance.

Type

(\(N_{bins}\)) numpy.ndarray

default_query_args

The default query arguments are {'mode': 'ball', 'r_max': self.r_max}.

nbins

The number of bins in the histogram

Type

int

plot

Plot complex correlation function.

Parameters

ax (matplotlib.axes.Axes, optional) – Axis to plot on. If None, make a new figure and axis. (Default value = None)

Returns

Axis with the plot.

Return type

(matplotlib.axes.Axes)

class freud.density.GaussianDensity

Bases: freud.util._Compute

Computes the density of a system on a grid.

Replaces particle positions with a Gaussian blur and calculates the contribution from each to the proscribed grid based upon the distance of the grid cell from the center of the Gaussian. The resulting data is a regular grid of particle densities that can be used in standard algorithms requiring evenly spaced point, such as Fast Fourier Transforms. The dimensions of the image (grid) are set in the constructor, and can either be set equally for all dimensions or for each dimension independently.

Parameters
  • width (int or list or tuple) – The number of bins to make the image in each direction (identical in all dimensions if a single integer value is provided).

  • r_max (float) – Distance over which to blur.

  • sigma (float) – Sigma parameter for Gaussian.

property box

Box used in the calculation.

Type

freud.box.Box

compute

Calculates the Gaussian blur for the specified points.

Parameters

system – Any object that is a valid argument to freud.locality.NeighborQuery.from_system.

property density

The image grid with the Gaussian density.

Type

(\(w_x\), \(w_y\), \(w_z\)) numpy.ndarray

plot

Plot Gaussian Density.

Parameters

ax (matplotlib.axes.Axes, optional) – Axis to plot on. If None, make a new figure and axis. (Default value = None)

Returns

Axis with the plot.

Return type

(matplotlib.axes.Axes)

r_max

Distance over which to blur.

Type

float

sigma

Sigma parameter for Gaussian.

Type

float

width

The number of bins to make the image in each direction (identical in all dimensions if a single integer value is provided).

Type

int or list or tuple

class freud.density.LocalDensity

Bases: freud.locality._PairCompute

Computes the local density around a particle.

The density of the local environment is computed and averaged for a given set of query points in a sea of data points. Providing the same points calculates them against themselves. Computing the local density results in an array listing the value of the local density around each query point. Also available is the number of neighbors for each query point, giving the user the ability to count the number of particles in that region. Note that the computed density is essentially a number density (that allows for fractional values as described below). If your particles have a specific volume, you can compute a volume density by simply multiplying the output by the volume of the particles.

In order to provide sufficiently smooth data, data points can be fractionally counted towards the density. Rather than perform compute-intensive area (volume) overlap calculations to determine the exact amount of overlap area (volume), the LocalDensity class performs a simple linear interpolation relative to the centers of the data points. Specifically, a point is counted as one neighbor of a given query point if it is entirely contained within the r_max, half of a neighbor if the distance to its center is exactly r_max, and zero if its center is a distance greater than or equal to r_max + diameter from the query point’s center. Graphically, this looks like:

_images/density.png
Parameters
  • r_max (float) – Maximum distance over which to calculate the density.

  • diameter (float) – Diameter of particle circumsphere.

property box

Box used in the calculation.

Type

freud.box.Box

compute

Calculates the local density for the specified points.

Parameters
default_query_args

The default query arguments are {'mode': 'ball', 'r_max': self.r_max + 0.5*self.diameter}.

property density

Density of points per query point.

Type

(\(N_{points}\)) numpy.ndarray

diameter

Diameter of particle circumsphere.

Type

float

property num_neighbors

Number of neighbor points for each query point.

Type

(\(N_{points}\)) numpy.ndarray

r_max

Maximum distance over which to calculate the density.

Type

float

class freud.density.RDF

Bases: freud.locality._SpatialHistogram1D

Computes the RDF \(g \left( r \right)\) for supplied data.

Note that the RDF is defined strictly according to the pair correlation function, i.e.

\[g(r) = V\frac{N-1}{N} \langle \delta(r) \rangle\]

In the thermodynamic limit, the fraction tends to unity and the limiting behavior of \(\lim_{r \to \infty} g(r)=1\) is recovered. However, for very small systems the long range behavior of the radial distribution will instead tend to \(\frac{N-1}{N}\). If you are analyzing a very small system but wish to recover the more familiar behavior, you may use the normalize flag to enforce this requirement upon construction of this object. Note that this will have little to no effect on larger systems (for example, for systems of 100 particles the RDF will differ by 1%).

Note

2D: freud.density.RDF properly handles 2D boxes. The points must be passed in as [x, y, 0].

Parameters
  • bins (unsigned int) – The number of bins in the RDF.

  • r_max (float) – Maximum interparticle distance to include in the calculation.

  • r_min (float, optional) – Minimum interparticle distance to include in the calculation (Default value = 0).

  • normalize (bool, optional) – Scale the RDF values by \(\frac{N_{query\_points}+1}{N_{query\_points}+1}\). This argument primarily exists to deal with standard RDF calculations where no special query_points or neighbors are provided, but where the number of query_points is small enough that the long-ranged limit of \(g(r)\) deviates significantly from \(1\). It should not be used if query_points is provided as a different set of points, or if unusual query arguments are provided to compute(), specifically if :code`exclude_ii` is set to False. This normalization is not meaningful in such cases and will simply convolute the data.

bin_centers

The centers of each bin in the histogram.

Type

\((N_{bins}, )\) numpy.ndarray

property bin_counts

The bin counts in the histogram.

Type

numpy.ndarray

bin_edges

The edges of each bin in the histogram. Is one element larger becauseeach bin has a lower and upper bound.

Type

\((N_{bins}+1, )\) numpy.ndarray

bounds

A tuple indicating upper and lower bounds of the histogram.

Type

tuple

property box

The box object used in the last computation.

Type

freud.box.Box

compute

Calculates the RDF and adds to the current RDF histogram.

Parameters
  • system – Any object that is a valid argument to freud.locality.NeighborQuery.from_system.

  • query_points ((\(N_{query\_points}\), 3) numpy.ndarray, optional) – Query points used to calculate the RDF. Uses the system’s points if None (Default value = None).

  • neighbors (freud.locality.NeighborList or dict, optional) –

    Either a NeighborList of neighbor pairs to use in the calculation, or a dictionary of query arguments (Default value: None).

  • reset (bool) – Whether to erase the previously computed values before adding the new computation; if False, will accumulate data (Default value: True).

default_query_args

The default query arguments are {'mode': 'ball', 'r_max': self.r_max}.

property n_r

Histogram of cumulative bin_counts values. More precisely, n_r[i] is the average number of points contained within a ball of radius R[i]+dr/2 centered at a given query_point averaged over all query_points in the last call to compute().

Type

(\(N_{bins}\),) numpy.ndarray

nbins

The number of bins in the histogram

Type

int

plot

Plot radial distribution function.

Parameters

ax (matplotlib.axes.Axes, optional) – Axis to plot on. If None, make a new figure and axis. (Default value = None)

Returns

Axis with the plot.

Return type

(matplotlib.axes.Axes)

property rdf

Histogram of RDF values.

Type

(\(N_{bins}\),) numpy.ndarray