Order Module¶
Overview
freud.order.CubaticOrderParameter |
Compute the cubatic order parameter [HajiAkbari2015] for a system of particles using simulated annealing instead of Newton-Raphson root finding. |
freud.order.NematicOrderParameter |
Compute the nematic order parameter for a system of particles. |
freud.order.HexOrderParameter |
Calculates the \(k\)-atic order parameter for each particle in the system. |
freud.order.TransOrderParameter |
Compute the translational order parameter for each particle. |
freud.order.LocalQl |
Compute the local Steinhardt [Steinhardt1983] rotationally invariant \(Q_l\) order parameter for a set of points. |
freud.order.LocalQlNear |
A variant of the LocalQl class that performs its average over nearest neighbor particles as determined by an instance of freud.locality.NeighborList . |
freud.order.LocalWl |
Compute the local Steinhardt [Steinhardt1983] rotationally invariant \(W_l\) order parameter for a set of points. |
freud.order.LocalWlNear |
A variant of the LocalWl class that performs its average over nearest neighbor particles as determined by an instance of freud.locality.NeighborList . |
freud.order.SolLiq |
Computes dot products of \(Q_{lm}\) between particles and uses these for clustering. |
freud.order.SolLiqNear |
A variant of the SolLiq class that performs its average over nearest neighbor particles as determined by an instance of freud.locality.NeighborList . |
freud.order.BondOrder |
Deprecated Compute the bond order diagram for the system of particles. |
freud.order.LocalDescriptors |
Deprecated Compute a set of descriptors (a numerical “fingerprint”) of a particle’s local environment. |
freud.order.MatchEnv |
Deprecated Clusters particles according to whether their local environments match or not, according to various shape matching metrics. |
freud.order.Pairing2D |
Deprecated Compute pairs for the system of particles. |
freud.order.AngularSeparation |
Deprecated Calculates the minimum angles of separation between particles and references. |
Details
The 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.
CubaticOrderParameter
(t_initial, t_final, scale, n_replicates, seed)¶ Compute the cubatic order parameter [HajiAkbari2015] for a system of particles using simulated annealing instead of Newton-Raphson root finding.
Module author: Eric Harper <harperic@umich.edu>
Parameters: Variables: - t_initial (float) – The value of the initial temperature.
- t_final (float) – The value of the final temperature.
- scale (float) – The scale
- cubatic_order_parameter (float) – The cubatic order parameter.
- orientation (\(\left(4 \right)\)
numpy.ndarray
) – The quaternion of global orientation. - particle_order_parameter (
numpy.ndarray
) – Cubatic order parameter. - particle_tensor (\(\left(N_{particles}, 3, 3, 3, 3 \right)\)
numpy.ndarray
) – Rank 5 tensor corresponding to each individual particle orientation. - global_tensor (\(\left(3, 3, 3, 3 \right)\)
numpy.ndarray
) – Rank 4 tensor corresponding to global orientation. - cubatic_tensor (\(\left(3, 3, 3, 3 \right)\)
numpy.ndarray
) – Rank 4 cubatic tensor. - gen_r4_tensor (\(\left(3, 3, 3, 3 \right)\)
numpy.ndarray
) – Rank 4 tensor corresponding to each individual particle orientation.
-
compute
¶ Calculates the per-particle and global order parameter.
Parameters: orientations ((\(N_{particles}\), 4) numpy.ndarray
) – Orientations as angles to use in computation.
-
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 (\(\left(3 \right)\)
numpy.ndarray
) – The nematic director of a single particle in the reference state (without any rotation applied).Variables: - nematic_order_parameter (float) – Nematic order parameter.
- director (\(\left(3 \right)\)
numpy.ndarray
) – The average nematic director. - particle_tensor (\(\left(N_{particles}, 3, 3 \right)\)
numpy.ndarray
) – One 3x3 matrix per-particle corresponding to each individual particle orientation. - nematic_tensor (\(\left(3, 3 \right)\)
numpy.ndarray
) – 3x3 matrix corresponding to the average particle orientation.
-
compute
¶ Calculates the per-particle and global order parameter.
Parameters: orientations (\(\left(N_{particles}, 4 \right)\) numpy.ndarray
) – Orientations to calculate the 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:
freud.order.HexOrderParameter
properly handles 2D boxes. The points must be passed in as[x, y, 0]
. Failing to set z=0 will lead to undefined behavior.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).
Variables: - psi (\(\left(N_{particles} \right)\)
numpy.ndarray
) – Order parameter. - box (
freud.box.Box
) – Box used in the calculation. - num_particles (unsigned int) – Number of particles.
- K (unsigned int) – Symmetry of the order parameter.
-
compute
¶ Calculates the correlation function and adds to the current histogram.
Parameters: - box (
freud.box.Box
) – Simulation box. - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
) – Neighborlist to use to find bonds.
- box (
-
class
freud.order.
TransOrderParameter
(rmax, k, n)¶ Compute the translational order parameter for each particle.
Module author: Michael Engel <engelmm@umich.edu>
Parameters: Variables: - d_r (\(\left(N_{particles}\right)\)
numpy.ndarray
) – Reference to the last computed translational order array. - box (
freud.box.Box
) – Box used in the calculation. - num_particles (unsigned int) – Number of particles.
-
compute
¶ Calculates the local descriptors.
Parameters: - box (
freud.box.Box
) – Simulation box. - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
) – Neighborlist to use to find bonds.
- box (
- d_r (\(\left(N_{particles}\right)\)
-
class
freud.order.
LocalQl
(box, rmax, l, rmin)¶ Compute the local Steinhardt [Steinhardt1983] rotationally invariant \(Q_l\) 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}))\). The particles included in the sum are determined by the rmax argument to the constructor.
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 }\).
The
computeAve()
method provides access to a variant of this parameter that performs a average over the first and second shell combined [Lechner2008]. 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
computeNorm()
andcomputeAveNorm()
methods provide normalized versions ofcompute()
andcomputeAve()
, where the normalization is performed by dividing by the average \(Q_{lm}\) values over all particles.Module author: Xiyu Du <xiyudu@umich.edu>
Module author: Vyas Ramasubramani <vramasub@umich.edu>
Parameters: - box (
freud.box.Box
) – Simulation box. - rmax (float) – Cutoff radius for the local order parameter. Values near the first minimum 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.
Variables: - box (
freud.box.Box
) – Box used in the calculation. - num_particles (unsigned int) – Number of particles.
- Ql (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(Q_l\) for each particle (filled with NaN for particles with no neighbors). - ave_Ql (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(\bar{Q_l}\) for each particle (filled with NaN for particles with no neighbors). - norm_Ql (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(Q_l\) for each particle normalized by the value over all particles (filled with NaN for particles with no neighbors). - ave_norm_Ql (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(\bar{Q_l}\) for each particle normalized by the value over all particles (filled with NaN for particles with no neighbors).
-
compute
¶ Compute the order parameter.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeAve
¶ Compute the order parameter over two nearest neighbor shells.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeAveNorm
¶ Compute the order parameter over two nearest neighbor shells normalized by the average spherical harmonic value over all the particles.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeNorm
¶ Compute the order parameter normalized by the average spherical harmonic value over all the particles.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
setBox
¶ Reset the simulation box.
Parameters: box ( freud.box.Box
) – Simulation box.
- box (
-
class
freud.order.
LocalQlNear
(box, rmax, l, kn)¶ A variant of the
LocalQl
class that performs its average over nearest neighbor particles as determined by an instance offreud.locality.NeighborList
. The number of included neighbors is determined by the kn parameter to the constructor.Module author: Xiyu Du <xiyudu@umich.edu>
Module author: Vyas Ramasubramani <vramasub@umich.edu>
Parameters: - box (
freud.box.Box
) – Simulation box. - rmax (float) – Cutoff radius for the local order parameter. Values near the first minimum 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.
Variables: - box (
freud.box.Box
) – Box used in the calculation. - num_particles (unsigned int) – Number of particles.
- Ql (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(Q_l\) for each particle (filled with NaN for particles with no neighbors). - ave_Ql (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(\bar{Q_l}\) for each particle (filled with NaN for particles with no neighbors). - norm_Ql (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(Q_l\) for each particle normalized by the value over all particles (filled with NaN for particles with no neighbors). - ave_norm_Ql (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(\bar{Q_l}\) for each particle normalized by the value over all particles (filled with NaN for particles with no neighbors).
-
compute
¶ Compute the order parameter.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeAve
¶ Compute the order parameter over two nearest neighbor shells.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeAveNorm
¶ Compute the order parameter over two nearest neighbor shells normalized by the average spherical harmonic value over all the particles.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeNorm
¶ Compute the order parameter normalized by the average spherical harmonic value over all the particles.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
setBox
¶ Reset the simulation box.
Parameters: box ( freud.box.Box
) – Simulation box.
- box (
-
class
freud.order.
LocalWl
(box, rmax, l)¶ Compute the local Steinhardt [Steinhardt1983] rotationally invariant \(W_l\) order parameter for a set of points.
Implements the local rotationally invariant \(W_l\) order parameter described by Steinhardt. For a particle i, we calculate the average \(W_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}))\). The particles included in the sum are determined by the rmax argument to the constructor.
The \(W_l\) is then defined as a weighted average over the \(\overline{Q}_{lm}(i)\) values using Wigner 3j symbols (CG coefficients). The resulting combination is rotationally (i.e. frame) invariant.
The
computeAve()
method provides access to a variant of this parameter that performs a average over the first and second shell combined [Lechner2008]. 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
computeNorm()
andcomputeAveNorm()
methods provide normalized versions ofcompute()
andcomputeAve()
, where the normalization is performed by dividing by the average \(Q_{lm}\) values over all particles.Module author: Xiyu Du <xiyudu@umich.edu>
Module author: Vyas Ramasubramani <vramasub@umich.edu>
Parameters: - box (
freud.box.Box
) – Simulation box. - rmax (float) – Cutoff radius for the local order parameter. Values near the first minimum of the RDF are recommended.
- l (unsigned int) – Spherical harmonic quantum number l. Must be a positive number
- rmin (float) – Lower bound for computing the local order parameter. Allows looking at, for instance, only the second shell, or some other arbitrary RDF region.
Variables: - box (
freud.box.Box
) – Box used in the calculation. - num_particles (unsigned int) – Number of particles.
- Wl (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(W_l\) for each particle (filled with NaN for particles with no neighbors). - ave_Wl (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(\bar{W}_l\) for each particle (filled with NaN for particles with no neighbors). - norm_Wl (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(W_l\) for each particle normalized by the value over all particles (filled with NaN for particles with no neighbors). - ave_norm_Wl (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(\bar{W}_l\) for each particle normalized by the value over all particles (filled with NaN for particles with no neighbors).
-
compute
¶ Compute the order parameter.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeAve
¶ Compute the order parameter over two nearest neighbor shells.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeAveNorm
¶ Compute the order parameter over two nearest neighbor shells normalized by the average spherical harmonic value over all the particles.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeNorm
¶ Compute the order parameter normalized by the average spherical harmonic value over all the particles.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
setBox
¶ Reset the simulation box.
Parameters: box ( freud.box.Box
) – Simulation box.
- box (
-
class
freud.order.
LocalWlNear
(box, rmax, l, kn)¶ A variant of the
LocalWl
class that performs its average over nearest neighbor particles as determined by an instance offreud.locality.NeighborList
. The number of included neighbors is determined by the kn parameter to the constructor.Module author: Xiyu Du <xiyudu@umich.edu>
Module author: Vyas Ramasubramani <vramasub@umich.edu>
Parameters: - box (
freud.box.Box
) – Simulation box. - rmax (float) – Cutoff radius for the local order parameter. Values near the first minimum 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.
Variables: - box (
freud.box.Box
) – Box used in the calculation. - num_particles (unsigned int) – Number of particles.
- Wl (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(W_l\) for each particle (filled with NaN for particles with no neighbors). - ave_Wl (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(\bar{W}_l\) for each particle (filled with NaN for particles with no neighbors). - norm_Wl (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(W_l\) for each particle normalized by the value over all particles (filled with NaN for particles with no neighbors). - ave_norm_Wl (\(\left(N_{particles}\right)\)
numpy.ndarray
) – The last computed \(\bar{W}_l\) for each particle normalized by the value over all particles (filled with NaN for particles with no neighbors).
-
compute
¶ Compute the order parameter.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeAve
¶ Compute the order parameter over two nearest neighbor shells.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeAveNorm
¶ Compute the order parameter over two nearest neighbor shells normalized by the average spherical harmonic value over all the particles.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeNorm
¶ Compute the order parameter normalized by the average spherical harmonic value over all the particles.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
setBox
¶ Reset the simulation box.
Parameters: box ( freud.box.Box
) – Simulation box.
- box (
-
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>
- Args:
- box (
freud.box.Box
): - Simulation box.
- rmax (float):
- Cutoff radius for the local order parameter. Values near first minimum 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 is generally good for FCC or BCC structures).
- l (unsigned int):
- Choose spherical harmonic \(Q_l\). Must be positive and even.
- box (
- Attributes:
- box (
freud.box.Box
): - Box used in the calculation.
- largest_cluster_size (unsigned int):
- The largest cluster size. Must call a compute method first.
- cluster_sizes (unsigned int):
- The sizes of all clusters.
- largest_cluster_size (unsigned int):
- The largest cluster size. Must call a compute method first.
- Ql_mi (\(\left(N_{particles}\right)\)
numpy.ndarray
): - The last computed \(Q_{lmi}\) for each particle.
- clusters (\(\left(N_{particles}\right)\)
numpy.ndarray
): - The last computed set of solid-like cluster indices for each particle.
- num_connections (\(\left(N_{particles}\right)\)
numpy.ndarray
): - The number of connections per particle.
- Ql_dot_ij (\(\left(N_{particles}\right)\)
numpy.ndarray
): - Reference to the qldot_ij values.
- num_particles (unsigned int):
- Number of particles.
- box (
-
compute
¶ Compute the solid-liquid order parameter.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeSolLiqNoNorm
¶ Compute the solid-liquid order parameter without normalizing the dot product.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
computeSolLiqVariant
¶ Compute a variant of the solid-liquid order parameter.
This variant method places a minimum threshold on the number of solid-like bonds a particle must have to be considered solid-like for clustering purposes.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
, optional) – Neighborlist to use to find bonds (Default value = None).
- points ((\(N_{particles}\), 3)
-
class
freud.order.
SolLiqNear
(box, rmax, Qthreshold, Sthreshold, l)¶ SolLiqNear(box, rmax, Qthreshold, Sthreshold, l, kn=12) A variant of the
SolLiq
class that performs its averageover nearest neighbor particles as determined by an instance of
freud.locality.NeighborList
. The number of included neighbors is determined by the kn parameter to the constructor.Module author: Richmond Newman <newmanrs@umich.edu>
- Args:
- box (
freud.box.Box
): - Simulation box.
- rmax (float):
- Cutoff radius for the local order parameter. Values near the first minimum 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 is 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.
- box (
- Attributes:
- box (
freud.box.Box
): - Box used in the calculation.
- largest_cluster_size (unsigned int):
- The largest cluster size. Must call a compute method first.
- cluster_sizes (unsigned int):
- The sizes of all clusters.
- largest_cluster_size (unsigned int):
- The largest cluster size. Must call a compute method first.
- Ql_mi (\(\left(N_{particles}\right)\)
numpy.ndarray
): - The last computed \(Q_{lmi}\) for each particle.
- clusters (\(\left(N_{particles}\right)\)
numpy.ndarray
): - The last computed set of solid-like cluster indices for each particle.
- num_connections (\(\left(N_{particles}\right)\)
numpy.ndarray
): - The number of connections per particle.
- Ql_dot_ij (\(\left(N_{particles}\right)\)
numpy.ndarray
): - Reference to the qldot_ij values.
- num_particles (unsigned int):
- Number of particles.
- box (
-
compute
¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
) – Neighborlist to use to find bonds.
- points ((\(N_{particles}\), 3)
-
computeSolLiqNoNorm
¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
) – Neighborlist to use to find bonds.
- points ((\(N_{particles}\), 3)
-
computeSolLiqVariant
¶ Compute the local rotationally invariant \(Q_l\) order parameter.
Parameters: - points ((\(N_{particles}\), 3)
numpy.ndarray
) – Points to calculate the order parameter. - nlist (
freud.locality.NeighborList
) – Neighborlist to use to find bonds.
- points ((\(N_{particles}\), 3)
Deprecated Classes¶
The below functions have all either been deprecated or moved to the Environment Module module
Bond Order¶
-
class
freud.order.
BondOrder
(rmax, k, n, nBinsT, nBinsP)¶ Deprecated Compute the bond order diagram for the system of particles.
Note
This class is only retained for backwards compatibility. Please use
freud.environment.BondOrder
instead.Deprecated since version 0.8.2: Use
freud.environment.BondOrder
instead.
Local Descriptors¶
-
class
freud.order.
LocalDescriptors
(box, nNeigh, lmax, rmax)¶ Deprecated Compute a set of descriptors (a numerical “fingerprint”) of a particle’s local environment.
Note
This class is only retained for backwards compatibility. Please use
freud.environment.LocalDescriptors
instead.Deprecated since version 0.8.2: Use
freud.environment.LocalDescriptors
instead.
Environment Matching¶
-
class
freud.order.
MatchEnv
(box, rmax, k)¶ Deprecated Clusters particles according to whether their local environments match or not, according to various shape matching metrics.
Note
This class is only retained for backwards compatibility. Please use
freud.environment.MatchEnv
instead.Deprecated since version 0.8.2: Use
freud.environment.MatchEnv
instead.
Pairing¶
Note
Pairing2D is deprecated and is replaced with Bond Module.
-
class
freud.order.
Pairing2D
(rmax, k, compDotTol)¶ Deprecated Compute pairs for the system of particles. .. note:
This class is only retained for backwards compatibility. Please use :py:mod:`freud.bond` instead.
Deprecated since version 0.8.2: Use
freud.bond
instead.
Angular Separation¶
-
class
freud.order.
AngularSeparation
(box, rmax, n)¶ Deprecated Calculates the minimum angles of separation between particles and references.
Note
This class is only retained for backwards compatibility. Please use
freud.environment.AngularSeparation
instead.Deprecated since version 0.8.2: Use
freud.environment.AngularSeparation
instead.