Voronoi Module¶
Overview
Compute the Voronoi tessellation of a 2D or 3D system using qhull. |
Details
The freud.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 ifbuff >= L/2
whereL
is the longest side of the simulation box. For dense systems with particles filling the entire simulation volume, a smaller value forbuff
is acceptable. 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. If either the polytopes or volumes lists that are computed is different from the size of the array of positions used in thefreud.voronoi.Voronoi.compute()
method, try recomputing using a larger buffer width.- Parameters
box (
freud.box.Box
) – Simulation box.buff (float) – Buffer width.
- Variables
buffer (float) – Buffer width.
nlist (
NeighborList
) – Returns a weighted neighbor list. 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.polytopes (list[
numpy.ndarray
]) – List of arrays, each containing Voronoi polytope vertices.volumes ((\(\left(N_{cells} \right)\))
numpy.ndarray
) – Returns an array of volumes (areas in 2D) corresponding to Voronoi cells.
-
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
ingetNeighbors()
. 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 the volumes attribute.
-
getNeighbors
¶ Get well-sorted neighbors from cumulative Voronoi shells for each particle by specifying
numShells
.Must call
computeNeighbors()
before this method.- Parameters
numShells (int) – Number of neighbor shells.
-
plot
¶ Plot Voronoi diagram.
- Parameters
ax (
matplotlib.axes.Axes
) – Axis to plot on. IfNone
, make a new figure and axis. (Default value =None
)- Returns
Axis with the plot.
- Return type