Voro++
cell.hh
Go to the documentation of this file.
00001 // Voro++, a 3D cell-based Voronoi library
00002 //
00003 // Author   : Chris H. Rycroft (LBL / UC Berkeley)
00004 // Email    : chr@alum.mit.edu
00005 // Date     : August 30th 2011
00006 
00007 /** \file cell.hh
00008  * \brief Header file for the voronoicell and related classes. */
00009 
00010 #ifndef VOROPP_CELL_HH
00011 #define VOROPP_CELL_HH
00012 
00013 #include <cstdio>
00014 #include <cstdlib>
00015 #include <cmath>
00016 #include <vector>
00017 using namespace std;
00018 
00019 #include "config.hh"
00020 #include "common.hh"
00021 
00022 namespace voro {
00023 
00024 /** \brief A class representing a single Voronoi cell.
00025  *
00026  * This class represents a single Voronoi cell, as a collection of vertices
00027  * that are connected by edges. The class contains routines for initializing
00028  * the Voronoi cell to be simple shapes such as a box, tetrahedron, or octahedron.
00029  * It the contains routines for recomputing the cell based on cutting it
00030  * by a plane, which forms the key routine for the Voronoi cell computation.
00031  * It contains numerous routine for computing statistics about the Voronoi cell,
00032  * and it can output the cell in several formats.
00033  *
00034  * This class is not intended for direct use, but forms the base of the
00035  * voronoicell and voronoicell_neighbor classes, which extend it based on
00036  * whether neighboring particle ID information needs to be tracked. */
00037 class voronoicell_base {
00038         public:
00039                 /** This holds the current size of the arrays ed and nu, which
00040                  * hold the vertex information. If more vertices are created
00041                  * than can fit in this array, then it is dynamically extended
00042                  * using the add_memory_vertices routine. */
00043                 int current_vertices;
00044                 /** This holds the current maximum allowed order of a vertex,
00045                  * which sets the size of the mem, mep, and mec arrays. If a
00046                  * vertex is created with more vertices than this, the arrays
00047                  * are dynamically extended using the add_memory_vorder routine.
00048                  */
00049                 int current_vertex_order;
00050                 /** This sets the size of the main delete stack. */
00051                 int current_delete_size;
00052                 /** This sets the size of the auxiliary delete stack. */
00053                 int current_delete2_size;
00054                 /** This sets the total number of vertices in the current cell.
00055                  */
00056                 int p;
00057                 /** This is the index of particular point in the cell, which is
00058                  * used to start the tracing routines for plane intersection
00059                  * and cutting. These routines will work starting from any
00060                  * point, but it's often most efficient to start from the last
00061                  * point considered, since in many cases, the cell construction
00062                  * algorithm may consider many planes with similar vectors
00063                  * concurrently. */
00064                 int up;
00065                 /** This is a two dimensional array that holds information
00066                  * about the edge connections of the vertices that make up the
00067                  * cell. The two dimensional array is not allocated in the
00068                  * usual method. To account for the fact the different vertices
00069                  * have different orders, and thus require different amounts of
00070                  * storage, the elements of ed[i] point to one-dimensional
00071                  * arrays in the mep[] array of different sizes.
00072                  *
00073                  * More specifically, if vertex i has order m, then ed[i]
00074                  * points to a one-dimensional array in mep[m] that has 2*m+1
00075                  * entries. The first m elements hold the neighboring edges, so
00076                  * that the jth edge of vertex i is held in ed[i][j]. The next
00077                  * m elements hold a table of relations which is redundant but
00078                  * helps speed up the computation. It satisfies the relation
00079                  * ed[ed[i][j]][ed[i][m+j]]=i. The final entry holds a back
00080                  * pointer, so that ed[i+2*m]=i. The back pointers are used
00081                  * when rearranging the memory. */
00082                 int **ed;
00083                 /** This array holds the order of the vertices in the Voronoi
00084                  * cell. This array is dynamically allocated, with its current
00085                  * size held by current_vertices. */
00086                 int *nu;
00087                 /** This in an array with size 3*current_vertices for holding
00088                  * the positions of the vertices. */
00089                 double *pts;
00090                 voronoicell_base();
00091                 ~voronoicell_base();
00092                 void init_base(double xmin,double xmax,double ymin,double ymax,double zmin,double zmax);
00093                 void init_octahedron_base(double l);
00094                 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);
00095                 void translate(double x,double y,double z);
00096                 void draw_pov(double x,double y,double z,FILE *fp=stdout);
00097                 /** Outputs the cell in POV-Ray format, using cylinders for edges
00098                  * and spheres for vertices, to a given file.
00099                  * \param[in] (x,y,z) a displacement to add to the cell's
00100                  *                    position.
00101                  * \param[in] filename the name of the file to write to. */
00102                 inline void draw_pov(double x,double y,double z,const char *filename) {
00103                         FILE *fp(fopen(filename,"w"));
00104                         if(fp==NULL) voro_fatal_error("Unable to open file",VOROPP_FILE_ERROR);
00105                         draw_pov(x,y,z,fp);
00106                         fclose(fp);
00107                 };
00108                 void draw_pov_mesh(double x,double y,double z,FILE *fp=stdout);
00109                 /** Outputs the cell in POV-Ray format as a mesh2 object to a
00110                  * given file.
00111                  * \param[in] (x,y,z) a displacement to add to the cell's
00112                  *                    position.
00113                  * \param[in] filename the name of the file to write to. */
00114                 inline void draw_pov_mesh(double x,double y,double z,const char *filename) {
00115                         FILE *fp(fopen(filename,"w"));
00116                         if(fp==NULL) voro_fatal_error("Unable to open file",VOROPP_FILE_ERROR);
00117                         draw_pov_mesh(x,y,z,fp);
00118                         fclose(fp);
00119                 }
00120                 void draw_gnuplot(double x,double y,double z,FILE *fp=stdout);
00121                 /** Outputs the cell in Gnuplot format a given file.
00122                  * \param[in] (x,y,z) a displacement to add to the cell's
00123                  *                    position.
00124                  * \param[in] filename the name of the file to write to. */
00125                 inline void draw_gnuplot(double x,double y,double z,const char *filename) {
00126                         FILE *fp(fopen(filename,"w"));
00127                         if(fp==NULL) voro_fatal_error("Unable to open file",VOROPP_FILE_ERROR);
00128                         draw_gnuplot(x,y,z,fp);
00129                         fclose(fp);
00130                 }
00131                 double volume();
00132                 double max_radius_squared();
00133                 double total_edge_distance();
00134                 double surface_area();
00135                 void centroid(double &cx,double &cy,double &cz);
00136                 int number_of_faces();
00137                 int number_of_edges();
00138                 void vertex_orders(vector<int> &v);
00139                 void output_vertex_orders(FILE *fp=stdout);
00140                 void vertices(vector<double> &v);
00141                 void output_vertices(FILE *fp=stdout);
00142                 void vertices(double x,double y,double z,vector<double> &v);
00143                 void output_vertices(double x,double y,double z,FILE *fp=stdout);
00144                 void face_areas(vector<double> &v);
00145                 /** Outputs the areas of the faces.
00146                  * \param[in] fp the file handle to write to. */
00147                 inline void output_face_areas(FILE *fp=stdout) {
00148                         vector<double> v;face_areas(v);
00149                         voro_print_vector(v,fp);
00150                 }
00151                 void face_orders(vector<int> &v);
00152                 /** Outputs a list of the number of sides of each face.
00153                  * \param[in] fp the file handle to write to. */
00154                 inline void output_face_orders(FILE *fp=stdout) {
00155                         vector<int> v;face_orders(v);
00156                         voro_print_vector(v,fp);
00157                 }
00158                 void face_freq_table(vector<int> &v);
00159                 /** Outputs a */
00160                 inline void output_face_freq_table(FILE *fp=stdout) {
00161                         vector<int> v;face_freq_table(v);
00162                         voro_print_vector(v,fp);
00163                 }
00164                 void face_vertices(vector<int> &v);
00165                 /** Outputs the */
00166                 inline void output_face_vertices(FILE *fp=stdout) {
00167                         vector<int> v;face_vertices(v);
00168                         voro_print_face_vertices(v,fp);
00169                 }
00170                 void face_perimeters(vector<double> &v);
00171                 /** Outputs a list of the perimeters of each face.
00172                  * \param[in] fp the file handle to write to. */
00173                 inline void output_face_perimeters(FILE *fp=stdout) {
00174                         vector<double> v;face_perimeters(v);
00175                         voro_print_vector(v,fp);
00176                 }
00177                 void normals(vector<double> &v);
00178                 /** Outputs a list of the perimeters of each face.
00179                  * \param[in] fp the file handle to write to. */
00180                 inline void output_normals(FILE *fp=stdout) {
00181                         vector<double> v;normals(v);
00182                         voro_print_positions(v,fp);
00183                 }
00184                 /** Outputs a custom string of information about the Voronoi
00185                  * cell to a file. It assumes the cell is at (0,0,0) and has a
00186                  * the default_radius associated with it.
00187                  * \param[in] format the custom format string to use.
00188                  * \param[in] fp the file handle to write to. */
00189                 inline void output_custom(const char *format,FILE *fp=stdout) {output_custom(format,0,0,0,0,default_radius,fp);}
00190                 void output_custom(const char *format,int i,double x,double y,double z,double r,FILE *fp=stdout);
00191                 template<class vc_class>
00192                 bool nplane(vc_class &vc,double x,double y,double z,double rs,int p_id);
00193                 bool plane_intersects(double x,double y,double z,double rs);
00194                 bool plane_intersects_guess(double x,double y,double z,double rs);
00195                 void construct_relations();
00196                 void check_relations();
00197                 void check_duplicates();
00198                 void print_edges();
00199                 /** Returns a list of IDs of neighboring particles
00200                  * corresponding to each face.
00201                  * \param[out] v a reference to a vector in which to return the
00202                  *               results. If no neighbor information is
00203                  *               available, a blank vector is returned. */
00204                 virtual void neighbors(vector<int> &v) {v.clear();}
00205                 /** This is a virtual function that is overridden by a routine
00206                  * to print a list of IDs of neighboring particles
00207                  * corresponding to each face. By default, when no neighbor
00208                  * information is available, the routine does nothing.
00209                  * \param[in] fp the file handle to write to. */
00210                 virtual void output_neighbors(FILE *fp=stdout) {}
00211                 /** This a virtual function that is overridden by a routine to
00212                  * print the neighboring particle IDs for a given vertex. By
00213                  * default, when no neighbor information is available, the
00214                  * routine does nothing.
00215                  * \param[in] i the vertex to consider. */
00216                 virtual void print_edges_neighbors(int i) {};
00217                 /** This is a simple inline function for picking out the index
00218                  * of the next edge counterclockwise at the current vertex.
00219                  * \param[in] a the index of an edge of the current vertex.
00220                  * \param[in] p the number of the vertex.
00221                  * \return 0 if a=nu[p]-1, or a+1 otherwise. */
00222                 inline int cycle_up(int a,int p) {return a==nu[p]-1?0:a+1;}
00223                 /** This is a simple inline function for picking out the index
00224                  * of the next edge clockwise from the current vertex.
00225                  * \param[in] a the index of an edge of the current vertex.
00226                  * \param[in] p the number of the vertex.
00227                  * \return nu[p]-1 if a=0, or a-1 otherwise. */
00228                 inline int cycle_down(int a,int p) {return a==0?nu[p]-1:a-1;}
00229         protected:
00230                 /** This a one dimensional array that holds the current sizes
00231                  * of the memory allocations for them mep array.*/
00232                 int *mem;
00233                 /** This is a one dimensional array that holds the current
00234                  * number of vertices of order p that are stored in the mep[p]
00235                  * array. */
00236                 int *mec;
00237                 /** This is a two dimensional array for holding the information
00238                  * about the edges of the Voronoi cell. mep[p] is a
00239                  * one-dimensional array for holding the edge information about
00240                  * all vertices of order p, with each vertex holding 2*p+1
00241                  * integers of information. The total number of vertices held
00242                  * on mep[p] is stored in mem[p]. If the space runs out, the
00243                  * code allocates more using the add_memory() routine. */
00244                 int **mep;
00245                 inline void reset_edges();
00246                 template<class vc_class>
00247                 void check_memory_for_copy(vc_class &vc,voronoicell_base* vb);
00248                 void copy(voronoicell_base* vb);
00249         private:
00250                 /** This is the delete stack, used to store the vertices which
00251                  * are going to be deleted during the plane cutting procedure.
00252                  */
00253                 int *ds,*stacke;
00254                 /** This is the auxiliary delete stack, which has size set by
00255                  * current_delete2_size. */
00256                 int *ds2,*stacke2;
00257                 /** This stores the current memory allocation for the marginal
00258                  * cases. */
00259                 int current_marginal;
00260                 /** This stores the total number of marginal points which are
00261                  * currently in the buffer. */
00262                 int n_marg;
00263                 /** This array contains a list of the marginal points, and also
00264                  * the outcomes of the marginal tests. */
00265                 int *marg;
00266                 /** The x coordinate of the normal vector to the test plane. */
00267                 double px;
00268                 /** The y coordinate of the normal vector to the test plane. */
00269                 double py;
00270                 /** The z coordinate of the normal vector to the test plane. */
00271                 double pz;
00272                 /** The magnitude of the normal vector to the test plane. */
00273                 double prsq;
00274                 template<class vc_class>
00275                 void add_memory(vc_class &vc,int i,int *stackp2);
00276                 template<class vc_class>
00277                 void add_memory_vertices(vc_class &vc);
00278                 template<class vc_class>
00279                 void add_memory_vorder(vc_class &vc);
00280                 void add_memory_ds(int *&stackp);
00281                 void add_memory_ds2(int *&stackp2);
00282                 template<class vc_class>
00283                 inline bool collapse_order1(vc_class &vc);
00284                 template<class vc_class>
00285                 inline bool collapse_order2(vc_class &vc);
00286                 template<class vc_class>
00287                 inline bool delete_connection(vc_class &vc,int j,int k,bool hand);
00288                 template<class vc_class>
00289                 inline bool search_for_outside_edge(vc_class &vc,int &up);
00290                 template<class vc_class>
00291                 inline void add_to_stack(vc_class &vc,int lp,int *&stackp2);
00292                 inline bool plane_intersects_track(double x,double y,double z,double rs,double g);
00293                 inline void normals_search(vector<double> &v,int i,int j,int k);
00294                 inline bool search_edge(int l,int &m,int &k);
00295                 inline int m_test(int n,double &ans);
00296                 int check_marginal(int n,double &ans);
00297                 friend class voronoicell;
00298                 friend class voronoicell_neighbor;
00299 };
00300 
00301 /** \brief Extension of the voronoicell_base class to represent a Voronoi
00302  * cell without neighbor information.
00303  *
00304  * This class is an extension of the voronoicell_base class, in cases when
00305  * is not necessary to track the IDs of neighboring particles associated
00306  * with each face of the Voronoi cell. */
00307 class voronoicell : public voronoicell_base {
00308         public:
00309                 using voronoicell_base::nplane;
00310                 /** Copies the information from another voronoicell class into
00311                  * this class, extending memory allocation if necessary.
00312                  * \param[in] c the class to copy. */
00313                 inline void operator=(voronoicell &c) {
00314                         voronoicell_base* vb((voronoicell_base*) &c);
00315                         check_memory_for_copy(*this,vb);copy(vb);
00316                 }
00317                 /** Cuts a Voronoi cell using by the plane corresponding to the
00318                  * perpendicular bisector of a particle.
00319                  * \param[in] (x,y,z) the position of the particle.
00320                  * \param[in] rsq the modulus squared of the vector.
00321                  * \param[in] p_id the plane ID, ignored for this case where no
00322                  *                 neighbor tracking is enabled.
00323                  * \return False if the plane cut deleted the cell entirely,
00324                  *         true otherwise. */
00325                 inline bool nplane(double x,double y,double z,double rsq,int p_id) {
00326                         return nplane(*this,x,y,z,rsq,0);
00327                 }
00328                 /** Cuts a Voronoi cell using by the plane corresponding to the
00329                  * perpendicular bisector of a particle.
00330                  * \param[in] (x,y,z) the position of the particle.
00331                  * \param[in] p_id the plane ID, ignored for this case where no
00332                  *                 neighbor tracking is enabled.
00333                  * \return False if the plane cut deleted the cell entirely,
00334                  *         true otherwise. */
00335                 inline bool nplane(double x,double y,double z,int p_id) {
00336                         double rsq=x*x+y*y+z*z;
00337                         return nplane(*this,x,y,z,rsq,0);
00338                 }
00339                 /** Cuts a Voronoi cell using by the plane corresponding to the
00340                  * perpendicular bisector of a particle.
00341                  * \param[in] (x,y,z) the position of the particle.
00342                  * \param[in] rsq the modulus squared of the vector.
00343                  * \return False if the plane cut deleted the cell entirely,
00344                  *         true otherwise. */
00345                 inline bool plane(double x,double y,double z,double rsq) {
00346                         return nplane(*this,x,y,z,rsq,0);
00347                 }
00348                 /** Cuts a Voronoi cell using by the plane corresponding to the
00349                  * perpendicular bisector of a particle.
00350                  * \param[in] (x,y,z) the position of the particle.
00351                  * \return False if the plane cut deleted the cell entirely,
00352                  *         true otherwise. */
00353                 inline bool plane(double x,double y,double z) {
00354                         double rsq=x*x+y*y+z*z;
00355                         return nplane(*this,x,y,z,rsq,0);
00356                 }
00357                 /** Initializes the Voronoi cell to be rectangular box with the
00358                  * given dimensions.
00359                  * \param[in] (xmin,xmax) the minimum and maximum x coordinates.
00360                  * \param[in] (ymin,ymax) the minimum and maximum y coordinates.
00361                  * \param[in] (zmin,zmax) the minimum and maximum z coordinates. */
00362                 inline void init(double xmin,double xmax,double ymin,double ymax,double zmin,double zmax) {
00363                         init_base(xmin,xmax,ymin,ymax,zmin,zmax);
00364                 }
00365                 /** Initializes the cell to be an octahedron with vertices at
00366                  * (l,0,0), (-l,0,0), (0,l,0), (0,-l,0), (0,0,l), and (0,0,-l).
00367                  * \param[in] l a parameter setting the size of the octahedron.
00368                  */
00369                 inline void init_octahedron(double l) {
00370                         init_octahedron_base(l);
00371                 }
00372                 /** Initializes the cell to be a tetrahedron.
00373                  * \param[in] (x0,y0,z0) the coordinates of the first vertex.
00374                  * \param[in] (x1,y1,z1) the coordinates of the second vertex.
00375                  * \param[in] (x2,y2,z2) the coordinates of the third vertex.
00376                  * \param[in] (x3,y3,z3) the coordinates of the fourth vertex.
00377                  */
00378                 inline void init_tetrahedron(double x0,double y0,double z0,double x1,double y1,double z1,double x2,double y2,double z2,double x3,double y3,double z3) {
00379                         init_tetrahedron_base(x0,y0,z0,x1,y1,z1,x2,y2,z2,x3,y3,z3);
00380                 }
00381         private:
00382                 inline void n_allocate(int i,int m) {};
00383                 inline void n_add_memory_vertices(int i) {};
00384                 inline void n_add_memory_vorder(int i) {};
00385                 inline void n_set_pointer(int p,int n) {};
00386                 inline void n_copy(int a,int b,int c,int d) {};
00387                 inline void n_set(int a,int b,int c) {};
00388                 inline void n_set_aux1(int k) {};
00389                 inline void n_copy_aux1(int a,int b) {};
00390                 inline void n_copy_aux1_shift(int a,int b) {};
00391                 inline void n_set_aux2_copy(int a,int b) {};
00392                 inline void n_copy_pointer(int a,int b) {};
00393                 inline void n_set_to_aux1(int j) {};
00394                 inline void n_set_to_aux2(int j) {};
00395                 inline void n_allocate_aux1(int i) {};
00396                 inline void n_switch_to_aux1(int i) {};
00397                 inline void n_copy_to_aux1(int i,int m) {};
00398                 inline void n_set_to_aux1_offset(int k,int m) {};
00399                 inline void n_neighbors(vector<int> &v) {v.clear();};
00400                 friend class voronoicell_base;
00401 };
00402 
00403 /** \brief Extension of the voronoicell_base class to represent a Voronoi cell
00404  * with neighbor information.
00405  *
00406  * This class is an extension of the voronoicell_base class, in cases when the
00407  * IDs of neighboring particles associated with each face of the Voronoi cell.
00408  * It contains additional data structures mne and ne for storing this
00409  * information. */
00410 class voronoicell_neighbor : public voronoicell_base {
00411         public:
00412                 using voronoicell_base::nplane;
00413                 /** This two dimensional array holds the neighbor information
00414                  * associated with each vertex. mne[p] is a one dimensional
00415                  * array which holds all of the neighbor information for
00416                  * vertices of order p. */
00417                 int **mne;
00418                 /** This is a two dimensional array that holds the neighbor
00419                  * information associated with each vertex. ne[i] points to a
00420                  * one-dimensional array in mne[nu[i]]. ne[i][j] holds the
00421                  * neighbor information associated with the jth edge of vertex
00422                  * i. It is set to the ID number of the plane that made the
00423                  * face that is clockwise from the jth edge. */
00424                 int **ne;
00425                 voronoicell_neighbor();
00426                 ~voronoicell_neighbor();
00427                 void operator=(voronoicell &c);
00428                 void operator=(voronoicell_neighbor &c);
00429                 /** Cuts the Voronoi cell by a particle whose center is at a
00430                  * separation of (x,y,z) from the cell center. The value of rsq
00431                  * should be initially set to \f$x^2+y^2+z^2\f$.
00432                  * \param[in] (x,y,z) the normal vector to the plane.
00433                  * \param[in] rsq the distance along this vector of the plane.
00434                  * \param[in] p_id the plane ID (for neighbor tracking only).
00435                  * \return False if the plane cut deleted the cell entirely,
00436                  * true otherwise. */
00437                 inline bool nplane(double x,double y,double z,double rsq,int p_id) {
00438                         return nplane(*this,x,y,z,rsq,p_id);
00439                 }
00440                 /** This routine calculates the modulus squared of the vector
00441                  * before passing it to the main nplane() routine with full
00442                  * arguments.
00443                  * \param[in] (x,y,z) the vector to cut the cell by.
00444                  * \param[in] p_id the plane ID (for neighbor tracking only).
00445                  * \return False if the plane cut deleted the cell entirely,
00446                  *         true otherwise. */
00447                 inline bool nplane(double x,double y,double z,int p_id) {
00448                         double rsq=x*x+y*y+z*z;
00449                         return nplane(*this,x,y,z,rsq,p_id);
00450                 }
00451                 /** This version of the plane routine just makes up the plane
00452                  * ID to be zero. It will only be referenced if neighbor
00453                  * tracking is enabled.
00454                  * \param[in] (x,y,z) the vector to cut the cell by.
00455                  * \param[in] rsq the modulus squared of the vector.
00456                  * \return False if the plane cut deleted the cell entirely,
00457                  *         true otherwise. */
00458                 inline bool plane(double x,double y,double z,double rsq) {
00459                         return nplane(*this,x,y,z,rsq,0);
00460                 }
00461                 /** Cuts a Voronoi cell using the influence of a particle at
00462                  * (x,y,z), first calculating the modulus squared of this
00463                  * vector before passing it to the main nplane() routine. Zero
00464                  * is supplied as the plane ID, which will be ignored unless
00465                  * neighbor tracking is enabled.
00466                  * \param[in] (x,y,z) the vector to cut the cell by.
00467                  * \return False if the plane cut deleted the cell entirely,
00468                  *         true otherwise. */
00469                 inline bool plane(double x,double y,double z) {
00470                         double rsq=x*x+y*y+z*z;
00471                         return nplane(*this,x,y,z,rsq,0);
00472                 }
00473                 void init(double xmin,double xmax,double ymin,double ymax,double zmin,double zmax);
00474                 void init_octahedron(double l);
00475                 void init_tetrahedron(double x0,double y0,double z0,double x1,double y1,double z1,double x2,double y2,double z2,double x3,double y3,double z3);
00476                 void check_facets();
00477                 virtual void neighbors(vector<int> &v);
00478                 virtual void print_edges_neighbors(int i);
00479                 virtual void output_neighbors(FILE *fp=stdout) {
00480                         vector<int> v;neighbors(v);
00481                         voro_print_vector(v,fp);
00482                 }
00483         private:
00484                 int *paux1;
00485                 int *paux2;
00486                 inline void n_allocate(int i,int m) {mne[i]=new int[m*i];}
00487                 inline void n_add_memory_vertices(int i) {
00488                         int **pp(new int*[i]);
00489                         for(int j=0;j<current_vertices;j++) pp[j]=ne[j];
00490                         delete [] ne;ne=pp;
00491                 }
00492                 inline void n_add_memory_vorder(int i) {
00493                         int **p2(new int*[i]);
00494                         for(int j=0;j<current_vertex_order;j++) p2[j]=mne[j];
00495                         delete [] mne;mne=p2;
00496                 }
00497                 inline void n_set_pointer(int p,int n) {
00498                         ne[p]=mne[n]+n*mec[n];
00499                 }
00500                 inline void n_copy(int a,int b,int c,int d) {ne[a][b]=ne[c][d];}
00501                 inline void n_set(int a,int b,int c) {ne[a][b]=c;}
00502                 inline void n_set_aux1(int k) {paux1=mne[k]+k*mec[k];}
00503                 inline void n_copy_aux1(int a,int b) {paux1[b]=ne[a][b];}
00504                 inline void n_copy_aux1_shift(int a,int b) {paux1[b]=ne[a][b+1];}
00505                 inline void n_set_aux2_copy(int a,int b) {
00506                         paux2=mne[b]+b*mec[b];
00507                         for(int i=0;i<b;i++) ne[a][i]=paux2[i];
00508                 }
00509                 inline void n_copy_pointer(int a,int b) {ne[a]=ne[b];}
00510                 inline void n_set_to_aux1(int j) {ne[j]=paux1;}
00511                 inline void n_set_to_aux2(int j) {ne[j]=paux2;}
00512                 inline void n_allocate_aux1(int i) {paux1=new int[i*mem[i]];}
00513                 inline void n_switch_to_aux1(int i) {delete [] mne[i];mne[i]=paux1;}
00514                 inline void n_copy_to_aux1(int i,int m) {paux1[m]=mne[i][m];}
00515                 inline void n_set_to_aux1_offset(int k,int m) {ne[k]=paux1+m;}
00516                 friend class voronoicell_base;
00517 };
00518 
00519 }
00520 
00521 #endif