Voro++
Public Member Functions | Data Fields | Protected Member Functions | Protected Attributes | Friends
voro::voronoicell_base Class Reference

A class representing a single Voronoi cell. More...

#include <cell.hh>

+ Inheritance diagram for voro::voronoicell_base:

Public Member Functions

 voronoicell_base ()
 
 ~voronoicell_base ()
 
void init_base (double xmin, double xmax, double ymin, double ymax, double zmin, double zmax)
 
void init_octahedron_base (double l)
 
void init_tetrahedron_base (double x0, double y0, double z0, double x1, double y1, double z1, double x2, double y2, double z2, double x3, double y3, double z3)
 
void translate (double x, double y, double z)
 
void draw_pov (double x, double y, double z, FILE *fp=stdout)
 
void draw_pov (double x, double y, double z, const char *filename)
 
void draw_pov_mesh (double x, double y, double z, FILE *fp=stdout)
 
void draw_pov_mesh (double x, double y, double z, const char *filename)
 
void draw_gnuplot (double x, double y, double z, FILE *fp=stdout)
 
void draw_gnuplot (double x, double y, double z, const char *filename)
 
double volume ()
 
double max_radius_squared ()
 
double total_edge_distance ()
 
double surface_area ()
 
void centroid (double &cx, double &cy, double &cz)
 
int number_of_faces ()
 
int number_of_edges ()
 
void vertex_orders (std::vector< int > &v)
 
void output_vertex_orders (FILE *fp=stdout)
 
void vertices (std::vector< double > &v)
 
void output_vertices (FILE *fp=stdout)
 
void vertices (double x, double y, double z, std::vector< double > &v)
 
void output_vertices (double x, double y, double z, FILE *fp=stdout)
 
void face_areas (std::vector< double > &v)
 
void output_face_areas (FILE *fp=stdout)
 
void face_orders (std::vector< int > &v)
 
void output_face_orders (FILE *fp=stdout)
 
void face_freq_table (std::vector< int > &v)
 
void output_face_freq_table (FILE *fp=stdout)
 
void face_vertices (std::vector< int > &v)
 
void output_face_vertices (FILE *fp=stdout)
 
void face_perimeters (std::vector< double > &v)
 
void output_face_perimeters (FILE *fp=stdout)
 
void normals (std::vector< double > &v)
 
void output_normals (FILE *fp=stdout)
 
void output_custom (const char *format, FILE *fp=stdout)
 
void output_custom (const char *format, int i, double x, double y, double z, double r, FILE *fp=stdout)
 
template<class vc_class >
bool nplane (vc_class &vc, double x, double y, double z, double rsq, int p_id)
 
bool plane_intersects (double x, double y, double z, double rsq)
 
bool plane_intersects_guess (double x, double y, double z, double rsq)
 
void construct_relations ()
 
void check_relations ()
 
void check_duplicates ()
 
void print_edges ()
 
virtual void neighbors (std::vector< int > &v)
 
virtual void output_neighbors (FILE *fp=stdout)
 
virtual void print_edges_neighbors (int i)
 
int cycle_up (int a, int p)
 
int cycle_down (int a, int p)
 

Data Fields

int current_vertices
 
int current_vertex_order
 
int current_delete_size
 
int current_delete2_size
 
int p
 
int up
 
int ** ed
 
int * nu
 
double * pts
 

Protected Member Functions

void reset_edges ()
 
template<class vc_class >
void check_memory_for_copy (vc_class &vc, voronoicell_base *vb)
 
void copy (voronoicell_base *vb)
 

Protected Attributes

int * mem
 
int * mec
 
int ** mep
 

Friends

class voronoicell
 
class voronoicell_neighbor
 

Detailed Description

This class represents a single Voronoi cell, as a collection of vertices that are connected by edges. The class contains routines for initializing the Voronoi cell to be simple shapes such as a box, tetrahedron, or octahedron. It the contains routines for recomputing the cell based on cutting it by a plane, which forms the key routine for the Voronoi cell computation. It contains numerous routine for computing statistics about the Voronoi cell, and it can output the cell in several formats.

