Voro++
|
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=safe_fopen(filename,"w"); 00104 draw_pov(x,y,z,fp); 00105 fclose(fp); 00106 }; 00107 void draw_pov_mesh(double x,double y,double z,FILE *fp=stdout); 00108 /** Outputs the cell in POV-Ray format as a mesh2 object to a 00109 * given file. 00110 * \param[in] (x,y,z) a displacement to add to the cell's 00111 * position. 00112 * \param[in] filename the name of the file to write to. */ 00113 inline void draw_pov_mesh(double x,double y,double z,const char *filename) { 00114 FILE *fp=safe_fopen(filename,"w"); 00115 draw_pov_mesh(x,y,z,fp); 00116 fclose(fp); 00117 } 00118 void draw_gnuplot(double x,double y,double z,FILE *fp=stdout); 00119 /** Outputs the cell in Gnuplot format a given file. 00120 * \param[in] (x,y,z) a displacement to add to the cell's 00121 * position. 00122 * \param[in] filename the name of the file to write to. */ 00123 inline void draw_gnuplot(double x,double y,double z,const char *filename) { 00124 FILE *fp=safe_fopen(filename,"w"); 00125 draw_gnuplot(x,y,z,fp); 00126 fclose(fp); 00127 } 00128 double volume(); 00129 double max_radius_squared(); 00130 double total_edge_distance(); 00131 double surface_area(); 00132 void centroid(double &cx,double &cy,double &cz); 00133 int number_of_faces(); 00134 int number_of_edges(); 00135 void vertex_orders(vector<int> &v); 00136 void output_vertex_orders(FILE *fp=stdout); 00137 void vertices(vector<double> &v); 00138 void output_vertices(FILE *fp=stdout); 00139 void vertices(double x,double y,double z,vector<double> &v); 00140 void output_vertices(double x,double y,double z,FILE *fp=stdout); 00141 void face_areas(vector<double> &v); 00142 /** Outputs the areas of the faces. 00143 * \param[in] fp the file handle to write to. */ 00144 inline void output_face_areas(FILE *fp=stdout) { 00145 vector<double> v;face_areas(v); 00146 voro_print_vector(v,fp); 00147 } 00148 void face_orders(vector<int> &v); 00149 /** Outputs a list of the number of sides of each face. 00150 * \param[in] fp the file handle to write to. */ 00151 inline void output_face_orders(FILE *fp=stdout) { 00152 vector<int> v;face_orders(v); 00153 voro_print_vector(v,fp); 00154 } 00155 void face_freq_table(vector<int> &v); 00156 /** Outputs a */ 00157 inline void output_face_freq_table(FILE *fp=stdout) { 00158 vector<int> v;face_freq_table(v); 00159 voro_print_vector(v,fp); 00160 } 00161 void face_vertices(vector<int> &v); 00162 /** Outputs the */ 00163 inline void output_face_vertices(FILE *fp=stdout) { 00164 vector<int> v;face_vertices(v); 00165 voro_print_face_vertices(v,fp); 00166 } 00167 void face_perimeters(vector<double> &v); 00168 /** Outputs a list of the perimeters of each face. 00169 * \param[in] fp the file handle to write to. */ 00170 inline void output_face_perimeters(FILE *fp=stdout) { 00171 vector<double> v;face_perimeters(v); 00172 voro_print_vector(v,fp); 00173 } 00174 void normals(vector<double> &v); 00175 /** Outputs a list of the perimeters of each face. 00176 * \param[in] fp the file handle to write to. */ 00177 inline void output_normals(FILE *fp=stdout) { 00178 vector<double> v;normals(v); 00179 voro_print_positions(v,fp); 00180 } 00181 /** Outputs a custom string of information about the Voronoi 00182 * cell to a file. It assumes the cell is at (0,0,0) and has a 00183 * the default_radius associated with it. 00184 * \param[in] format the custom format string to use. 00185 * \param[in] fp the file handle to write to. */ 00186 inline void output_custom(const char *format,FILE *fp=stdout) {output_custom(format,0,0,0,0,default_radius,fp);} 00187 void output_custom(const char *format,int i,double x,double y,double z,double r,FILE *fp=stdout); 00188 template<class vc_class> 00189 bool nplane(vc_class &vc,double x,double y,double z,double rs,int p_id); 00190 bool plane_intersects(double x,double y,double z,double rs); 00191 bool plane_intersects_guess(double x,double y,double z,double rs); 00192 void construct_relations(); 00193 void check_relations(); 00194 void check_duplicates(); 00195 void print_edges(); 00196 /** Returns a list of IDs of neighboring particles 00197 * corresponding to each face. 00198 * \param[out] v a reference to a vector in which to return the 00199 * results. If no neighbor information is 00200 * available, a blank vector is returned. */ 00201 virtual void neighbors(vector<int> &v) {v.clear();} 00202 /** This is a virtual function that is overridden by a routine 00203 * to print a list of IDs of neighboring particles 00204 * corresponding to each face. By default, when no neighbor 00205 * information is available, the routine does nothing. 00206 * \param[in] fp the file handle to write to. */ 00207 virtual void output_neighbors(FILE *fp=stdout) {} 00208 /** This a virtual function that is overridden by a routine to 00209 * print the neighboring particle IDs for a given vertex. By 00210 * default, when no neighbor information is available, the 00211 * routine does nothing. 00212 * \param[in] i the vertex to consider. */ 00213 virtual void print_edges_neighbors(int i) {}; 00214 /** This is a simple inline function for picking out the index 00215 * of the next edge counterclockwise at the current vertex. 00216 * \param[in] a the index of an edge of the current vertex. 00217 * \param[in] p the number of the vertex. 00218 * \return 0 if a=nu[p]-1, or a+1 otherwise. */ 00219 inline int cycle_up(int a,int p) {return a==nu[p]-1?0:a+1;} 00220 /** This is a simple inline function for picking out the index 00221 * of the next edge clockwise from the current vertex. 00222 * \param[in] a the index of an edge of the current vertex. 00223 * \param[in] p the number of the vertex. 00224 * \return nu[p]-1 if a=0, or a-1 otherwise. */ 00225 inline int cycle_down(int a,int p) {return a==0?nu[p]-1:a-1;} 00226 protected: 00227 /** This a one dimensional array that holds the current sizes 00228 * of the memory allocations for them mep array.*/ 00229 int *mem; 00230 /** This is a one dimensional array that holds the current 00231 * number of vertices of order p that are stored in the mep[p] 00232 * array. */ 00233 int *mec; 00234 /** This is a two dimensional array for holding the information 00235 * about the edges of the Voronoi cell. mep[p] is a 00236 * one-dimensional array for holding the edge information about 00237 * all vertices of order p, with each vertex holding 2*p+1 00238 * integers of information. The total number of vertices held 00239 * on mep[p] is stored in mem[p]. If the space runs out, the 00240 * code allocates more using the add_memory() routine. */ 00241 int **mep; 00242 inline void reset_edges(); 00243 template<class vc_class> 00244 void check_memory_for_copy(vc_class &vc,voronoicell_base* vb); 00245 void copy(voronoicell_base* vb); 00246 private: 00247 /** This is the delete stack, used to store the vertices which 00248 * are going to be deleted during the plane cutting procedure. 00249 */ 00250 int *ds,*stacke; 00251 /** This is the auxiliary delete stack, which has size set by 00252 * current_delete2_size. */ 00253 int *ds2,*stacke2; 00254 /** This stores the current memory allocation for the marginal 00255 * cases. */ 00256 int current_marginal; 00257 /** This stores the total number of marginal points which are 00258 * currently in the buffer. */ 00259 int n_marg; 00260 /** This array contains a list of the marginal points, and also 00261 * the outcomes of the marginal tests. */ 00262 int *marg; 00263 /** The x coordinate of the normal vector to the test plane. */ 00264 double px; 00265 /** The y coordinate of the normal vector to the test plane. */ 00266 double py; 00267 /** The z coordinate of the normal vector to the test plane. */ 00268 double pz; 00269 /** The magnitude of the normal vector to the test plane. */ 00270 double prsq; 00271 template<class vc_class> 00272 void add_memory(vc_class &vc,int i,int *stackp2); 00273 template<class vc_class> 00274 void add_memory_vertices(vc_class &vc); 00275 template<class vc_class> 00276 void add_memory_vorder(vc_class &vc); 00277 void add_memory_ds(int *&stackp); 00278 void add_memory_ds2(int *&stackp2); 00279 template<class vc_class> 00280 inline bool collapse_order1(vc_class &vc); 00281 template<class vc_class> 00282 inline bool collapse_order2(vc_class &vc); 00283 template<class vc_class> 00284 inline bool delete_connection(vc_class &vc,int j,int k,bool hand); 00285 template<class vc_class> 00286 inline bool search_for_outside_edge(vc_class &vc,int &up); 00287 template<class vc_class> 00288 inline void add_to_stack(vc_class &vc,int lp,int *&stackp2); 00289 inline bool plane_intersects_track(double x,double y,double z,double rs,double g); 00290 inline void normals_search(vector<double> &v,int i,int j,int k); 00291 inline bool search_edge(int l,int &m,int &k); 00292 inline int m_test(int n,double &ans); 00293 int check_marginal(int n,double &ans); 00294 friend class voronoicell; 00295 friend class voronoicell_neighbor; 00296 }; 00297 00298 /** \brief Extension of the voronoicell_base class to represent a Voronoi 00299 * cell without neighbor information. 00300 * 00301 * This class is an extension of the voronoicell_base class, in cases when 00302 * is not necessary to track the IDs of neighboring particles associated 00303 * with each face of the Voronoi cell. */ 00304 class voronoicell : public voronoicell_base { 00305 public: 00306 using voronoicell_base::nplane; 00307 /** Copies the information from another voronoicell class into 00308 * this class, extending memory allocation if necessary. 00309 * \param[in] c the class to copy. */ 00310 inline void operator=(voronoicell &c) { 00311 voronoicell_base* vb((voronoicell_base*) &c); 00312 check_memory_for_copy(*this,vb);copy(vb); 00313 } 00314 /** Cuts a Voronoi cell using by the plane corresponding to the 00315 * perpendicular bisector of a particle. 00316 * \param[in] (x,y,z) the position of the particle. 00317 * \param[in] rsq the modulus squared of the vector. 00318 * \param[in] p_id the plane ID, ignored for this case where no 00319 * neighbor tracking is enabled. 00320 * \return False if the plane cut deleted the cell entirely, 00321 * true otherwise. */ 00322 inline bool nplane(double x,double y,double z,double rsq,int p_id) { 00323 return nplane(*this,x,y,z,rsq,0); 00324 } 00325 /** Cuts a Voronoi cell using by the plane corresponding to the 00326 * perpendicular bisector of a particle. 00327 * \param[in] (x,y,z) the position of the particle. 00328 * \param[in] p_id the plane ID, ignored for this case where no 00329 * neighbor tracking is enabled. 00330 * \return False if the plane cut deleted the cell entirely, 00331 * true otherwise. */ 00332 inline bool nplane(double x,double y,double z,int p_id) { 00333 double rsq=x*x+y*y+z*z; 00334 return nplane(*this,x,y,z,rsq,0); 00335 } 00336 /** Cuts a Voronoi cell using by the plane corresponding to the 00337 * perpendicular bisector of a particle. 00338 * \param[in] (x,y,z) the position of the particle. 00339 * \param[in] rsq the modulus squared of the vector. 00340 * \return False if the plane cut deleted the cell entirely, 00341 * true otherwise. */ 00342 inline bool plane(double x,double y,double z,double rsq) { 00343 return nplane(*this,x,y,z,rsq,0); 00344 } 00345 /** Cuts a Voronoi cell using by the plane corresponding to the 00346 * perpendicular bisector of a particle. 00347 * \param[in] (x,y,z) the position of the particle. 00348 * \return False if the plane cut deleted the cell entirely, 00349 * true otherwise. */ 00350 inline bool plane(double x,double y,double z) { 00351 double rsq=x*x+y*y+z*z; 00352 return nplane(*this,x,y,z,rsq,0); 00353 } 00354 /** Initializes the Voronoi cell to be rectangular box with the 00355 * given dimensions. 00356 * \param[in] (xmin,xmax) the minimum and maximum x coordinates. 00357 * \param[in] (ymin,ymax) the minimum and maximum y coordinates. 00358 * \param[in] (zmin,zmax) the minimum and maximum z coordinates. */ 00359 inline void init(double xmin,double xmax,double ymin,double ymax,double zmin,double zmax) { 00360 init_base(xmin,xmax,ymin,ymax,zmin,zmax); 00361 } 00362 /** Initializes the cell to be an octahedron with vertices at 00363 * (l,0,0), (-l,0,0), (0,l,0), (0,-l,0), (0,0,l), and (0,0,-l). 00364 * \param[in] l a parameter setting the size of the octahedron. 00365 */ 00366 inline void init_octahedron(double l) { 00367 init_octahedron_base(l); 00368 } 00369 /** Initializes the cell to be a tetrahedron. 00370 * \param[in] (x0,y0,z0) the coordinates of the first vertex. 00371 * \param[in] (x1,y1,z1) the coordinates of the second vertex. 00372 * \param[in] (x2,y2,z2) the coordinates of the third vertex. 00373 * \param[in] (x3,y3,z3) the coordinates of the fourth vertex. 00374 */ 00375 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) { 00376 init_tetrahedron_base(x0,y0,z0,x1,y1,z1,x2,y2,z2,x3,y3,z3); 00377 } 00378 private: 00379 inline void n_allocate(int i,int m) {}; 00380 inline void n_add_memory_vertices(int i) {}; 00381 inline void n_add_memory_vorder(int i) {}; 00382 inline void n_set_pointer(int p,int n) {}; 00383 inline void n_copy(int a,int b,int c,int d) {}; 00384 inline void n_set(int a,int b,int c) {}; 00385 inline void n_set_aux1(int k) {}; 00386 inline void n_copy_aux1(int a,int b) {}; 00387 inline void n_copy_aux1_shift(int a,int b) {}; 00388 inline void n_set_aux2_copy(int a,int b) {}; 00389 inline void n_copy_pointer(int a,int b) {}; 00390 inline void n_set_to_aux1(int j) {}; 00391 inline void n_set_to_aux2(int j) {}; 00392 inline void n_allocate_aux1(int i) {}; 00393 inline void n_switch_to_aux1(int i) {}; 00394 inline void n_copy_to_aux1(int i,int m) {}; 00395 inline void n_set_to_aux1_offset(int k,int m) {}; 00396 inline void n_neighbors(vector<int> &v) {v.clear();}; 00397 friend class voronoicell_base; 00398 }; 00399 00400 /** \brief Extension of the voronoicell_base class to represent a Voronoi cell 00401 * with neighbor information. 00402 * 00403 * This class is an extension of the voronoicell_base class, in cases when the 00404 * IDs of neighboring particles associated with each face of the Voronoi cell. 00405 * It contains additional data structures mne and ne for storing this 00406 * information. */ 00407 class voronoicell_neighbor : public voronoicell_base { 00408 public: 00409 using voronoicell_base::nplane; 00410 /** This two dimensional array holds the neighbor information 00411 * associated with each vertex. mne[p] is a one dimensional 00412 * array which holds all of the neighbor information for 00413 * vertices of order p. */ 00414 int **mne; 00415 /** This is a two dimensional array that holds the neighbor 00416 * information associated with each vertex. ne[i] points to a 00417 * one-dimensional array in mne[nu[i]]. ne[i][j] holds the 00418 * neighbor information associated with the jth edge of vertex 00419 * i. It is set to the ID number of the plane that made the 00420 * face that is clockwise from the jth edge. */ 00421 int **ne; 00422 voronoicell_neighbor(); 00423 ~voronoicell_neighbor(); 00424 void operator=(voronoicell &c); 00425 void operator=(voronoicell_neighbor &c); 00426 /** Cuts the Voronoi cell by a particle whose center is at a 00427 * separation of (x,y,z) from the cell center. The value of rsq 00428 * should be initially set to \f$x^2+y^2+z^2\f$. 00429 * \param[in] (x,y,z) the normal vector to the plane. 00430 * \param[in] rsq the distance along this vector of the plane. 00431 * \param[in] p_id the plane ID (for neighbor tracking only). 00432 * \return False if the plane cut deleted the cell entirely, 00433 * true otherwise. */ 00434 inline bool nplane(double x,double y,double z,double rsq,int p_id) { 00435 return nplane(*this,x,y,z,rsq,p_id); 00436 } 00437 /** This routine calculates the modulus squared of the vector 00438 * before passing it to the main nplane() routine with full 00439 * arguments. 00440 * \param[in] (x,y,z) the vector to cut the cell by. 00441 * \param[in] p_id the plane ID (for neighbor tracking only). 00442 * \return False if the plane cut deleted the cell entirely, 00443 * true otherwise. */ 00444 inline bool nplane(double x,double y,double z,int p_id) { 00445 double rsq=x*x+y*y+z*z; 00446 return nplane(*this,x,y,z,rsq,p_id); 00447 } 00448 /** This version of the plane routine just makes up the plane 00449 * ID to be zero. It will only be referenced if neighbor 00450 * tracking is enabled. 00451 * \param[in] (x,y,z) the vector to cut the cell by. 00452 * \param[in] rsq the modulus squared of the vector. 00453 * \return False if the plane cut deleted the cell entirely, 00454 * true otherwise. */ 00455 inline bool plane(double x,double y,double z,double rsq) { 00456 return nplane(*this,x,y,z,rsq,0); 00457 } 00458 /** Cuts a Voronoi cell using the influence of a particle at 00459 * (x,y,z), first calculating the modulus squared of this 00460 * vector before passing it to the main nplane() routine. Zero 00461 * is supplied as the plane ID, which will be ignored unless 00462 * neighbor tracking is enabled. 00463 * \param[in] (x,y,z) the vector to cut the cell by. 00464 * \return False if the plane cut deleted the cell entirely, 00465 * true otherwise. */ 00466 inline bool plane(double x,double y,double z) { 00467 double rsq=x*x+y*y+z*z; 00468 return nplane(*this,x,y,z,rsq,0); 00469 } 00470 void init(double xmin,double xmax,double ymin,double ymax,double zmin,double zmax); 00471 void init_octahedron(double l); 00472 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); 00473 void check_facets(); 00474 virtual void neighbors(vector<int> &v); 00475 virtual void print_edges_neighbors(int i); 00476 virtual void output_neighbors(FILE *fp=stdout) { 00477 vector<int> v;neighbors(v); 00478 voro_print_vector(v,fp); 00479 } 00480 private: 00481 int *paux1; 00482 int *paux2; 00483 inline void n_allocate(int i,int m) {mne[i]=new int[m*i];} 00484 inline void n_add_memory_vertices(int i) { 00485 int **pp=new int*[i]; 00486 for(int j=0;j<current_vertices;j++) pp[j]=ne[j]; 00487 delete [] ne;ne=pp; 00488 } 00489 inline void n_add_memory_vorder(int i) { 00490 int **p2=new int*[i]; 00491 for(int j=0;j<current_vertex_order;j++) p2[j]=mne[j]; 00492 delete [] mne;mne=p2; 00493 } 00494 inline void n_set_pointer(int p,int n) { 00495 ne[p]=mne[n]+n*mec[n]; 00496 } 00497 inline void n_copy(int a,int b,int c,int d) {ne[a][b]=ne[c][d];} 00498 inline void n_set(int a,int b,int c) {ne[a][b]=c;} 00499 inline void n_set_aux1(int k) {paux1=mne[k]+k*mec[k];} 00500 inline void n_copy_aux1(int a,int b) {paux1[b]=ne[a][b];} 00501 inline void n_copy_aux1_shift(int a,int b) {paux1[b]=ne[a][b+1];} 00502 inline void n_set_aux2_copy(int a,int b) { 00503 paux2=mne[b]+b*mec[b]; 00504 for(int i=0;i<b;i++) ne[a][i]=paux2[i]; 00505 } 00506 inline void n_copy_pointer(int a,int b) {ne[a]=ne[b];} 00507 inline void n_set_to_aux1(int j) {ne[j]=paux1;} 00508 inline void n_set_to_aux2(int j) {ne[j]=paux2;} 00509 inline void n_allocate_aux1(int i) {paux1=new int[i*mem[i]];} 00510 inline void n_switch_to_aux1(int i) {delete [] mne[i];mne[i]=paux1;} 00511 inline void n_copy_to_aux1(int i,int m) {paux1[m]=mne[i][m];} 00512 inline void n_set_to_aux1_offset(int k,int m) {ne[k]=paux1+m;} 00513 friend class voronoicell_base; 00514 }; 00515 00516 } 00517 00518 #endif