Cluster Module

Cluster Functions

class freud.cluster.Cluster(box, rcut)

Finds clusters in a set of points.

Given a set of coordinates and a cutoff, freud.cluster.Cluster will determine all of the clusters of points that are made up of points that are closer than the cutoff. Clusters are labelled from 0 to the number of clusters-1 and an index array is returned where cluster_idx[i] is the cluster index in which particle i is found. By the definition of a cluster, points that are not within the cutoff of another point end up in their own 1-particle cluster.

Identifying micelles is one primary use-case for finding clusters. This operation is somewhat different, though. In a cluster of points, each and every point belongs to one and only one cluster. However, because a string of points belongs to a polymer, that single polymer may be present in more than one cluster. To handle this situation, an optional layer is presented on top of the cluster_idx array. Given a key value per particle (i.e. the polymer id), the computeClusterMembership function will process cluster_idx with the key values in mind and provide a list of keys that are present in each cluster.

Module author: Joshua Anderson <joaander@umich.edu>

Parameters:

Note

2D: freud.cluster.Cluster properly handles 2D boxes. The points must be passed in as [x, y, 0]. Failing to set z=0 will lead to undefined behavior.

box

Return the stored freud Box.

cluster_idx

Returns 1D array of Cluster idx for each particle.

cluster_keys

Returns the keys contained in each cluster.

computeClusterMembership(self, keys)

Compute the clusters with key membership.

Loops over all particles and adds them to a list of sets. Each set contains all the keys that are part of that cluster.

Get the computed list with getClusterKeys().

Parameters:keys (numpy.ndarray, shape=(\(N_{particles}\)), dtype= numpy.uint32) – Membership keys, one for each particle
computeClusters(self, points, nlist=None, box=None)

Compute the clusters for the given set of points.

Parameters:
getBox(self)

Return the stored freud Box.

Returns:freud Box
Return type:freud.box.Box
getClusterIdx(self)

Returns 1D array of Cluster idx for each particle

Returns:1D array of cluster idx
Return type:numpy.ndarray, shape=(\(N_{particles}\)), dtype= numpy.uint32
getClusterKeys(self)

Returns the keys contained in each cluster.

Returns:list of lists of each key contained in clusters
Return type:list
getNumClusters(self)

Returns the number of clusters.

Returns:number of clusters
Return type:int
getNumParticles(self)

Returns the number of particles.

Returns:number of particles
Return type:int
num_clusters

Returns the number of clusters.

num_particles

Returns the number of particles.

class freud.cluster.ClusterProperties(box)

Routines for computing properties of point clusters.

Given a set of points and cluster ids (from Cluster, or another source), ClusterProperties determines the following properties for each cluster:

  • Center of mass
  • Gyration tensor

The computed center of mass for each cluster (properly handling periodic boundary conditions) can be accessed with getClusterCOM(). This returns a numpy.ndarray, shape= \(\left(N_{clusters}, 3 \right)\).

The \(3 \times 3\) gyration tensor \(G\) can be accessed with getClusterG(). This returns a numpy.ndarray, shape= \(\left(N_{clusters} \times 3 \times 3\right)\). The tensor is symmetric for each cluster.

Module author: Joshua Anderson <joaander@umich.edu>

Parameters:box (freud.box.Box) – simulation box
box

Return the stored freud Box.

cluster_COM

Returns the center of mass of the last computed cluster.

cluster_G

Returns the cluster \(G\) tensors computed by the last call to computeProperties(). computeProperties.

cluster_sizes

Returns the cluster sizes computed by the last call to computeProperties(). computeProperties.

computeProperties(self, points, cluster_idx, box=None)

Compute properties of the point clusters.

Loops over all points in the given array and determines the center of mass of the cluster as well as the \(G\) tensor. These can be accessed after the call to ~.computeProperties() with getClusterCOM() and getClusterG().

Parameters:
  • points (numpy.ndarray, shape=(\(N_{particles}\), 3), dtype= numpy.float32) – Positions of the particles making up the clusters
  • cluster_idx (numpy.ndarray, shape=(\(N_{particles}\)), dtype= numpy.uint32) – List of cluster indexes for each particle
  • box (freud.box.Box) – simulation box
getBox(self)

Return the stored freud.box.Box object.

Returns:freud Box
Return type:freud.box.Box
getClusterCOM(self)

Returns the center of mass of the last computed cluster.

Returns:numpy array of cluster center of mass coordinates \(\left(x,y,z\right)\)
Return type:numpy.ndarray, shape=(\(N_{clusters}\), 3), dtype= numpy.float32
getClusterG(self)

Returns the cluster \(G\) tensors computed by the last call to computeProperties().

Returns:list of gyration tensors for each cluster
Return type:numpy.ndarray, shape=(\(N_{clusters}\), 3, 3), dtype= numpy.float32
getClusterSizes(self)

Returns the cluster sizes computed by the last call to computeProperties(). computeProperties.

Returns:sizes of each cluster
Return type:numpy.ndarray, shape=(\(N_{clusters}\)), dtype= numpy.uint32
getNumClusters(self)

Count the number of clusters found in the last call to computeProperties()

Returns:number of clusters
Return type:int
num_clusters

Returns the number of clusters.