frustum.cc – Using the cone wall object to create a frustum

Voronoi tessellation in a frustum

This code makes use of the conical wall object to carry out a Voronoi tessellation for particles in a circular frustum. A container class is created using non-periodic boundary conditions, where the z coordinate ranges from zero to one. On line 23, a wall_cone class is constructed, with the first three arguments specifying the apex position at (0,0,2), and the next three coordinates specifying the axis direction along (0,0,-1). The final argument sets the angle of the cone to be atan(0.5). The intersection of this cone with the container is a frustum, where the base at z=0 has radius 1, and the top at z=1 has radius 0.5.

In lines 28 to 32, a rectangular grid of particles is added to container, only allowing those points that return true to the point_inside() routine. The remainder of the code outputs the particle positions and Voronoi cells in Gnuplot and POV-Ray formats.

The image above was rendered in POV-Ray using the header file “frustum.pov”. The cells mesh together to form the frustum, although some discrepancies are visible between them. This is due to the fact that the wall_cone class approximates the conical wall surface with a single plane cut, in the same way as for the cylinder example. In lines 44 to 48 the accuracy of this approximation is tested, by comparing the sum of Voronoi volumes, to the exact frustum volume:

Exact frustum volume : 1.8326
Voronoi cell volume  : 1.85132
Difference           : 0.0187281

This is an error of 1.02%. Increasing the number of inserted particles or improving the accuracy of the wall_cone class would reduce this.

Movies

Code listing

 1: // Frustum example code
 2: //
 3: // Author   : Chris H. Rycroft (LBL / UC Berkeley)
 4: // Email    : chr@alum.mit.edu
 5: // Date     : August 30th 2011
 6: 
 7: #include "voro++.hh"
 8: using namespace voro;
 9: 
10: const double pi=3.1415926535897932384626433832795;
11: 
12: int main() {
13:         int i=0;
14:         double x,y,z,evol,vvol;
15: 
16:         // Create a container with the geometry given above, and make it
17:         // non-periodic in each of the three coordinates. Allocate space for
18:         // eight particles within each computational block.
19:         container con(-1.2,1.2,-1.2,1.2,0,1,14,14,7,
20:                         false,false,false,8);
21: 
22:         // Add a cylindrical wall to the container
23:         wall_cone cone(0,0,2,0,0,-1,atan(0.5));
24:         con.add_wall(cone);
25: 
26:         // Place particles in a regular grid within the frustum, for points
27:         // which are within the wall boundaries
28:         for(z=0.1;z<1;z+=0.2) for(y=-0.85;y<1;y+=0.2) for(x=-0.95;x<1;x+=0.2) {
29:                 if (con.point_inside(x,y,z)) {
30:                         con.put(i,x,y,z);i++;
31:                 }
32:         }
33: 
34:         // Output the particle positions and Voronoi cells in Gnuplot format
35:         con.draw_particles("frustum_p.gnu");
36:         con.draw_cells_gnuplot("frustum_v.gnu");
37: 
38:         // Output the particle positions and Voronoi cells in POV-Ray format
39:         con.draw_particles_pov("frustum_p.pov");
40:         con.draw_cells_pov("frustum_v.pov");
41: 
42:         // Compute the volume of the Voronoi cells and compare it to the
43:         // exact frustum volume
44:         evol=pi*1*(0.5*0.5+0.5*1+1*1)/3;
45:         vvol=con.sum_cell_volumes();
46:         printf("Exact frustum volume : %g\n"
47:                "Voronoi cell volume  : %g\n"
48:                "Difference           : %g\n",evol,vvol,vvol-evol);
49: }