Bugfix in search_for_outside_edge routine.
[voro++.git] / branches / exact / examples / interface / polygons.cc
blob205a6e4e1219a37031cf7feadffd4a88a292506c
1 // Direct C++ interface example code
2 //
3 // Author : Chris H. Rycroft (LBL / UC Berkeley)
4 // Email : chr@alum.mit.edu
5 // Date : August 30th 2011
7 #include "voro++.hh"
8 using namespace voro;
10 void draw_polygon(FILE *fp,vector<int> &f_vert,vector<double> &v,int j);
12 int main() {
13 unsigned int i,j;
14 int id,nx,ny,nz;
15 double x,y,z;
16 voronoicell_neighbor c;
17 vector<int> neigh,f_vert;
18 vector<double> v;
20 // Create a pre-container class to import the input file and guess the
21 // best computational grid size to use.
22 pre_container pcon(-3,3,-3,3,0,6,false,false,false);
23 pcon.import("pack_six_cube");
24 pcon.guess_optimal(nx,ny,nz);
26 // Set up the container class and import the particles from the
27 // pre-container
28 container con(-3,3,-3,3,0,6,nx,ny,nz,false,false,false,8);
29 pcon.setup(con);
31 // Open the output files
32 FILE *fp4=safe_fopen("polygons4_v.pov","w"),
33 *fp5=safe_fopen("polygons5_v.pov","w"),
34 *fp6=safe_fopen("polygons6_v.pov","w");
36 // Loop over all particles in the container and compute each Voronoi
37 // cell
38 c_loop_all cl(con);
39 if(cl.start()) do if(con.compute_cell(c,cl)) {
40 cl.pos(x,y,z);id=cl.pid();
42 // Gather information about the computed Voronoi cell
43 c.neighbors(neigh);
44 c.face_vertices(f_vert);
45 c.vertices(x,y,z,v);
47 // Loop over all faces of the Voronoi cell
48 for(i=0,j=0;i<neigh.size();i++) {
50 // Draw all quadrilaterals, pentagons, and hexagons.
51 // Skip if the neighbor information is smaller than
52 // this particle's ID, to avoid double counting. This
53 // also removes faces that touch the walls, since the
54 // neighbor information is set to negative numbers for
55 // these cases.
56 if(neigh[i]>id) {
57 switch(f_vert[j]) {
58 case 4: draw_polygon(fp4,f_vert,v,j);
59 break;
60 case 5: draw_polygon(fp5,f_vert,v,j);
61 break;
62 case 6: draw_polygon(fp6,f_vert,v,j);
66 // Skip to the next entry in the face vertex list
67 j+=f_vert[j]+1;
69 } while (cl.inc());
71 // Close the output files
72 fclose(fp4);
73 fclose(fp5);
74 fclose(fp6);
76 // Draw the outline of the domain
77 con.draw_domain_pov("polygons_d.pov");
80 void draw_polygon(FILE *fp,vector<int> &f_vert,vector<double> &v,int j) {
81 static char s[6][128];
82 int k,l,n=f_vert[j];
84 // Create POV-Ray vector strings for each of the vertices
85 for(k=0;k<n;k++) {
86 l=3*f_vert[j+k+1];
87 sprintf(s[k],"<%g,%g,%g>",v[l],v[l+1],v[l+2]);
90 // Draw the interior of the polygon
91 fputs("union{\n",fp);
92 for(k=2;k<n;k++) fprintf(fp,"\ttriangle{%s,%s,%s}\n",s[0],s[k-1],s[k]);
93 fputs("\ttexture{t1}\n}\n",fp);
95 // Draw the outline of the polygon
96 fputs("union{\n",fp);
97 for(k=0;k<n;k++) {
98 l=(k+1)%n;
99 fprintf(fp,"\tcylinder{%s,%s,r}\n\tsphere{%s,r}\n",
100 s[k],s[l],s[l]);
102 fputs("\ttexture{t2}\n}\n",fp);