branches/exact/src/v_compute.hh

   1 // Voro++, a 3D cell-based Voronoi library
   2 //
   3 // Author   : Chris H. Rycroft (LBL / UC Berkeley)
   4 // Email    : chr@alum.mit.edu
   5 // Date     : August 30th 2011
   6
   7 /** \file v_compute.hh
   8  * \brief Header file for the voro_compute template and related classes. */
   9
  10 #ifndef VOROPP_V_COMPUTE_HH
  11 #define VOROPP_V_COMPUTE_HH
  12
  13 #include <cstdio>
  14 #include <cstdlib>
  15 #include <cmath>
  16 #include <vector>
  17 using namespace std;
  18
  19 #include "config.hh"
  20 #include "worklist.hh"
  21 #include "cell.hh"
  22
  23 namespace voro {
  24
  25 /** \brief Structure for holding information about a particle.
  26  *
  27  * This small structure holds information about a single particle, and is used
  28  * by several of the routines in the voro_compute template for passing
  29  * information by reference between functions. */
  30 struct particle_record {
  31         /** The index of the block that the particle is within. */
  32         int ijk;
  33         /** The number of particle within its block. */
  34         int l;
  35         /** The x-index of the block. */
  36         int di;
  37         /** The y-index of the block. */
  38         int dj;
  39         /** The z-index of the block. */
  40         int dk;
  41 };
  42
  43 /** \brief Template for carrying out Voronoi cell computations. */
  44 template <class c_class>
  45 class voro_compute {
  46         public:
  47                 /** A reference to the container class on which to carry out*/
  48                 c_class &con;
  49                 /** The size of an internal computational block in the x
  50                  * direction. */
  51                 const double boxx;
  52                 /** The size of an internal computational block in the y
  53                  * direction. */
  54                 const double boxy;
  55                 /** The size of an internal computational block in the z
  56                  * direction. */
  57                 const double boxz;
  58                 /** The inverse box length in the x direction, set to
  59                  * nx/(bx-ax). */
  60                 const double xsp;
  61                 /** The inverse box length in the y direction, set to
  62                  * ny/(by-ay). */
  63                 const double ysp;
  64                 /** The inverse box length in the z direction, set to
  65                  * nz/(bz-az). */
  66                 const double zsp;
  67                 /** The number of boxes in the x direction for the searching mask. */
  68                 const int hx;
  69                 /** The number of boxes in the y direction for the searching mask. */
  70                 const int hy;
  71                 /** The number of boxes in the z direction for the searching mask. */
  72                 const int hz;
  73                 /** A constant, set to the value of hx multiplied by hy, which
  74                  * is used in the routines which step through mask boxes in
  75                  * sequence. */
  76                 const int hxy;
  77                 /** A constant, set to the value of hx*hy*hz, which is used in
  78                  * the routines which step through mask boxes in sequence. */
  79                 const int hxyz;
  80                 /** The number of floating point entries to store for each
  81                  * particle. */
  82                 const int ps;
  83                 /** This array holds the numerical IDs of each particle in each
  84                  * computational box. */
  85                 int **id;
  86                 /** A two dimensional array holding particle positions. For the
  87                  * derived container_poly class, this also holds particle
  88                  * radii. */
  89                 double **p;
  90                 /** An array holding the number of particles within each
  91                  * computational box of the container. */
  92                 int *co;
  93                 voro_compute(c_class &con_,int hx_,int hy_,int hz_);
  94                 /** The class destructor frees the dynamically allocated memory
  95                  * for the mask and queue. */
  96                 ~voro_compute() {
  97                         delete [] qu;
  98                         delete [] mask;
  99                 }
 100                 template<class v_cell>
 101                 bool compute_cell(v_cell &c,int ijk,int s,int ci,int cj,int ck);
 102                 void find_voronoi_cell(double x,double y,double z,int ci,int cj,int ck,int ijk,particle_record &w,double &mrs);
 103         private:
 104                 /** A constant set to boxx*boxx+boxy*boxy+boxz*boxz, which is
 105                  * frequently used in the computation. */
 106                 const double bxsq;
 107                 /** This sets the current value being used to mark tested blocks
 108                  * in the mask. */
 109                 unsigned int mv;
 110                 /** The current size of the search list. */
 111                 int qu_size;
 112                 /** A pointer to the array of worklists. */
 113                 const unsigned int *wl;
 114                 /** An pointer to the array holding the minimum distances
 115                  * associated with the worklists. */
 116                 double *mrad;
 117                 /** This array is used during the cell computation to determine
 118                  * which blocks have been considered. */
 119                 unsigned int *mask;
 120                 /** An array is used to store the queue of blocks to test
 121                  * during the Voronoi cell computation. */
 122                 int *qu;
 123                 /** A pointer to the end of the queue array, used to determine
 124                  * when the queue is full. */
 125                 int *qu_l;
 126                 template<class v_cell>
 127                 bool corner_test(v_cell &c,double xl,double yl,double zl,double xh,double yh,double zh);
 128                 template<class v_cell>
 129                 inline bool edge_x_test(v_cell &c,double x0,double yl,double zl,double x1,double yh,double zh);
 130                 template<class v_cell>
 131                 inline bool edge_y_test(v_cell &c,double xl,double y0,double zl,double xh,double y1,double zh);
 132                 template<class v_cell>
 133                 inline bool edge_z_test(v_cell &c,double xl,double yl,double z0,double xh,double yh,double z1);
 134                 template<class v_cell>
 135                 inline bool face_x_test(v_cell &c,double xl,double y0,double z0,double y1,double z1);
 136                 template<class v_cell>
 137                 inline bool face_y_test(v_cell &c,double x0,double yl,double z0,double x1,double z1);
 138                 template<class v_cell>
 139                 inline bool face_z_test(v_cell &c,double x0,double y0,double zl,double x1,double y1);
 140                 bool compute_min_max_radius(int di,int dj,int dk,double fx,double fy,double fz,double gx,double gy,double gz,double& crs,double mrs);
 141                 bool compute_min_radius(int di,int dj,int dk,double fx,double fy,double fz,double mrs);
 142                 inline void add_to_mask(int ei,int ej,int ek,int *&qu_e);
 143                 inline void scan_bits_mask_add(unsigned int q,unsigned int *mijk,int ei,int ej,int ek,int *&qu_e);
 144                 inline void scan_all(int ijk,double x,double y,double z,int di,int dj,int dk,particle_record &w,double &mrs);
 145                 void add_list_memory(int*& qu_s,int*& qu_e);
 146                 /** Resets the mask in cases where the mask counter wraps
 147                  * around. */
 148                 inline void reset_mask() {
 149                         for(unsigned int *mp(mask);mp<mask+hxyz;mp++) *mp=0;
 150                 }
 151 };
 152
 153 }
 154
 155 #endif