This class is not intended for direct use, but forms the base of the voronoicell and voronoicell_neighbor classes, which extend it based on whether neighboring particle ID information needs to be tracked.

Definition at line 33 of file cell.hh.

Constructor & Destructor Documentation

voro::voronoicell_base::voronoicell_base ( )

Constructs a Voronoi cell and sets up the initial memory.

Definition at line 20 of file cell.cc.

voro::voronoicell_base::~voronoicell_base ( )

The voronoicell destructor deallocates all the dynamic memory.

Definition at line 43 of file cell.cc.

Member Function Documentation

void voro::voronoicell_base::centroid ( double &  cx,
double &  cy,
double &  cz 
)

Calculates the centroid of the Voronoi cell, by decomposing the cell into tetrahedra extending outward from the zeroth vertex.

Parameters
[out](cx,cy,cz)references to floating point numbers in which to pass back the centroid vector.

Definition at line 1408 of file cell.cc.

void voro::voronoicell_base::check_duplicates ( )

This routine checks for any two vertices that are connected by more than one edge. The plane algorithm is designed so that this should not happen, so any occurrences are most likely errors. Note that the routine is O(p), so running it every time the plane routine is called will result in a significant slowdown.

Definition at line 342 of file cell.cc.

template<class vc_class >
template void voro::voronoicell_base::check_memory_for_copy ( vc_class &  vc,
voronoicell_base vb 
)
protected

Ensures that enough memory is allocated prior to carrying out a copy.

Parameters
[in]vca reference to the specialized version of the calling class.
[in]vba pointered to the class to be copied.

Definition at line 56 of file cell.cc.

void voro::voronoicell_base::check_relations ( )

Checks that the relational table of the Voronoi cell is accurate, and prints out any errors. This algorithm is O(p), so running it every time the plane routine is called will result in a significant slowdown.

Definition at line 331 of file cell.cc.

void voro::voronoicell_base::construct_relations ( )

Constructs the relational table if the edges have been specified.

Definition at line 349 of file cell.cc.

void voro::voronoicell_base::copy ( voronoicell_base vb)
protected

Copies the vertex and edge information from another class. The routine assumes that enough memory is available for the copy.

Parameters
[in]vba pointer to the class to copy.

Definition at line 65 of file cell.cc.

int voro::voronoicell_base::cycle_down ( int  a,
int  p 
)
inline

This is a simple inline function for picking out the index of the next edge clockwise from the current vertex.

Parameters
[in]athe index of an edge of the current vertex.
[in]pthe number of the vertex.
Returns
nu[p]-1 if a=0, or a-1 otherwise.

Definition at line 221 of file cell.hh.

int voro::voronoicell_base::cycle_up ( int  a,
int  p 
)
inline

This is a simple inline function for picking out the index of the next edge counterclockwise at the current vertex.

Parameters
[in]athe index of an edge of the current vertex.
[in]pthe number of the vertex.
Returns
0 if a=nu[p]-1, or a+1 otherwise.

Definition at line 215 of file cell.hh.

void voro::voronoicell_base::draw_gnuplot ( double  x,
double  y,
double  z,
FILE *  fp = stdout 
)

Outputs the edges of the Voronoi cell in gnuplot format to an output stream.

Parameters
[in](x,y,z)a displacement vector to be added to the cell's position.
[in]fpa file handle to write to.

Definition at line 1507 of file cell.cc.

void voro::voronoicell_base::draw_gnuplot ( double  x,
double  y,
double  z,
const char *  filename 
)
inline

Outputs the cell in Gnuplot format a given file.

Parameters
[in](x,y,z)a displacement to add to the cell's position.
[in]filenamethe name of the file to write to.

Definition at line 119 of file cell.hh.

void voro::voronoicell_base::draw_pov ( double  x,
double  y,
double  z,
FILE *  fp = stdout 
)

Outputs the edges of the Voronoi cell in POV-Ray format to an open file stream, displacing the cell by given vector.

Parameters
[in](x,y,z)a displacement vector to be added to the cell's position.
[in]fpa file handle to write to.

Definition at line 1487 of file cell.cc.

