Order Module¶
The order module contains functions which compute order parameters for the whole system or individual particles.
Bond Order¶
-
class
freud.order.
BondOrder
(rmax, k, n, nBinsT, nBinsP)¶ Compute the bond order diagram for the system of particles.
Available modes of calculation:
- If
mode='bod'
(Bond Order Diagram, default): Create the 2D histogram containing the number of bonds formed through the surface of a unit sphere based on the azimuthal \(\left( \theta \right)\) and polar \(\left( \phi \right)\) angles. - If
mode='lbod'
(Local Bond Order Diagram): Create the 2D histogram containing the number of bonds formed, rotated into the local orientation of the central particle, through the surface of a unit sphere based on the azimuthal \(\left( \theta \right)\) and polar \(\left( \phi \right)\) angles. - If
mode='obcd'
(Orientation Bond Correlation Diagram): Create the 2D histogram containing the number of bonds formed, rotated by the rotation that takes the orientation of neighboring particle j to the orientation of each particle i, through the surface of a unit sphere based on the azimuthal \(\left( \theta \right)\) and polar \(\left( \phi \right)\) angles. - If
mode='oocd'
(Orientation Orientation Correlation Diagram): Create the 2D histogram containing the directors of neighboring particles (\(\hat{z}\) rotated by their quaternion), rotated into the local orientation of the central particle, through the surface of a unit sphere based on the azimuthal \(\left( \theta \right)\) and polar \(\left( \phi \right)\) angles.
Module author: Erin Teich <erteich@umich.edu>
Parameters: - r_max (float) – distance over which to calculate
- k (unsigned int) – order parameter i. to be removed
- n (unsigned int) – number of neighbors to find
- n_bins_t (unsigned int) – number of theta bins
- n_bins_p (unsigned int) – number of phi bins
-
accumulate
(self, box, ref_points, ref_orientations, points, orientations, str mode='bod', nlist=None)¶ Calculates the correlation function and adds to the current histogram.
Parameters: - box (
freud.box.Box
) – simulation box - ref_points (
numpy.ndarray
, shape=(\(N_{particles}\), 3), dtype=numpy.float32
) – reference points to calculate the local density - ref_orientations (
numpy.ndarray
, shape=(\(N_{particles}\), 4), dtype=numpy.float32
) – orientations to use in computation - points (
numpy.ndarray
, shape=(\(N_{particles}\), 3), dtype=numpy.float32
) – points to calculate the local density - orientations (
numpy.ndarray
, shape=(\(N_{particles}\), 4), dtype=numpy.float32
) – orientations to use in computation - mode (str) – mode to calc bond order. “bod”, “lbod”, “obcd”, and “oocd”
- nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- box (
-
bond_order
¶ Bond order.
-
box
¶ Box used in the calculation.
-
compute
(self, box, ref_points, ref_orientations, points, orientations, mode='bod', nlist=None)¶ Calculates the bond order histogram. Will overwrite the current histogram.
Parameters: - box (
freud.box.Box
) – simulation box - ref_points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – reference points to calculate the local density - ref_orientations (
numpy.ndarray
, shape= \(\left(N_{particles}, 4\right)\), dtype=numpy.float32
) – orientations to use in computation - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – points to calculate the local density - orientations (
numpy.ndarray
, shape= \(\left(N_{particles}, 4 \right)\), dtype=numpy.float32
) – orientations to use in computation - mode (str) – mode to calc bond order. “bod”, “lbod”, “obcd”, and “oocd”
- nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- box (
-
getBondOrder
(self)¶ Get the bond order.
Returns: bond order Return type: numpy.ndarray
, shape= \(\left(N_{\phi}, N_{\theta} \right)\), dtype=numpy.float32
-
getBox
(self)¶ Get the box used in the calculation.
Returns: freud Box Return type: freud.box.Box
-
getNBinsPhi
(self)¶ Get the number of bins in the Phi-dimension of histogram.
Returns: \(N_{\phi}\) Return type: unsigned int
-
getNBinsTheta
(self)¶ Get the number of bins in the Theta-dimension of histogram.
Returns: \(N_{\theta}\) Return type: unsigned int
-
getPhi
(self)¶ Returns: values of bin centers for Phi Return type: numpy.ndarray
, shape= \(\left(N_{\phi} \right)\), dtype=numpy.float32
-
getTheta
(self)¶ Returns: values of bin centers for Theta Return type: numpy.ndarray
, shape= \(\left(N_{\theta} \right)\), dtype=numpy.float32
-
reduceBondOrder
(self)¶ Reduces the histogram in the values over N processors to a single histogram. This is called automatically by
freud.order.BondOrder.getBondOrder()
.
-
resetBondOrder
(self)¶ resets the values of the bond order in memory
- If
Order Parameters¶
Order parameters take bond order data and interpret it in some way to quantify the degree of order in a system. This is often done through computing spherical harmonics of the bond order diagram, which are the spherical analogue of Fourier Transforms.
Cubatic Order Parameter¶
-
class
freud.order.
CubaticOrderParameter
(t_initial, t_final, scale, n_replicates, seed)¶ Compute the cubatic order parameter [Cit1] for a system of particles using simulated annealing instead of Newton-Raphson root finding.
Module author: Eric Harper <harperic@umich.edu>
Parameters: -
compute
(self, orientations)¶ Calculates the per-particle and global order parameter.
Parameters: - box (
freud.box.Box
) – simulation box - orientations (
numpy.ndarray
, shape= \(\left(N_{particles}, 4 \right)\), dtype=numpy.float32
) – orientations to calculate the order parameter
- box (
-
get_cubatic_tensor
(self)¶ Returns: Rank 4 tensor corresponding to each individual particle orientation Return type: numpy.ndarray
, shape= \(\left(3, 3, 3, 3 \right)\), dtype=numpy.float32
-
get_gen_r4_tensor
(self)¶ Returns: Rank 4 tensor corresponding to each individual particle orientation Return type: numpy.ndarray
, shape= \(\left(3, 3, 3, 3 \right)\), dtype=numpy.float32
-
get_global_tensor
(self)¶ Returns: Rank 4 tensor corresponding to each individual particle orientation Return type: numpy.ndarray
, shape= \(\left(3, 3, 3, 3 \right)\), dtype=numpy.float32
-
get_orientation
(self)¶ Returns: orientation of global orientation Return type: numpy.ndarray
, shape= \(\left(4 \right)\), dtype=numpy.float32
-
get_particle_tensor
(self)¶ Returns: Rank 4 tensor corresponding to each individual particle orientation Return type: numpy.ndarray
, shape= \(\left(N_{particles}, 3, 3, 3, 3 \right)\), dtype=numpy.float32
-
Nematic Order Parameter¶
-
class
freud.order.
NematicOrderParameter
(u)¶ Compute the nematic order parameter for a system of particles.
Module author: Jens Glaser <jsglaser@umich.edu>
New in version 0.7.0.
Parameters: u ( numpy.ndarray
, shape= \(\left(3 \right)\), dtype=numpy.float32
) – The nematic director of a single particle in the reference state (without any rotation applied)-
compute
(self, orientations)¶ Calculates the per-particle and global order parameter.
Parameters: orientations ( numpy.ndarray
, shape= \(\left(N_{particles}, 4 \right)\), dtype=numpy.float32
) – orientations to calculate the order parameter
-
get_director
(self)¶ The director (eigenvector corresponding to the order parameter).
Returns: The average nematic director Return type: numpy.ndarray
, shape= \(\left(3 \right)\), dtype=numpy.float32
-
get_nematic_order_parameter
(self)¶ The nematic order parameter.
Returns: Nematic order parameter Return type: float
-
get_nematic_tensor
(self)¶ The nematic Q tensor.
Returns: 3x3 matrix corresponding to the average particle orientation Return type: numpy.ndarray
, shape= \(\left(3, 3 \right)\), dtype=numpy.float32
-
get_particle_tensor
(self)¶ The full per-particle tensor of orientation information.
Returns: 3x3 matrix corresponding to each individual particle orientation Return type: numpy.ndarray
, shape= \(\left(N_{particles}, 3, 3 \right)\), dtype=numpy.float32
-
Hexatic Order Parameter¶
-
class
freud.order.
HexOrderParameter
(rmax, k, n)¶ Calculates the \(k\)-atic order parameter for each particle in the system.
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 while the parameter \(n\) governs the number of neighbors of particle \(i\) to average over. \(\phi_{ij}\) is the angle between the vector \(r_{ij}\) and \(\left( 1,0 \right)\)
Note
2D: This calculation is defined for 2D systems only. However, particle positions are still required to be passed in as
[x, y, 0]
.Module author: Eric Harper <harperic@umich.edu>
Parameters: - rmax (float) – +/- r distance to search for neighbors
- k (unsigned int) – symmetry of order parameter (\(k=6\) is hexatic)
- n (unsigned int) – number of neighbors (\(n=k\) if \(n\) not specified)
-
box
¶ Get the box used in the calculation.
-
compute
(self, box, points, nlist=None)¶ Calculates the correlation function and adds to the current histogram.
Parameters: - box (
freud.box.Box
) – simulation box - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- box (
-
getBox
(self)¶ Get the box used in the calculation.
Returns: freud Box Return type: freud.box.Box
-
getK
(self)¶ Get the symmetry of the order parameter.
Returns: \(k\) Return type: unsigned int
-
getNP
(self)¶ Get the number of particles.
Returns: \(N_{particles}\) Return type: unsigned int
-
getPsi
(self)¶ Get the order parameter.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles} \right)\), dtype=numpy.complex64
-
k
¶ Symmetry of the order parameter.
-
num_particles
¶ Get the number of particles.
-
psi
¶ Order parameter.
Local Descriptors¶
-
class
freud.order.
LocalDescriptors
(box, nNeigh, lmax, rmax)¶ Compute a set of descriptors (a numerical “fingerprint”) of a particle’s local environment.
Module author: Matthew Spellings <mspells@umich.edu>
Parameters: - num_neighbors (unsigned int) – Maximum number of neighbors to compute descriptors for
- lmax (unsigned int) – Maximum spherical harmonic \(l\) to consider
- rmax (float) – Initial guess of the maximum radius to looks for neighbors
- negative_m (bool) – True if we should also calculate \(Y_{lm}\) for negative \(m\)
-
compute
(self, box, unsigned int num_neighbors, points_ref, points=None, orientations=None, mode='neighborhood', nlist=None)¶ Calculates the local descriptors of bonds from a set of source points to a set of destination points.
Parameters: - num_neighbors – Number of neighbors to compute with
- points_ref (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – source points to calculate the order parameter - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – destination points to calculate the order parameter - orientations (
numpy.ndarray
, shape= \(\left(N_{particles}, 4 \right)\), dtype=numpy.float32
or None) – Orientation of each reference point - mode (str) – Orientation mode to use for environments, either ‘neighborhood’ to use the orientation of the local neighborhood, ‘particle_local’ to use the given particle orientations, or ‘global’ to not rotate environments
- nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
-
computeNList
(self, box, points_ref, points=None)¶ Compute the neighbor list for bonds from a set of source points to a set of destination points.
Parameters: - num_neighbors – Number of neighbors to compute with
- points_ref (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – source points to calculate the order parameter - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – destination points to calculate the order parameter
-
getLMax
(self)¶ Get the maximum spherical harmonic \(l\) to calculate for.
Returns: \(l\) Return type: unsigned int
-
getNP
(self)¶ Get the number of particles.
Returns: \(N_{particles}\) Return type: unsigned int
-
getNSphs
(self)¶ Get the number of neighbors.
Returns: \(N_{neighbors}\) Return type: unsigned int
-
getSph
(self)¶ Get a reference to the last computed spherical harmonic array.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{bonds}, \text{SphWidth} \right)\), dtype=numpy.complex64
-
l_max
¶ Get the maximum spherical harmonic \(l\) to calculate for.
-
num_neighbors
¶ Get the number of neighbors.
-
num_particles
¶ Get the number of particles.
-
r_max
¶ Get the cutoff radius.
-
sph
¶ A reference to the last computed spherical harmonic array.
Translational Order Parameter¶
-
class
freud.order.
TransOrderParameter
(rmax, k, n)¶ Compute the translational order parameter for each particle.
Module author: Michael Engel <engelmm@umich.edu>
Parameters: -
box
¶ Get the box used in the calculation.
-
compute
(self, box, points, nlist=None)¶ Calculates the local descriptors.
Parameters: - box (
freud.box.Box
) – simulation box - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- box (
-
d_r
¶ Get a reference to the last computed spherical harmonic array.
-
getBox
(self)¶ Get the box used in the calculation.
Returns: freud Box Return type: freud.box.Box
-
getDr
(self)¶ Get a reference to the last computed spherical harmonic array.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.complex64
-
getNP
(self)¶ Get the number of particles.
Returns: \(N_{particles}\) Return type: unsigned int
-
num_particles
¶ Get the number of particles.
-
Local \(Q_l\)¶
-
class
freud.order.
LocalQl
(box, rmax, l, rmin)¶ LocalQl(box, rmax, l, rmin=0)
Compute the local Steinhardt rotationally invariant \(Q_l\) [Cit4] order parameter for a set of points.
Implements the local rotationally invariant \(Q_l\) order parameter described by Steinhardt. For a particle i, we calculate the average \(Q_l\) 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}))\)
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 more details see PJ Steinhardt (1983) (DOI: 10.1103/PhysRevB.28.784)
Added first/second shell combined average \(Q_l\) order parameter for a set of points:
- Variation of the Steinhardt \(Q_l\) order parameter
- For a particle i, we calculate the average \(Q_l\) by summing the spherical harmonics between particle i and its neighbors j and the neighbors k of neighbor j in a local region
Module author: Xiyu Du <xiyudu@umich.edu>
Parameters: - box (
freud.box.Box
) – simulation box - rmax (float) – Cutoff radius for the local order parameter. Values near first minima of the RDF are recommended
- l (unsigned int) – Spherical harmonic quantum number l. Must be a positive number
- rmin (float) – can look at only the second shell or some arbitrary RDF region
-
Ql
¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
-
ave_Ql
¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
-
ave_norm_Ql
¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
-
box
¶ Get the box used in the calculation.
-
compute
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeAve
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeAveNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3 \right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
getAveQl
(self)¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
-
getBox
(self)¶ Get the box used in the calculation.
Returns: freud Box Return type: freud.box.Box
-
getNP
(self)¶ Get the number of particles.
Returns: \(N_p\) Return type: unsigned int
-
getQl
(self)¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
-
getQlAveNorm
(self)¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
-
getQlNorm
(self)¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
-
norm_Ql
¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
-
num_particles
¶ Get the number of particles.
-
setBox
(self, box)¶ Reset the simulation box.
Parameters: box ( freud.box.Box
) – simulation box
Nearest Neighbors Local \(Q_l\)¶
-
class
freud.order.
LocalQlNear
(box, rmax, l, kn)¶ LocalQlNear(box, rmax, l, kn=12)
Compute the local Steinhardt rotationally invariant \(Q_l\) order parameter [Cit4] for a set of points.
Implements the local rotationally invariant \(Q_l\) order parameter described by Steinhardt. For a particle i, we calculate the average \(Q_l\) 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}))\)
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 more details see PJ Steinhardt (1983) (DOI: 10.1103/PhysRevB.28.784)
Added first/second shell combined average \(Q_l\) order parameter for a set of points:
- Variation of the Steinhardt \(Q_l\) order parameter
- For a particle i, we calculate the average \(Q_l\) by summing the spherical harmonics between particle i and its neighbors j and the neighbors k of neighbor j in a local region
Module author: Xiyu Du <xiyudu@umich.edu>
Parameters: - box (
freud.box.Box
) – simulation box - rmax (float) – Cutoff radius for the local order parameter. Values near first minima of the RDF are recommended
- l (unsigned int) – Spherical harmonic quantum number l. Must be a positive number
- kn (unsigned int) – number of nearest neighbors. must be a positive integer
-
compute
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeAve
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeAveNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
Local \(W_l\)¶
-
class
freud.order.
LocalWl
(box, rmax, l)¶ LocalWl(box, rmax, l)
Compute the local Steinhardt rotationally invariant \(W_l\) order parameter [Cit4] for a set of points.
Implements the local rotationally invariant \(W_l\) order parameter described by Steinhardt that can aid in distinguishing between FCC, HCP, and BCC.
For more details see PJ Steinhardt (1983) (DOI: 10.1103/PhysRevB.28.784)
Added first/second shell combined average \(W_l\) order parameter for a set of points:
- Variation of the Steinhardt \(W_l\) order parameter
- For a particle i, we calculate the average \(W_l\) by summing the spherical harmonics between particle i and its neighbors j and the neighbors k of neighbor j in a local region
Module author: Xiyu Du <xiyudu@umich.edu>
Parameters: - box (
freud.box.Box
) – simulation box - rmax (float) – Cutoff radius for the local order parameter. Values near first minima of the RDF are recommended
- l (unsigned int) – Spherical harmonic quantum number l. Must be a positive number
-
Ql
¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
-
Wl
¶ Get a reference to the last computed \(W_l\) for each particle. Returns NaN instead of \(W_l\) for particles with no neighbors.
-
ave_Wl
¶ Get a reference to the last computed \(W_l\) for each particle. Returns NaN instead of \(W_l\) for particles with no neighbors.
-
ave_norm_Wl
¶ Get a reference to the last computed \(W_l\) for each particle. Returns NaN instead of \(W_l\) for particles with no neighbors.
-
box
¶ Get the box used in the calculation.
-
compute
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeAve
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeAveNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
getAveWl
(self)¶ Get a reference to the last computed \(W_l\) for each particle. Returns NaN instead of \(W_l\) for particles with no neighbors.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
-
getBox
(self)¶ Get the box used in the calculation.
Returns: freud Box Return type: freud.box.Box
-
getNP
(self)¶ Get the number of particles.
Returns: \(N_{particles}\) Return type: unsigned int
-
getQl
(self)¶ Get a reference to the last computed \(Q_l\) for each particle. Returns NaN instead of \(Q_l\) for particles with no neighbors.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
-
getWl
(self)¶ Get a reference to the last computed \(W_l\) for each particle. Returns NaN instead of \(W_l\) for particles with no neighbors.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.complex64
-
getWlAveNorm
(self)¶ Get a reference to the last computed \(W_l\) for each particle. Returns NaN instead of \(W_l\) for particles with no neighbors.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
-
getWlNorm
(self)¶ Get a reference to the last computed \(W_l\) for each particle. Returns NaN instead of \(W_l\) for particles with no neighbors.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
-
norm_Wl
¶ Get a reference to the last computed \(W_l\) for each particle. Returns NaN instead of \(W_l\) for particles with no neighbors.
-
num_particles
¶ Get the number of particles.
-
setBox
(self, box)¶ Reset the simulation box.
Parameters: box ( freud.box.Box
) – simulation box
Nearest Neighbors Local \(W_l\)¶
-
class
freud.order.
LocalWlNear
(box, rmax, l, kn)¶ LocalWlNear(box, rmax, l, kn=12)
Compute the local Steinhardt rotationally invariant \(W_l\) order parameter [Cit4] for a set of points.
Implements the local rotationally invariant \(W_l\) order parameter described by Steinhardt that can aid in distinguishing between FCC, HCP, and BCC.
For more details see PJ Steinhardt (1983) (DOI: 10.1103/PhysRevB.28.784)
Added first/second shell combined average \(W_l\) order parameter for a set of points:
- Variation of the Steinhardt \(W_l\) order parameter
- For a particle i, we calculate the average \(W_l\) by summing the spherical harmonics between particle i and its neighbors j and the neighbors k of neighbor j in a local region
Module author: Xiyu Du <xiyudu@umich.edu>
Parameters: - box (
freud.box.Box
) – simulation box - rmax (float) – Cutoff radius for the local order parameter. Values near first minima of the RDF are recommended
- l (unsigned int) – Spherical harmonic quantum number l. Must be a positive number
- kn (unsigned int) – Number of nearest neighbors. Must be a positive number
-
compute
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeAve
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeAveNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
Solid-Liquid Order Parameter¶
-
class
freud.order.
SolLiq
(box, rmax, Qthreshold, Sthreshold, l)¶ SolLiq(box, rmax, Qthreshold, Sthreshold, l)
Computes dot products of \(Q_{lm}\) between particles and uses these for clustering.
Module author: Richmond Newman <newmanrs@umich.edu>
Parameters: - box (
freud.box.Box
) – simulation box - rmax (float) – Cutoff radius for the local order parameter. Values near first minima of the RDF are recommended
- Qthreshold (float) – Value of dot product threshold when evaluating \(Q_{lm}^*(i) Q_{lm}(j)\) to determine if a neighbor pair is a solid-like bond. (For \(l=6\), 0.7 generally good for FCC or BCC structures)
- Sthreshold (unsigned int) – Minimum required number of adjacent solid-link bonds for a particle to be considered solid-like for clustering. (For \(l=6\), 6-8 generally good for FCC or BCC structures)
- l (unsigned int) – Choose spherical harmonic \(Q_l\). Must be positive and even.
-
Ql_dot_ij
¶ Get a reference to the number of connections per particle.
-
Ql_mi
¶ Get a reference to the last computed \(Q_{lmi}\) for each particle.
-
box
¶ Get the box used in the calculation.
-
cluster_sizes
¶ Return the sizes of all clusters.
-
clusters
¶ Get a reference to the last computed set of solid-like cluster indices for each particle.
-
compute
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeSolLiqNoNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeSolLiqVariant
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
getBox
(self)¶ Get the box used in the calculation.
Returns: freud Box Return type: freud.box.Box
-
getClusterSizes
(self)¶ Return the sizes of all clusters.
Returns: largest cluster size Return type: numpy.ndarray
, shape= \(\left(N_{clusters}\right)\), dtype=numpy.uint32
-
getClusters
(self)¶ Get a reference to the last computed set of solid-like cluster indices for each particle.
Returns: clusters Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.uint32
-
getLargestClusterSize
(self)¶ Returns the largest cluster size. Must call a compute method first.
Returns: largest cluster size Return type: unsigned int
-
getNP
(self)¶ Get the number of particles.
Returns: np Return type: unsigned int
-
getNumberOfConnections
(self)¶ Get a reference to the number of connections per particle.
Returns: clusters Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.uint32
-
getQldot_ij
(self)¶ Get a reference to the qldot_ij values.
Returns: largest cluster size Return type: numpy.ndarray
, shape= \(\left(N_{clusters}\right)\), dtype=numpy.complex64
-
getQlmi
(self)¶ Get a reference to the last computed \(Q_{lmi}\) for each particle.
Returns: order parameter Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.complex64
-
largest_cluster_size
¶ Returns the largest cluster size. Must call a compute method first.
-
num_connections
¶ Get a reference to the number of connections per particle.
-
num_particles
¶ Get the number of particles.
-
setBox
(self, box)¶ Reset the simulation box.
Parameters: box ( freud.box.Box
) – simulation box
- box (
Nearest Neighbors Solid-Liquid Order Parameter¶
-
class
freud.order.
SolLiqNear
(box, rmax, Qthreshold, Sthreshold, l)¶ SolLiqNear(box, rmax, Qthreshold, Sthreshold, l, kn=12)
Computes dot products of \(Q_{lm}\) between particles and uses these for clustering.
Module author: Richmond Newman <newmanrs@umich.edu>
Parameters: - box (
freud.box.Box
) – simulation box - rmax (float) – Cutoff radius for the local order parameter. Values near first minima of the RDF are recommended
- Qthreshold (float) – Value of dot product threshold when evaluating \(Q_{lm}^*(i) Q_{lm}(j)\) to determine if a neighbor pair is a solid-like bond. (For \(l=6\), 0.7 generally good for FCC or BCC structures)
- Sthreshold (unsigned int) – Minimum required number of adjacent solid-link bonds for a particle to be considered solid-like for clustering. (For \(l=6\), 6-8 generally good for FCC or BCC structures)
- l (unsigned int) – Choose spherical harmonic \(Q_l\). Must be positive and even.
- kn (unsigned int) – Number of nearest neighbors. Must be a positive number
-
compute
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeSolLiqNoNorm
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
computeSolLiqVariant
(self, points, nlist=None)¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – points to calculate the order parameter - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
- box (
Environment Matching¶
-
class
freud.order.
MatchEnv
(box, rmax, k)¶ Clusters particles according to whether their local environments match or not, according to various shape matching metrics.
Module author: Erin Teich <erteich@umich.edu>
Parameters: - box (
freud.box.Box
) – Simulation box - rmax (float) – Cutoff radius for cell list and clustering algorithm. Values near first minimum of the RDF are recommended.
- k (unsigned int) – Number of nearest neighbors taken to define the local environment of any given particle.
-
cluster
(self, points, threshold, hard_r=False, registration=False, global_search=False, env_nlist=None, nlist=None)¶ Determine clusters of particles with matching environments.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – particle positions - threshold (float) – maximum magnitude of the vector difference between two vectors, below which they are “matching”
- hard_r (bool) – If True, add all particles that fall within the threshold of m_rmaxsq to the environment
- registration (bool) – If True, first use brute force registration to orient one set of environment vectors with respect to the other set such that it minimizes the RMSD between the two sets.
- global_search (bool) – If True, do an exhaustive search wherein the environments of every single pair of particles in the simulation are compared. If False, only compare the environments of neighboring particles.
- nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find neighbors of every particle, to compare environments - env_nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find the environment of every particle
- points (
-
getClusters
(self)¶ Get a reference to the particles, indexed into clusters according to their matching local environments
Returns: clusters Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.uint32
-
getEnvironment
(self, i)¶ Returns the set of vectors defining the environment indexed by i.
Parameters: i (unsigned int) – environment index Returns: the array of vectors Return type: numpy.ndarray
, shape= \(\left(N_{neighbors}, 3\right)\), dtype=numpy.float32
-
getNP
(self)¶ Get the number of particles.
Returns: \(N_{particles}\) Return type: unsigned int
-
getNumClusters
(self)¶ Get the number of clusters.
Returns: \(N_{clusters}\) Return type: unsigned int
-
getTotEnvironment
(self)¶ Returns the entire m_Np by m_maxk by 3 matrix of all environments for all particles.
Returns: the array of vectors Return type: numpy.ndarray
, shape= \(\left(N_{particles}, N_{neighbors}, 3\right)\), dtype=numpy.float32
-
isSimilar
(self, refPoints1, refPoints2, threshold, registration=False)¶ Test if the motif provided by refPoints1 is similar to the motif provided by refPoints2.
Parameters: - refPoints1 (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – vectors that make up motif 1 - refPoints2 (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – vectors that make up motif 2 - threshold (float) – maximum magnitude of the vector difference between two vectors, below which they are considered “matching”
- registration (bool) – If true, first use brute force registration to orient one set of environment vectors with respect to the other set such that it minimizes the RMSD between the two sets.
Returns: a doublet that gives the rotated (or not) set of refPoints2, and the mapping between the vectors of refPoints1 and refPoints2 that will make them correspond to each other. empty if they do not correspond to each other.
Return type: tuple[(
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
), map[int, int]]- refPoints1 (
-
matchMotif
(self, points, refPoints, threshold, registration=False, nlist=None)¶ Determine clusters of particles that match the motif provided by refPoints.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – particle positions - refPoints (
numpy.ndarray
, shape= \(\left(N_{neighbors}, 3\right)\), dtype=numpy.float32
) – vectors that make up the motif against which we are matching - threshold (float) – maximum magnitude of the vector difference between two vectors, below which they are considered “matching”
- registration (bool) – If true, first use brute force registration to orient one set of environment vectors with respect to the other set such that it minimizes the RMSD between the two sets.
- nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- points (
-
minRMSDMotif
(self, points, refPoints, registration=False, nlist=None)¶ Rotate (if registration=True) and permute the environments of all particles to minimize their RMSD wrt the motif provided by refPoints.
Parameters: - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – particle positions - refPoints (
numpy.ndarray
, shape= \(\left(N_{neighbors}, 3\right)\), dtype=numpy.float32
) – vectors that make up the motif against which we are matching - registration (bool) – If true, first use brute force registration to orient one set of environment vectors with respect to the other set such that it minimizes the RMSD between the two sets.
- nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
Returns: vector of minimal RMSD values, one value per particle.
Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
- points (
-
minimizeRMSD
(self, refPoints1, refPoints2, registration=False)¶ Get the somewhat-optimal RMSD between the set of vectors refPoints1 and the set of vectors refPoints2.
Parameters: - refPoints1 (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – vectors that make up motif 1 - refPoints2 (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – vectors that make up motif 2 - registration (bool) – if true, first use brute force registration to orient one set of environment vectors with respect to the other set such that it minimizes the RMSD between the two sets
Returns: a triplet that gives the associated min_rmsd, rotated (or not) set of refPoints2, and the mapping between the vectors of refPoints1 and refPoints2 that somewhat minimizes the RMSD.
Return type: tuple[float, (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
), map[int, int]]- refPoints1 (
-
num_clusters
¶ Get the number of clusters.
-
num_particles
¶ Get the number of particles.
-
setBox
(self, box)¶ Reset the simulation box.
Parameters: box ( freud.box.Box
) – simulation box
-
tot_environment
¶ Returns the entire m_Np by m_maxk by 3 matrix of all environments for all particles.
- box (
Pairing¶
Note
This module is deprecated and is replaced with Bond Module.
-
class
freud.order.
Pairing2D
(rmax, k, compDotTol)¶ Compute pairs for the system of particles.
Module author: Eric Harper <harperic@umich.edu>
Parameters: -
box
¶ Get the box used in the calculation.
-
compute
(self, box, points, orientations, compOrientations, nlist=None)¶ Calculates the correlation function and adds to the current histogram.
Parameters: - box (
freud.box.Box
) – simulation box - points (
numpy.ndarray
, shape= \(\left(N_{particles}, 3\right)\), dtype=numpy.float32
) – reference points to calculate the local density - orientations (
numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
) – orientations to use in computation - compOrientations (
numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.float32
) – possible orientations to check for bonds - nlist (
freud.locality.NeighborList
) –freud.locality.NeighborList
object to use to find bonds
- box (
-
getBox
(self)¶ Get the box used in the calculation.
Returns: freud Box Return type: freud.box.Box
-
getMatch
(self)¶ Get the match.
Returns: match Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.uint32
-
getPair
(self)¶ Get the pair.
Returns: pair Return type: numpy.ndarray
, shape= \(\left(N_{particles}\right)\), dtype=numpy.uint32
-
match
¶ Match.
-
pair
¶ Pair.
-