poly_occ/src/poly_voro.cpp

#include "poly_voro.h"
#include "poly_point.h"
#include "poly_noeud.h"
#include "poly_face.h"
#include "poly_cellule.h"


#include "tpl_grille.h"


#include <algorithm>

#include "../voro++-0.4.6/src/voro++.hh"

Poly_Voro::Poly_Voro(std::vector<Poly_Point*> list_pnts,double dg): POLY_AFFICHE(),list_points(list_pnts),epsfusion(dg)
{
}

void Poly_Voro::construit(void)
{
  // Géométrie du conteneur
  double x_min= list_points[0]->get_x(), x_max= list_points[6]->get_x();
  double y_min= list_points[0]->get_y(), y_max= list_points[6]->get_y();
  double z_min= list_points[0]->get_z(), z_max= list_points[6]->get_z();
  double cvol = (x_max - x_min) * (y_max - y_min) * (x_max - x_min);

  // Nombre de blocs qui divise le container (blocs calculatoires)
  int n_x= 6, n_y= 6, n_z= 6;

  // Cré un container avec la géométrie données
  voro::container con(x_min, x_max, y_min, y_max, z_min, z_max, n_x, n_y, n_z, false, false, false, 8);

  // Ajoute les particules au container, on n'insère pas les 8 coins du cube
  double x,y,z;
  for(int i=8; i < list_points.size(); i++) {
      x = list_points[i]->get_x();
      y = list_points[i]->get_y();
      z = list_points[i]->get_z();
      con.put(i, x, y, z);
  }
  
  // Additionne les volumes et vérifi que ça concorde avec le volume du container
  double vvol = con.sum_cell_volumes();
  char message[2000];
  sprintf(message,"   Volume de l'echantillon = %lf",cvol);
  affiche(message);
  sprintf(message,"   Volume des cellules de Voronoi = %lf",vvol);
  affiche(message);
  
  
  voro::voronoicell voro_cell;
  voro::c_loop_all loop(con);
  
  // loop a travers les cellules
  if(loop.start()) do if(con.compute_cell(voro_cell,loop)) 
  {
    Poly_Cellule* cell = new Poly_Cellule;
    add_cell(cell);
    loop.pos(x,y,z); // position globale de la particule en cours
    
    // récupération des noeuds
    std::vector< double > vert_xyz;
    voro_cell.vertices(x, y, z, vert_xyz); // récupère tous les vertex de la cellule dans le système global
    
    for(int i=0; i < vert_xyz.size(); i+=3)
    {
      static int numeronoeud=0;
      numeronoeud++;
      Poly_Noeud* nd = new Poly_Noeud(numeronoeud, vert_xyz[i], vert_xyz[i+1], vert_xyz[i+2]);
      cell->add_noeud(nd);
    }
    
    // récupération des faces
    std::vector< int > face_vert;
    voro_cell.face_vertices(face_vert); // liste de vertex pour chaque face
    
    for(int i=0; i < face_vert.size(); i++)
    {
      Poly_Face* face = new Poly_Face;
      int nbVert = face_vert[i]; // première valeur = le nombre de vertex de la face
      for(int j=i+1; j <  i+nbVert+1; j++)
      {
        face->add_noeud(face_vert[j]);
      }
      cell->add_face(face);
      i += nbVert;
    }
  } while (loop.inc());
}

// destructor
Poly_Voro::~Poly_Voro()
{
    for (int i=0; i<get_nb_cell(); i++)
    {
        Poly_Cellule* c = list_cellules[i];
        delete c;
    }      
}