void voro::voronoicell_base::draw_pov ( double  x,
double  y,
double  z,
const char *  filename 
)
inline

Outputs the cell in POV-Ray format, using cylinders for edges and spheres for vertices, to a given file.

Parameters
[in](x,y,z)a displacement to add to the cell's position.
[in]filenamethe name of the file to write to.

Definition at line 98 of file cell.hh.

void voro::voronoicell_base::draw_pov_mesh ( double  x,
double  y,
double  z,
FILE *  fp = stdout 
)

Outputs the Voronoi cell in the POV mesh2 format, described in section 1.3.2.2 of the POV-Ray documentation. The mesh2 output consists of a list of vertex vectors, followed by a list of triangular faces. The routine also makes use of the optional inside_vector specification, which makes the mesh object solid, so the the POV-Ray Constructive Solid Geometry (CSG) can be applied.

Parameters
[in](x,y,z)a displacement vector to be added to the cell's position.
[in]fpa file handle to write to.

Definition at line 1542 of file cell.cc.

void voro::voronoicell_base::draw_pov_mesh ( double  x,
double  y,
double  z,
const char *  filename 
)
inline

Outputs the cell in POV-Ray format as a mesh2 object to a given file.

Parameters
[in](x,y,z)a displacement to add to the cell's position.
[in]filenamethe name of the file to write to.

Definition at line 109 of file cell.hh.

void voro::voronoicell_base::face_areas ( std::vector< double > &  v)

Calculates the areas of each face of the Voronoi cell and prints the results to an output stream.

Parameters
[out]vthe vector to store the results in.

Definition at line 1336 of file cell.cc.

void voro::voronoicell_base::face_freq_table ( std::vector< int > &  v)

Computes the number of edges that each face has and outputs a frequency table of the results.

Parameters
[out]vthe vector to store the results in.

Definition at line 1891 of file cell.cc.

void voro::voronoicell_base::face_orders ( std::vector< int > &  v)

Outputs a list of the number of edges in each face.

Parameters
[out]vthe vector to store the results in.

Definition at line 1866 of file cell.cc.

void voro::voronoicell_base::face_perimeters ( std::vector< double > &  v)

This routine returns the perimeters of each face.

Parameters
[out]vthe vector to store the results in.

Definition at line 1807 of file cell.cc.

void voro::voronoicell_base::face_vertices ( std::vector< int > &  v)

For each face, this routine outputs a bracketed sequence of numbers containing a list of all the vertices that make up that face.

Parameters
[out]vthe vector to store the results in.

Definition at line 1839 of file cell.cc.

void voro::voronoicell_base::init_base ( double  xmin,
double  xmax,
double  ymin,
double  ymax,
double  zmin,
double  zmax 
)

Initializes a Voronoi cell as a rectangular box with the given dimensions.

Parameters
[in](xmin,xmax)the minimum and maximum x coordinates.
[in](ymin,ymax)the minimum and maximum y coordinates.
[in](zmin,zmax)the minimum and maximum z coordinates.

Definition at line 257 of file cell.cc.

void voro::voronoicell_base::init_octahedron_base ( double  l)

Initializes a Voronoi cell as a regular octahedron.

Parameters
[in]lThe distance from the octahedron center to a vertex. Six vertices are initialized at (-l,0,0), (l,0,0), (0,-l,0), (0,l,0), (0,0,-l), and (0,0,l).

Definition at line 286 of file cell.cc.

void voro::voronoicell_base::init_tetrahedron_base ( double  x0,
double  y0,
double  z0,
double  x1,
double  y1,
double  z1,
double  x2,
double  y2,
double  z2,
double  x3,
double  y3,
double  z3 
)

Initializes a Voronoi cell as a tetrahedron. It assumes that the normal to the face for the first three vertices points inside.

Parameters
(x0,y0,z0)a position vector for the first vertex.
(x1,y1,z1)a position vector for the second vertex.
(x2,y2,z2)a position vector for the third vertex.
(x3,y3,z3)a position vector for the fourth vertex.

Definition at line 312 of file cell.cc.

double voro::voronoicell_base::max_radius_squared ( )

