Order Module

Overview

freud.order.Cubatic

Compute the cubatic order parameter [HAG15] for a system of particles using simulated annealing instead of Newton-Raphson root finding.

freud.order.Nematic

Compute the nematic order parameter for a system of particles.

freud.order.Hexatic

Calculates the \(k\)-atic order parameter for 2D systems.

freud.order.Translational

Compute the translational order parameter for each particle.

freud.order.Steinhardt

Compute the local Steinhardt [Ste83] rotationally invariant \(Q_l\) \(W_l\) order parameter for a set of points.

freud.order.SolidLiquid

Identifies solid-like clusters using dot products of \(q_{lm}\).

freud.order.RotationalAutocorrelation

Calculates a measure of total rotational autocorrelation.

Details

The freud.order module contains functions which compute order parameters for the whole system or individual particles. Order parameters take bond order data and interpret it in some way to quantify the degree of order in a system using a scalar value. This is often done through computing spherical harmonics of the bond order diagram, which are the spherical analogue of Fourier Transforms.

class freud.order.Cubatic

Bases: freud.util._Compute

Compute the cubatic order parameter [HAG15] for a system of particles using simulated annealing instead of Newton-Raphson root finding.

Parameters
  • t_initial (float) – Starting temperature.

  • t_final (float) – Final temperature.

  • scale (float) – Scaling factor to reduce temperature.

  • n_replicates (unsigned int, optional) – Number of replicate simulated annealing runs. (Default value = 1).

  • seed (unsigned int, optional) – Random seed to use in calculations. If None, system time is used. (Default value = None).

compute

Calculates the per-particle and global order parameter.

Parameters

orientations ((\(N_{particles}\), 4) numpy.ndarray) – Orientations as angles to use in computation.

property cubatic_tensor

Rank 4 homogeneous tensor representing the optimal system-wide coordinates.

Type

\(\left(3, 3, 3, 3 \right)\) numpy.ndarray

property global_tensor

Rank 4 tensor corresponding to the global orientation. Computed from all orientations.

Type

\(\left(3, 3, 3, 3 \right)\) numpy.ndarray

property order

Cubatic order parameter of the system.

Type

float

property orientation

The quaternion of global orientation.

Type

\(\left(4 \right)\) numpy.ndarray

property particle_order

Order parameter.

Type

\(\left(N_{particles} \right)\) numpy.ndarray

scale

The scale.

Type

float

seed

Random seed to use in calculations.

Type

unsigned int

t_final

The value of the final temperature.

Type

float

t_initial

The value of the initial temperature.

Type

float

class freud.order.Hexatic

Bases: freud.locality._PairCompute

Calculates the \(k\)-atic order parameter for 2D systems.

The \(k\)-atic order parameter (called the hexatic order parameter for \(k = 6\)) is analogous to Steinhardt order parameters, and is used to measure order in the bonds of 2D systems.

The \(k\)-atic order parameter for a particle \(i\) and its \(n\) neighbors \(j\) is given by:

\(\psi_k \left( i \right) = \frac{1}{n} \sum_j^n e^{k i \phi_{ij}}\)

The parameter \(k\) governs the symmetry of the order parameter and typically matches the number of neighbors to be found for each particle. The quantity \(\phi_{ij}\) is the angle between the vector \(r_{ij}\) and \(\left( 1,0 \right)\).

Note

2D: freud.order.Hexatic is only defined for 2D systems. The points must be passed in as [x, y, 0].

Parameters

k (unsigned int, optional) – Symmetry of order parameter. (Default value = 6).

compute

Calculates the hexatic order parameter.

Parameters
default_query_args

The default query arguments are {'mode': 'nearest', 'num_neighbors': self.k}.

k

Symmetry of the order parameter.

Type

unsigned int

property particle_order

Order parameter.

Type

\(\left(N_{particles} \right)\) numpy.ndarray

class freud.order.Nematic

Bases: freud.util._Compute

Compute the nematic order parameter for a system of particles.

Parameters

u (\(\left(3 \right)\) numpy.ndarray) – The nematic director of a single particle in the reference state (without any rotation applied).

compute

Calculates the per-particle and global order parameter.

Parameters

orientations (\(\left(N_{particles}, 4 \right)\) numpy.ndarray) – Orientations to calculate the order parameter.

property director

The average nematic director.

Type

\(\left(3 \right)\) numpy.ndarray

property nematic_tensor

3x3 matrix corresponding to the average particle orientation.

Type

