Voronoi Module

Overview

freud.voronoi.Voronoi Compute the Voronoi tessellation of a 2D or 3D system using qhull.

Details

The voronoi module contains tools to characterize Voronoi cells of a system.

class freud.voronoi.Voronoi(box, buff)

Compute the Voronoi tessellation of a 2D or 3D system using qhull. This uses scipy.spatial.Voronoi, accounting for periodic boundary conditions.

Module author: Benjamin Schultz <baschult@umich.edu>

Module author: Yina Geng <yinageng@umich.edu>

Module author: Mayank Agrawal <amayank@umich.edu>

Module author: Bradley Dice <bdice@bradleydice.com>

Since qhull does not support periodic boundary conditions natively, we expand the box to include a portion of the particles’ periodic images. The buffer width is given by the parameter buff. The computation of Voronoi tessellations and neighbors is only guaranteed to be correct if buff >= L/2 where L is the longest side of the simulation box. For dense systems with particles filling the entire simulation volume, a smaller value for buff is acceptable.

Parameters:
compute

Compute Voronoi diagram.

Parameters:
  • positions ((\(N_{particles}\), 3) numpy.ndarray) – Points to calculate Voronoi diagram for.
  • box (freud.box.Box) – Simulation box (Default value = None).
  • buff (float) – Buffer distance within which to look for images (Default value = None).
computeNeighbors

Compute the neighbors of each particle based on the Voronoi tessellation. One can include neighbors from multiple Voronoi shells by specifying numShells in getNeighbors(). An example of computing neighbors from the first two Voronoi shells for a 2D mesh is shown below.

Retrieve the results with getNeighbors().

Example:

from freud import box, voronoi
import numpy as np
vor = voronoi.Voronoi(box.Box(5, 5, is2D=True))
pos = np.array([[0, 0, 0], [0, 1, 0], [0, 2, 0],
                [1, 0, 0], [1, 1, 0], [1, 2, 0],
                [2, 0, 0], [2, 1, 0], [2, 2, 0]], dtype=np.float32)
first_shell = vor.computeNeighbors(pos).getNeighbors(1)
second_shell = vor.computeNeighbors(pos).getNeighbors(2)
print('First shell:', first_shell)
print('Second shell:', second_shell)

Note

Input positions must be a 3D array. For 2D, set the z value to 0.

Parameters:
  • positions ((\(N_{particles}\), 3) numpy.ndarray) – Points to calculate Voronoi diagram for.
  • box (freud.box.Box) – Simulation box (Default value = None).
  • buff (float) – Buffer distance within which to look for images (Default value = None).
  • exclude_ii (bool, optional) – True if pairs of points with identical indices should be excluded (Default value = True).
computeVolumes

Computes volumes (areas in 2D) of Voronoi cells.

New in version 0.8.

Must call freud.voronoi.Voronoi.compute() before this method. Retrieve the results with freud.voronoi.Voronoi.getVolumes().

getBuffer

Returns the buffer width.

Returns:Buffer width.
Return type:float
getNeighborList

Returns a neighbor list object.

In the neighbor list, each neighbor pair has a weight value.

In 2D systems, the bond weight is the “ridge length” of the Voronoi boundary line between the neighboring particles.

In 3D systems, the bond weight is the “ridge area” of the Voronoi boundary polygon between the neighboring particles.

Returns:Neighbor list.
Return type:NeighborList
getNeighbors

Get numShells of neighbors for each particle

Must call computeNeighbors() before this method.

Parameters:numShells (int) – Number of neighbor shells.
getVolumes

Returns an array of volumes (areas in 2D) corresponding to Voronoi cells.

New in version 0.8.

Must call freud.voronoi.Voronoi.computeVolumes() before this method.

If the buffer width is too small, then some polytopes may not be closed (they may have a boundary at infinity), and these polytopes’ volumes/areas are excluded from the list.

The length of the list returned by this method should be the same as the array of positions used in the freud.voronoi.Voronoi.compute() method, if all the polytopes are closed. Otherwise try using a larger buffer width.

Returns:Voronoi polytope volumes/areas.
Return type:(\(\left(N_{cells} \right)\)) numpy.ndarray
getVoronoiPolytopes

Returns a list of polytope vertices corresponding to Voronoi cells.

If the buffer width is too small, then some polytopes may not be closed (they may have a boundary at infinity), and these polytopes’ vertices are excluded from the list.

The length of the list returned by this method should be the same as the array of positions used in the freud.voronoi.Voronoi.compute() method, if all the polytopes are closed. Otherwise try using a larger buffer width.

Returns:List of numpy.ndarray containing Voronoi polytope vertices.
Return type:list
setBox

Reset the simulation box.

Parameters:box (freud.box.Box) – Simulation box.
setBufferWidth

Reset the buffer width.

Parameters:buff (float) – Buffer width.