Computes the maximum radius squared of a vertex from the center of the cell. It can be used to determine when enough particles have been testing an all planes that could cut the cell have been considered.

Returns
The maximum radius squared of a vertex.

Definition at line 1454 of file cell.cc.

virtual void voro::voronoicell_base::neighbors ( std::vector< int > &  v)
inlinevirtual

Returns a list of IDs of neighboring particles corresponding to each face.

Parameters
[out]va reference to a vector in which to return the results. If no neighbor information is available, a blank vector is returned.

Reimplemented in voro::voronoicell_neighbor.

Definition at line 197 of file cell.hh.

void voro::voronoicell_base::normals ( std::vector< double > &  v)

For each face of the Voronoi cell, this routine prints the out the normal vector of the face, and scales it to the distance from the cell center to that plane.

Parameters
[out]vthe vector to store the results in.

Definition at line 1639 of file cell.cc.

template<class vc_class >
template bool voro::voronoicell_base::nplane ( vc_class &  vc,
double  x,
double  y,
double  z,
double  rsq,
int  p_id 
)

Cuts the Voronoi cell by a particle whose center is at a separation of (x,y,z) from the cell center. The value of rsq should be initially set to $x^2+y^2+z^2$.

Parameters
[in]vca reference to the specialized version of the calling class.
[in](x,y,z)the normal vector to the plane.
[in]rsqthe distance along this vector of the plane.
[in]p_idthe plane ID (for neighbor tracking only).
Returns
False if the plane cut deleted the cell entirely, true otherwise.

Definition at line 404 of file cell.cc.

int voro::voronoicell_base::number_of_edges ( )

Counts the number of edges of the Voronoi cell.

Returns
the number of edges.

Definition at line 2007 of file cell.cc.

int voro::voronoicell_base::number_of_faces ( )

Returns the number of faces of a computed Voronoi cell.

Returns
The number of faces.

Definition at line 1721 of file cell.cc.

void voro::voronoicell_base::output_custom ( const char *  format,
FILE *  fp = stdout 
)
inline

Outputs a custom string of information about the Voronoi cell to a file. It assumes the cell is at (0,0,0) and has a the default_radius associated with it.

Parameters
[in]formatthe custom format string to use.
[in]fpthe file handle to write to.

Definition at line 182 of file cell.hh.

void voro::voronoicell_base::output_custom ( const char *  format,
int  i,
double  x,
double  y,
double  z,
double  r,
FILE *  fp = stdout 
)

Outputs a custom string of information about the Voronoi cell. The string of information follows a similar style as the C printf command, and detailed information about its format is available at http://math.lbl.gov/voro++/doc/custom.html.

Parameters
[in]formatthe custom string to print.
[in]ithe ID of the particle associated with this Voronoi cell.
[in](x,y,z)the position of the particle associated with this Voronoi cell.
[in]ra radius associated with the particle.
[in]fpthe file handle to write to.

Definition at line 2023 of file cell.cc.

void voro::voronoicell_base::output_face_areas ( FILE *  fp = stdout)
inline

Outputs the areas of the faces.

Parameters
[in]fpthe file handle to write to.

Definition at line 140 of file cell.hh.

void voro::voronoicell_base::output_face_freq_table ( FILE *  fp = stdout)
inline

Outputs a

Definition at line 153 of file cell.hh.

void voro::voronoicell_base::output_face_orders ( FILE *  fp = stdout)
inline

Outputs a list of the number of sides of each face.

Parameters
[in]fpthe file handle to write to.

Definition at line 147 of file cell.hh.

void voro::voronoicell_base::output_face_perimeters ( FILE *  fp = stdout)
inline

Outputs a list of the perimeters of each face.

Parameters
[in]fpthe file handle to write to.

Definition at line 166 of file cell.hh.

void voro::voronoicell_base::output_face_vertices ( FILE *  fp = stdout)
inline

Outputs the

Definition at line 159 of file cell.hh.

virtual void voro::voronoicell_base::output_neighbors ( FILE *  fp = stdout)
inlinevirtual