\(\left(3, 3 \right)\) numpy.ndarray

property order

Nematic order parameter of the system.

Type

float

property particle_tensor

\(\left(N_{particles}, 3, 3 \right)\) numpy.ndarray: One 3x3 matrix per-particle corresponding to each individual particle orientation.

u

The normalized reference director (the normalized vector provided on construction).

Type

\(\left(3 \right)\) numpy.ndarray

class freud.order.RotationalAutocorrelation

Bases: freud.util._Compute

Calculates a measure of total rotational autocorrelation.

For any calculation of rotational correlations of extended (i.e. non-point) particles, encoding the symmetries of these particles is crucial to appropriately determining correlations. For systems of anisotropic particles in three dimensions, representing such equivalence can be quite mathematically challenging. This calculation is based on the hyperspherical harmonics as laid out in [KDvAG19]. Generalizations of spherical harmonics to four dimensions, hyperspherical harmonics provide a natural basis for periodic functions in 4 dimensions just as harmonic functions (sines and cosines) or spherical harmonics do in lower dimensions. The idea behind this calculation is to embed orientation quaternions into a 4-dimensional space and then use hyperspherical harmonics to find correlations in a symmetry-aware fashion.

The choice of the hyperspherical harmonic parameter \(l\) determines the symmetry of the functions. The output is not a correlation function, but rather a scalar value that measures total system orientational correlation with an initial state. As such, the output can be treated as an order parameter measuring degrees of rotational (de)correlation. For analysis of a trajectory, the compute call needs to be done at each trajectory frame.

Parameters

l (int) – Order of the hyperspherical harmonic. Must be a positive, even integer.

compute

Calculates the rotational autocorrelation function for a single frame.

Parameters
  • ref_orientations ((\(N_{orientations}\), 4) numpy.ndarray) – Orientations for the initial frame.

  • orientations ((\(N_{orientations}\), 4) numpy.ndarray) – Orientations for the frame of interest.

l

The azimuthal quantum number, which defines the order of the hyperspherical harmonic.

Type

int

property order

Autocorrelation of the system.

Type

float

property particle_order

Rotational autocorrelation values calculated for each orientation.

Type

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

class freud.order.SolidLiquid

Bases: freud.locality._PairCompute

Identifies solid-like clusters using dot products of \(q_{lm}\).

The solid-liquid order parameter [tW95][FHND10] uses a Steinhardt-like approach to identify solid-like particles. First, a bond parameter \(q_l(i, j)\) is computed for each neighbor bond.

If normalize_q is true (default), the bond parameter is given by \(q_l(i, j) = \frac{\sum_{m=-l}^{l} \text{Re}~q_{lm}(i) q_{lm}^*(j)} {\sqrt{\sum_{m=-l}^{l} \lvert q_{lm}(i) \rvert^2} \sqrt{\sum_{m=-l}^{l} \lvert q_{lm}(j) \rvert^2}}\)

If normalize_q is false, then the denominator of the above expression is left out.

Next, the bonds are filtered to keep only “solid-like” bonds with \(q_l(i, j)\) above a cutoff value \(q_{threshold}\).

If a particle has more than \(S_{threshold}\) solid-like bonds, then the particle is considered solid-like. Finally, solid-like particles are clustered.

Parameters
  • l (unsigned int) – Spherical harmonic quantum number l.

  • q_threshold (float) – Value of dot product threshold when evaluating \(q_l(i, j)\) to determine if a bond is solid-like. For \(l=6\), 0.7 is generally good for FCC or BCC structures [FHND10].

  • solid_threshold (unsigned int) – Minimum required number of adjacent solid-like bonds for a particle to be considered solid-like for clustering. For \(l=6\), 6-8 is generally good for FCC or BCC structures.

  • normalize_q (bool) – Whether to normalize the dot product (Default value = True).

property cluster_idx

\(\left(N_{particles}\right)\) numpy.ndarray: Solid-like cluster indices for each particle.

property cluster_sizes

The sizes of all clusters.

Type

\((N_{clusters}, )\) np.ndarray

compute

Compute the order parameter.

Parameters
default_query_args

No default query arguments.

l

Spherical harmonic quantum number l.

Type

unsigned int

property largest_cluster_size

The largest cluster size.

Type

unsigned int

property nlist

Neighbor list of solid-like bonds.

Type

freud.locality.NeighborList

normalize_q

Whether the dot product is normalized.

Type

bool

property num_connections

The number of solid-like bonds for each particle.

Type

