#include "voro++.hh"
using namespace voro;
#include <vector>
using namespace std;
const double x_min=-15,x_max=15;
const double y_min=-7,y_max=7;
const double z_min=-15,z_max=15;
const double Phi=0.5*(1+sqrt(5.0));
const double phi=0.5*(1-sqrt(5.0));
const int n_x=8,n_y=8,n_z=8;
const int wid=-10;
class wall_initial_shape : public wall {
public:
wall_initial_shape() {
v.init(-2,2,-2,2,-2,2);
v.nplane(0,Phi,1,wid);v.nplane(0,-Phi,1,wid);v.nplane(0,Phi,-1,wid);
v.nplane(0,-Phi,-1,wid);v.nplane(1,0,Phi,wid);v.nplane(-1,0,Phi,wid);
v.nplane(1,0,-Phi,wid);v.nplane(-1,0,-Phi,wid);v.nplane(Phi,1,0,wid);
v.nplane(-Phi,1,0,wid);v.nplane(Phi,-1,0,wid);v.nplane(-Phi,-1,0,wid);
};
bool point_inside(double x,double y,double z) {return true;}
bool cut_cell(voronoicell &c,double x,double y,double z) {
return true;
}
bool cut_cell(voronoicell_neighbor &c,double x,double y,double z) {
c=v;
return true;
}
private:
voronoicell_neighbor v;
};
bool has_dodec_sides(vector<int> &vi) {
for(unsigned int i=0;i<vi.size();i++) if(vi[i]==wid) return true;
return false;
}
int main() {
container con(x_min,x_max,y_min,y_max,z_min,z_max,n_x,n_y,n_z,
false,false,false,8);
wall_initial_shape(wis);
con.add_wall(wis);
con.import("pack_semicircle");
FILE *finside=safe_fopen("finite_sys_in.pov","w"),
*foutside=safe_fopen("finite_sys_out.pov","w");
double x,y,z;
vector<int> vi;
voronoicell_neighbor c;
c_loop_all cl(con);
if(cl.start()) do {
cl.pos(x,y,z);
if(y<0) continue;
if(con.compute_cell(c,cl)) {
c.neighbors(vi);
fprintf(has_dodec_sides(vi)?foutside:finside,
"sphere{<%g,%g,%g>,s}\n",x,y,z);
}
} while(cl.inc());
fclose(finside);
fclose(foutside);
}