This is a virtual function that is overridden by a routine to print a list of IDs of neighboring particles corresponding to each face. By default, when no neighbor information is available, the routine does nothing.

Parameters
[in]fpthe file handle to write to.

Reimplemented in voro::voronoicell_neighbor.

Definition at line 203 of file cell.hh.

void voro::voronoicell_base::output_normals ( FILE *  fp = stdout)
inline

Outputs a list of the perimeters of each face.

Parameters
[in]fpthe file handle to write to.

Definition at line 173 of file cell.hh.

void voro::voronoicell_base::output_vertex_orders ( FILE *  fp = stdout)

Outputs the vertex orders.

Parameters
[out]fpthe file handle to write to.

Definition at line 1751 of file cell.cc.

void voro::voronoicell_base::output_vertices ( FILE *  fp = stdout)

Outputs the vertex vectors using the local coordinate system.

Parameters
[out]fpthe file handle to write to.

Definition at line 1772 of file cell.cc.

void voro::voronoicell_base::output_vertices ( double  x,
double  y,
double  z,
FILE *  fp = stdout 
)

Outputs the vertex vectors using the global coordinate system.

Parameters
[out]fpthe file handle to write to.
[in](x,y,z)the position vector of the particle in the global coordinate system.

Definition at line 1798 of file cell.cc.

bool voro::voronoicell_base::plane_intersects ( double  x,
double  y,
double  z,
double  rsq 
)

This routine tests to see whether the cell intersects a plane by starting from the guess point up. If up intersects, then it immediately returns true. Otherwise, it calls the plane_intersects_track() routine.

Parameters
[in](x,y,z)the normal vector to the plane.
[in]rsqthe distance along this vector of the plane.
Returns
False if the plane does not intersect the plane, true if it does.

Definition at line 1920 of file cell.cc.

bool voro::voronoicell_base::plane_intersects_guess ( double  x,
double  y,
double  z,
double  rsq 
)

This routine tests to see if a cell intersects a plane. It first tests a random sample of approximately sqrt(p)/4 points. If any of those are intersect, then it immediately returns true. Otherwise, it takes the closest point and passes that to plane_intersect_track() routine.

Parameters
[in](x,y,z)the normal vector to the plane.
[in]rsqthe distance along this vector of the plane.
Returns
False if the plane does not intersect the plane, true if it does.

Definition at line 1933 of file cell.cc.

void voro::voronoicell_base::print_edges ( )

Prints the vertices, their edges, the relation table, and also notifies if any memory errors are visible.

Definition at line 2220 of file cell.cc.

virtual void voro::voronoicell_base::print_edges_neighbors ( int  i)
inlinevirtual

This a virtual function that is overridden by a routine to print the neighboring particle IDs for a given vertex. By default, when no neighbor information is available, the routine does nothing.

Parameters
[in]ithe vertex to consider.

Reimplemented in voro::voronoicell_neighbor.

Definition at line 209 of file cell.hh.

void voro::voronoicell_base::reset_edges ( )
inlineprotected

Several routines in the class that gather cell-based statistics internally track their progress by flipping edges to negative so that they know what parts of the cell have already been tested. This function resets them back to positive. When it is called, it assumes that every edge in the routine should have already been flipped to negative, and it bails out with an internal error if it encounters a positive edge.

Definition at line 1572 of file cell.cc.

double voro::voronoicell_base::surface_area ( )

Calculates the total surface area of the Voronoi cell.

Returns
The computed area.

Definition at line 1372 of file cell.cc.

double voro::voronoicell_base::total_edge_distance ( )

Calculates the total edge distance of the Voronoi cell.

Returns
A floating point number holding the calculated distance.

Definition at line 1468 of file cell.cc.

void voro::voronoicell_base::translate ( double  x,
double  y,
double  z 
)

Translates the vertices of the Voronoi cell by a given vector.

Parameters
[in](x,y,z)the coordinates of the vector.

Definition at line 105 of file cell.cc.

void voro::voronoicell_base::vertex_orders ( std::vector< int > &  v)

Returns a vector of the vertex orders.

Parameters
[out]vthe vector to store the results in.

Definition at line 1744 of file cell.cc.