\(\left(N_{particles}\right)\) numpy.ndarray

plot

Plot solid-like cluster distribution.

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)

q_threshold

Value of dot product threshold.

Type

float

property ql_ij

Bond dot products \(q_l(i, j)\). Indexed by the elements of self.nlist.

Type

\(\left(N_{bonds}\right)\) numpy.ndarray

solid_threshold

Value of number-of-bonds threshold.

Type

float

class freud.order.Steinhardt

Bases: freud.locality._PairCompute

Compute the local Steinhardt [Ste83] rotationally invariant \(Q_l\) \(W_l\) order parameter for a set of points.

Implements the local rotationally invariant \(q_l\) or \(w_l\) order parameter described by Steinhardt. For a particle i, we calculate the average order parameter by summing the spherical harmonics between particle \(i\) and its neighbors \(j\) in a local region: \(\overline{q}_{lm}(i) = \frac{1}{N_b} \displaystyle\sum_{j=1}^{N_b} Y_{lm}(\theta(\vec{r}_{ij}), \phi(\vec{r}_{ij}))\). The particles included in the sum are determined by the r_max argument to the constructor.

For \(q_l\), this is then combined in a rotationally invariant fashion to remove local orientational order as follows: \(q_l(i)=\sqrt{\frac{4\pi}{2l+1} \displaystyle\sum_{m=-l}^{l} |\overline{q}_{lm}|^2 }\).

For \(w_l\), it is then defined as a weighted average over the \(\overline{q}_{lm}(i)\) values using Wigner 3j symbols (Clebsch-Gordan coefficients). The resulting combination is rotationally (i.e. frame) invariant.

The average argument in the constructor provides access to a variant of this parameter that performs a average over the first and second shell combined [LD08]. To compute this parameter, we perform a second averaging over the first neighbor shell of the particle to implicitly include information about the second neighbor shell. This averaging is performed by replacing the value \(\overline{q}_{lm}(i)\) in the original definition by the average value of \(\overline{q}_{lm}(k)\) over all the \(k\) neighbors of particle \(i\) as well as itself.

The norm attribute argument provides normalized versions of the order parameter, where the normalization is performed by averaging the \(q_{lm}\) values over all particles before computing the order parameter of choice.

Parameters
  • l (unsigned int) – Spherical harmonic quantum number l.

  • average (bool, optional) – Determines whether to calculate the averaged Steinhardt order parameter. (Default value = False)

  • wl (bool, optional) – Determines whether to use the \(w_l\) version of the Steinhardt order parameter. (Default value = False)

  • weighted (bool, optional) – Determines whether to use neighbor weights in the computation of spherical harmonics over neighbors. If enabled and used with a Voronoi neighbor list, this results in the Minkowski Structure Metrics \(q'_l\). (Default value = False)

  • wl_normalize (bool, optional) – Determines whether to normalize the \(w_l\) version of the Steinhardt order parameter. (Default value = False)

average

Whether the the averaged Steinhardt order parameter was calculated.

Type

bool

compute

Compute the order parameter.

Parameters
default_query_args

No default query arguments.

l

Spherical harmonic quantum number l.

Type

unsigned int

property order

The system wide normalization of the \(q_l\) or \(w_l\) order parameter.

Type

float

property particle_order

Variant of the Steinhardt order parameter for each particle (filled with nan for particles with no neighbors).

Type

\(\left(N_{particles}\right)\) numpy.ndarray

plot

Plot order parameter distribution.

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 ql

\(\left(N_{particles}\right)\) numpy.ndarray: \(q_l\) Steinhardt order parameter for each particle (filled with nan for particles with no neighbors). This is always available, no matter which options are selected.

weighted

Whether neighbor weights were used in the computation of spherical harmonics over neighbors.

Type

bool

wl

Whether the \(W_l\) version of the Steinhardt order parameter was used.

Type

bool

class freud.order.Translational

Bases: freud.locality._PairCompute

Compute the translational order parameter for each particle.

Note

2D: freud.order.Translational is only defined for 2D systems. The points must be passed in as [x, y, 0].

Parameters

k (float, optional) – Symmetry of order parameter. (Default value = 6.0).

compute

Calculates the local descriptors.

Parameters
default_query_args

{‘mode’: ‘nearest’, ‘num_neighbors’: int(self.k)}.

Type

The default query arguments are

Type

code

k

Normalization of the order parameter.

Type

unsigned int

property particle_order

Order parameter.

Type

\(\left(N_{particles} \right)\) numpy.ndarray