void Poly_Voro::fusion_noeuds(void)
{
  TPL_MAP_ENTITE<Poly_Noeud*> lstnoeud;
  TPL_GRILLE<Poly_Noeud*>  grille;
  double x_min= list_points[0]->get_x(), x_max= list_points[6]->get_x();
  double y_min= list_points[0]->get_y(), y_max= list_points[6]->get_y();
  double z_min= list_points[0]->get_z(), z_max= list_points[6]->get_z();
  BOITE_3D boiteenglobante(x_min,y_min,z_min,x_max,y_max,z_max);
  boiteenglobante.change_grosseur(1.1);
  int nb=pow(list_points.size(),0.333333333333333333)+1;
  grille.initialiser(boiteenglobante.get_xmin(),boiteenglobante.get_ymin(),boiteenglobante.get_zmin(),boiteenglobante.get_xmax(),boiteenglobante.get_ymax(),boiteenglobante.get_zmax(),nb,nb,nb);
  
  double eps = epsfusion;
  
  
  for(int i=0; i < get_nb_cell(); i++)
  {
    Poly_Cellule* c = get_cell(i);
    for(int j=0; j < c->get_nb_noeud(); j++)
    {
      grille.inserer(c->get_noeud(j));
      lstnoeud.ajouter(c->get_noeud(j));
    }
  }
  
  
  TPL_MAP_ENTITE<Poly_Noeud*>::ITERATEUR it;
  
  for (Poly_Noeud* nd=lstnoeud.get_premier(it);nd!=NULL;nd=lstnoeud.get_suivant(it))
    if (nd->est_fusionne()==false)
        {
        double x=nd->get_x();
        double y=nd->get_y();
        double z=nd->get_z();
        TPL_MAP_ENTITE<Poly_Noeud*> liste_entite_trouve;
        grille.rechercher(x,y,z,eps,liste_entite_trouve);
        TPL_MAP_ENTITE<Poly_Noeud*>::ITERATEUR it2;
        for (Poly_Noeud* nd2=liste_entite_trouve.get_premier(it2);nd2!=NULL;nd2=liste_entite_trouve.get_suivant(it2))
                if (nd2->est_fusionne()==false)
                  {
                  double x2=nd2->get_x();
                  double y2=nd2->get_y();
                  double z2=nd2->get_z();
                  if (fabs(x-x2)<eps)
                    if (fabs(y-y2)<eps)
                      if (fabs(z-z2)<eps)
                          {
                          nd2->set_x(x);  
                          nd2->set_y(y);  
                          nd2->set_z(z);  
                          nd2->active_fusion();
                          nd2->change_maitre_fusion(nd);
                          }
                    
                  }
        
        nd->active_fusion();
        nd->change_maitre_fusion(nd);
        }
}


//
// Accessor et modifier
//

void Poly_Voro::add_cell(Poly_Cellule* element)
{ 
  list_cellules.push_back(element); 
}

Poly_Point* Poly_Voro::get_point(int num)
{ 
  return list_points[num]; 
}

Poly_Cellule* Poly_Voro::get_cell(int num)
{ 
  return list_cellules[num]; 
}

int Poly_Voro::get_nb_cell(void)
{ 
  return list_cellules.size();
}


Revision:	1007
Committed:	Mon Mar 25 16:36:48 2019 UTC (6 years, 1 month ago) by francois
File size:	5137 byte(s)
Log Message:	fusion controlee dans le generateur de polycristal
#	User	Rev	Content
1	francois	979	#include "poly_voro.h"
2			#include "poly_point.h"
3			#include "poly_noeud.h"
4			#include "poly_face.h"
5			#include "poly_cellule.h"
6
7	francois	1007
8			#include "tpl_grille.h"
9
10
11	francois	979	#include <algorithm>
12
13			#include "../voro++-0.4.6/src/voro++.hh"
14
15	francois	1007	Poly_Voro::Poly_Voro(std::vector<Poly_Point*> list_pnts,double dg): POLY_AFFICHE(),list_points(list_pnts),epsfusion(dg)
16			{
17			}
18
19			void Poly_Voro::construit(void)
20			{
21	francois	979	// Géométrie du conteneur
22			double x_min= list_points[0]->get_x(), x_max= list_points[6]->get_x();
23			double y_min= list_points[0]->get_y(), y_max= list_points[6]->get_y();
24			double z_min= list_points[0]->get_z(), z_max= list_points[6]->get_z();
25			double cvol = (x_max - x_min) * (y_max - y_min) * (x_max - x_min);
26
27			// Nombre de blocs qui divise le container (blocs calculatoires)
28			int n_x= 6, n_y= 6, n_z= 6;
29
30			// Cré un container avec la géométrie données
31			voro::container con(x_min, x_max, y_min, y_max, z_min, z_max, n_x, n_y, n_z, false, false, false, 8);
32
33			// Ajoute les particules au container, on n'insère pas les 8 coins du cube
34			double x,y,z;
35			for(int i=8; i < list_points.size(); i++) {
36			x = list_points[i]->get_x();
37			y = list_points[i]->get_y();
38			z = list_points[i]->get_z();
39			con.put(i, x, y, z);
40			}
41
42			// Additionne les volumes et vérifi que ça concorde avec le volume du container
43			double vvol = con.sum_cell_volumes();
44	francois	1007	char message[2000];
45			sprintf(message," Volume de l'echantillon = %lf",cvol);
46			affiche(message);
47			sprintf(message," Volume des cellules de Voronoi = %lf",vvol);
48			affiche(message);
49
50
51	francois	979	voro::voronoicell voro_cell;
52			voro::c_loop_all loop(con);
53
54			// loop a travers les cellules
55			if(loop.start()) do if(con.compute_cell(voro_cell,loop))
56			{
57			Poly_Cellule* cell = new Poly_Cellule;
58			add_cell(cell);
59			loop.pos(x,y,z); // position globale de la particule en cours
60
61			// récupération des noeuds
62			std::vector< double > vert_xyz;
63			voro_cell.vertices(x, y, z, vert_xyz); // récupère tous les vertex de la cellule dans le système global
64
65			for(int i=0; i < vert_xyz.size(); i+=3)
66			{
67	francois	1007	static int numeronoeud=0;
68			numeronoeud++;
69			Poly_Noeud* nd = new Poly_Noeud(numeronoeud, vert_xyz[i], vert_xyz[i+1], vert_xyz[i+2]);
70	francois	979	cell->add_noeud(nd);
71			}
72
73			// récupération des faces
74			std::vector< int > face_vert;
75			voro_cell.face_vertices(face_vert); // liste de vertex pour chaque face
76
77			for(int i=0; i < face_vert.size(); i++)
78			{
79			Poly_Face* face = new Poly_Face;
80			int nbVert = face_vert[i]; // première valeur = le nombre de vertex de la face
81			for(int j=i+1; j < i+nbVert+1; j++)
82			{
83			face->add_noeud(face_vert[j]);
84			}
85			cell->add_face(face);
86			i += nbVert;
87			}
88			} while (loop.inc());
89			}
90
91			// destructor
92			Poly_Voro::~Poly_Voro()
93			{
94			for (int i=0; i<get_nb_cell(); i++)
95			{
96			Poly_Cellule* c = list_cellules[i];
97			delete c;
98			}
99			}
100
101			void Poly_Voro::fusion_noeuds(void)
102			{
103	francois	1007	TPL_MAP_ENTITE<Poly_Noeud*> lstnoeud;
104			TPL_GRILLE<Poly_Noeud*> grille;
105			double x_min= list_points[0]->get_x(), x_max= list_points[6]->get_x();
106			double y_min= list_points[0]->get_y(), y_max= list_points[6]->get_y();
107			double z_min= list_points[0]->get_z(), z_max= list_points[6]->get_z();
108			BOITE_3D boiteenglobante(x_min,y_min,z_min,x_max,y_max,z_max);
109			boiteenglobante.change_grosseur(1.1);
110			int nb=pow(list_points.size(),0.333333333333333333)+1;
111			grille.initialiser(boiteenglobante.get_xmin(),boiteenglobante.get_ymin(),boiteenglobante.get_zmin(),boiteenglobante.get_xmax(),boiteenglobante.get_ymax(),boiteenglobante.get_zmax(),nb,nb,nb);
112	francois	979
113	francois	1007	double eps = epsfusion;
114	francois	979
115	francois	1007
116	francois	979	for(int i=0; i < get_nb_cell(); i++)
117			{
118			Poly_Cellule* c = get_cell(i);
119			for(int j=0; j < c->get_nb_noeud(); j++)
120			{
121	francois	1007	grille.inserer(c->get_noeud(j));
122			lstnoeud.ajouter(c->get_noeud(j));
123	francois	979	}
124			}
125
126
127	francois	1007	TPL_MAP_ENTITE<Poly_Noeud*>::ITERATEUR it;
128
129			for (Poly_Noeud* nd=lstnoeud.get_premier(it);nd!=NULL;nd=lstnoeud.get_suivant(it))
130			if (nd->est_fusionne()==false)
131			{
132			double x=nd->get_x();
133			double y=nd->get_y();
134			double z=nd->get_z();
135			TPL_MAP_ENTITE<Poly_Noeud*> liste_entite_trouve;
136			grille.rechercher(x,y,z,eps,liste_entite_trouve);
137			TPL_MAP_ENTITE<Poly_Noeud*>::ITERATEUR it2;
138			for (Poly_Noeud* nd2=liste_entite_trouve.get_premier(it2);nd2!=NULL;nd2=liste_entite_trouve.get_suivant(it2))
139			if (nd2->est_fusionne()==false)
140			{
141			double x2=nd2->get_x();
142			double y2=nd2->get_y();
143			double z2=nd2->get_z();
144			if (fabs(x-x2)<eps)
145			if (fabs(y-y2)<eps)
146			if (fabs(z-z2)<eps)
147			{
148			nd2->set_x(x);
149			nd2->set_y(y);
150			nd2->set_z(z);
151			nd2->active_fusion();
152			nd2->change_maitre_fusion(nd);
153			}
154
155			}
156
157			nd->active_fusion();
158			nd->change_maitre_fusion(nd);
159			}
160	francois	979	}
161
162
163			//
164			// Accessor et modifier
165			//
166
167			void Poly_Voro::add_cell(Poly_Cellule* element)
168			{
169			list_cellules.push_back(element);
170			}
171
172			Poly_Point* Poly_Voro::get_point(int num)
173			{
174			return list_points[num];
175			}
176
177			Poly_Cellule* Poly_Voro::get_cell(int num)
178			{
179			return list_cellules[num];
180			}
181
182			int Poly_Voro::get_nb_cell(void)
183			{
184			return list_cellules.size();
185	francois	1007	}
186
187