void voro::voronoicell_base::vertices ( std::vector< double > &  v)

Returns a vector of the vertex vectors using the local coordinate system.

Parameters
[out]vthe vector to store the results in.

Definition at line 1760 of file cell.cc.

void voro::voronoicell_base::vertices ( double  x,
double  y,
double  z,
std::vector< double > &  v 
)

Returns a vector of the vertex vectors in the global coordinate system.

Parameters
[out]vthe vector to store the results in.
[in](x,y,z)the position vector of the particle in the global coordinate system.

Definition at line 1784 of file cell.cc.

double voro::voronoicell_base::volume ( )

Calculates the volume of the Voronoi cell, by decomposing the cell into tetrahedra extending outward from the zeroth vertex, whose volumes are evaluated using a scalar triple product.

Returns
A floating point number holding the calculated volume.

Definition at line 1299 of file cell.cc.

Field Documentation

int voro::voronoicell_base::current_delete2_size

This sets the size of the auxiliary delete stack.

Definition at line 49 of file cell.hh.

int voro::voronoicell_base::current_delete_size

This sets the size of the main delete stack.

Definition at line 47 of file cell.hh.

int voro::voronoicell_base::current_vertex_order

This holds the current maximum allowed order of a vertex, which sets the size of the mem, mep, and mec arrays. If a vertex is created with more vertices than this, the arrays are dynamically extended using the add_memory_vorder routine.

Definition at line 45 of file cell.hh.

int voro::voronoicell_base::current_vertices

This holds the current size of the arrays ed and nu, which hold the vertex information. If more vertices are created than can fit in this array, then it is dynamically extended using the add_memory_vertices routine.

Definition at line 39 of file cell.hh.

int** voro::voronoicell_base::ed

This is a two dimensional array that holds information about the edge connections of the vertices that make up the cell. The two dimensional array is not allocated in the usual method. To account for the fact the different vertices have different orders, and thus require different amounts of storage, the elements of ed[i] point to one-dimensional arrays in the mep[] array of different sizes.

More specifically, if vertex i has order m, then ed[i] points to a one-dimensional array in mep[m] that has 2*m+1 entries. The first m elements hold the neighboring edges, so that the jth edge of vertex i is held in ed[i][j]. The next m elements hold a table of relations which is redundant but helps speed up the computation. It satisfies the relation ed[ed[i][j]][ed[i][m+j]]=i. The final entry holds a back pointer, so that ed[i+2*m]=i. The back pointers are used when rearranging the memory.

Definition at line 78 of file cell.hh.

int* voro::voronoicell_base::mec
protected

This is a one dimensional array that holds the current number of vertices of order p that are stored in the mep[p] array.

Definition at line 229 of file cell.hh.

int* voro::voronoicell_base::mem
protected

This a one dimensional array that holds the current sizes of the memory allocations for them mep array.

Definition at line 225 of file cell.hh.

int** voro::voronoicell_base::mep
protected

This is a two dimensional array for holding the information about the edges of the Voronoi cell. mep[p] is a one-dimensional array for holding the edge information about all vertices of order p, with each vertex holding 2*p+1 integers of information. The total number of vertices held on mep[p] is stored in mem[p]. If the space runs out, the code allocates more using the add_memory() routine.

Definition at line 237 of file cell.hh.

int* voro::voronoicell_base::nu

This array holds the order of the vertices in the Voronoi cell. This array is dynamically allocated, with its current size held by current_vertices.

Definition at line 82 of file cell.hh.

int voro::voronoicell_base::p

This sets the total number of vertices in the current cell.

Definition at line 52 of file cell.hh.

double* voro::voronoicell_base::pts

This in an array with size 3*current_vertices for holding the positions of the vertices.

Definition at line 85 of file cell.hh.

int voro::voronoicell_base::up

This is the index of particular point in the cell, which is used to start the tracing routines for plane intersection and cutting. These routines will work starting from any point, but it's often most efficient to start from the last point considered, since in many cases, the cell construction algorithm may consider many planes with similar vectors concurrently.

Definition at line 60 of file cell.hh.


The documentation for this class was generated from the following files: