ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/remailleur.cpp
Revision: 425
Committed: Tue Sep 24 22:01:46 2013 UTC (11 years, 11 months ago) by francois
File size: 54377 byte(s)
Log Message:
ajout d'une procedure d'analyse de qualité de maillage + amelioration de la sortie sur terminal des informations dans les mailleurs

File Contents

# Content
1 //---------------------------------------------------------------------------
2
3 #include "gestionversion.h"
4 #include <math.h>
5 #include "remailleur.h"
6 #include "mg_gestionnaire.h"
7 #include "mailleur3d.h"
8 #include "mailleur0d.h"
9 #include "mailleur1d.h"
10 #include "mailleur2d.h"
11 #include "ot_mathematique.h"
12 #include "mg_geometrie_outils.h"
13 //---------------------------------------------------------------------------
14 #pragma package(smart_init)
15
16
17
18
19
20 REMAILLEUR::REMAILLEUR(class MG_GESTIONNAIRE* g1,MG_GESTIONNAIRE* g2,MG_GEOMETRIE* geo1,MG_GEOMETRIE* geo2,class FCT_TAILLE* fct_taille1,FCT_TAILLE* fct_taille2,MG_MAILLAGE* mori,MG_MAILLAGE* mmodi,VCT_COMPARAISON_RESULTAT& cmp):gestorigine(g1),gestmodifie(g2),carteori(fct_taille1),cartemod(fct_taille2),rescmp(cmp),maiorigine(mori),maimodifie(mmodi),geoorigine(geo1),geomodifie(geo2)
21 {
22 int nbvolume=geomodifie->get_nb_mg_volume();
23 lsttrifront=new TPL_LISTE_ENTITE<MG_TRIANGLE*>[nbvolume];
24 }
25
26 REMAILLEUR::~REMAILLEUR()
27 {
28 //TPL_MAP_ENTITE<class MG_SEGMENT_FRONTIERE*>::ITERATEUR it;
29 //for (MG_SEGMENT_FRONTIERE* seg=lstsegfront.get_premier(it);seg!=NULL;seg=lstsegfront.get_suivant(it) ) delete seg;
30 delete octree;
31 delete [] lsttrifront;
32 }
33
34
35 void REMAILLEUR::maille(MG_GROUPE_TOPOLOGIQUE* mggt)
36 {
37 maille(20);
38 }
39
40 void REMAILLEUR::maille(int etape)
41 {
42 //etape 0 initialisation
43 initialise_octree();
44 //etape 1 Destruction des entites autour des entites disparues
45 if (etape<1) return;
46 affiche((char*)" Destruction autour des disparitions");
47 TPL_MAP_ENTITE<MG_ELEMENT_MAILLAGE*> lstdetruire;
48 int nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_DISPARUE);
49 for (int i=0;i<nb;i++)
50 {
51 MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_DISPARUE,i);
52 int dim=ele->get_dimension();
53 if (dim==0)
54 {
55 MG_SOMMET* som=(MG_SOMMET*)ele;
56 int nblien=som->get_lien_maillage()->get_nb();
57 for (int j=0;j<nblien;j++)
58 {
59 MG_NOEUD* no=(MG_NOEUD*)som->get_lien_maillage()->get(j);
60 if (maiorigine->get_mg_noeudid(no->get_id())!=NULL) lstdetruire.ajouter(no);
61 }
62 }
63 if (dim==1)
64 {
65 MG_ARETE* are=(MG_ARETE*)ele;
66 int nblien=are->get_lien_maillage()->get_nb();
67 for (int j=0;j<nblien;j++)
68 {
69 MG_SEGMENT* seg=(MG_SEGMENT*)are->get_lien_maillage()->get(j);
70 if (maiorigine->get_mg_segmentid(seg->get_id())!=NULL)
71 {
72 lstdetruire.ajouter(seg->get_noeud1());
73 lstdetruire.ajouter(seg->get_noeud2());
74 }
75 }
76 }
77 if (dim==2)
78 {
79 MG_FACE* face=(MG_FACE*)ele;
80 int nblien=face->get_lien_maillage()->get_nb();
81 for (int j=0;j<nblien;j++)
82 {
83 MG_TRIANGLE* tri=(MG_TRIANGLE*)face->get_lien_maillage()->get(j);
84 if (maiorigine->get_mg_triangleid(tri->get_id())!=NULL)
85 {
86 lstdetruire.ajouter(tri->get_noeud1());
87 lstdetruire.ajouter(tri->get_noeud2());
88 lstdetruire.ajouter(tri->get_noeud3());
89 }
90 }
91 }
92 }
93 TPL_MAP_ENTITE<MG_ELEMENT_MAILLAGE*>::ITERATEUR it;
94 for (MG_NOEUD* no=(MG_NOEUD*)lstdetruire.get_premier(it);no!=NULL;no=(MG_NOEUD*)lstdetruire.get_suivant(it))
95 {
96 double ereelle=0.;
97 int nbseg=no->get_lien_segment()->get_nb();
98 for (int i=0;i<nbseg;i++)
99 ereelle=ereelle+no->get_lien_segment()->get(i)->get_longueur();
100 ereelle=ereelle/nbseg;
101 detruit_noeud(no,2*ereelle);
102 }
103 //etape 2 : transfert des noeuds de sommets;
104 if (etape<2) return;
105 affiche((char*)" Transfert sommet");
106 nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE);
107 for (int i=0;i<nb;i++)
108 {
109 MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE,i);
110 int dim=ele->get_dimension();
111 if (dim==0)
112 {
113 MG_SOMMET* som=(MG_SOMMET*)ele;
114 int nblien=som->get_lien_maillage()->get_nb();
115 for (int j=0;j<nblien;j++)
116 {
117 MG_NOEUD* no=(MG_NOEUD*)som->get_lien_maillage()->get(j);
118 if ((maiorigine->get_mg_noeudid(no->get_id())!=NULL) && (no->get_nouveau_numero()==CONSERVE))
119 {
120 MG_SOMMET* nvsom=geomodifie->get_mg_sommetid(rescmp.get_liste_correspondance_modifie(som->get_id()));
121 transfert_noeud(no,nvsom);
122 }
123 }
124
125 }
126 }
127 // etape 3 maillage 0d
128 if (etape<3) return;
129 affiche((char*)" MAILLAGE 0D");
130 MAILLEUR0D m0d(maimodifie,geomodifie);
131 m0d.adapte();
132 // etape 4 destruction autour des nouveaux sommets
133 if (etape<4) return;
134 affiche((char*)" Destruction autour des nouveaux sommets");
135 lstdetruire.vide();
136 nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::MODIFIE_APPARUE);
137 for (int i=0;i<nb;i++)
138 {
139 MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::MODIFIE_APPARUE,i);
140 int dim=ele->get_dimension();
141 if (dim==0)
142 {
143 MG_SOMMET* som=(MG_SOMMET*)ele;
144 int nblien=som->get_lien_maillage()->get_nb();
145 for (int j=0;j<nblien;j++)
146 {
147 MG_NOEUD* no=(MG_NOEUD*)som->get_lien_maillage()->get(j);
148 if (maimodifie->get_mg_noeudid(no->get_id())!=NULL) lstdetruire.ajouter(no);
149 }
150 }
151 }
152 for (MG_ELEMENT_MAILLAGE* ele=lstdetruire.get_premier(it);ele!=NULL;ele=lstdetruire.get_suivant(it))
153 {
154 detruit_noeud(ele,-1,1);
155 }
156 // etape 5 transferts des segments d arete
157 if (etape<5) return;
158 /*affiche(" Transfert arete");
159 nb=rescmp.get_nb_liste_topologie(ORIGINE_CONSERVEE);
160 for (int i=0;i<nb;i++)
161 {
162 MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(ORIGINE_CONSERVEE,i);
163 int dim=ele->get_dimension();
164 if (dim==1)
165 {
166 MG_ARETE* are=(MG_ARETE*)ele;
167 int nblien=are->get_lien_maillage()->get_nb();
168 for (int j=0;j<nblien;j++)
169 {
170 MG_SEGMENT* seg=(MG_SEGMENT*)are->get_lien_maillage()->get(j);
171 if ((maiorigine->get_mg_segmentid(seg->get_id())!=NULL))
172 {
173 MG_NOEUD* n1=seg->get_noeud1();
174 MG_NOEUD* n2=seg->get_noeud2();
175 if (n1->get_nouveau_numero()==CONSERVE)
176 if (n2->get_nouveau_numero()==CONSERVE)
177 {
178 MG_ARETE* nvare=geomodifie->get_mg_areteid(rescmp.get_liste_correspondance(are->get_id()));
179 MG_NOEUD* nv1=get_liste_correspondance(n1);
180 MG_NOEUD* nv2=get_liste_correspondance(n2);
181 if (nv1==NULL) nv1=transfert_noeud(n1,nvare);
182 if (nv2==NULL) nv2=transfert_noeud(n2,nvare);
183 MG_SEGMENT *seg=new MG_SEGMENT(nvare,nv1,nv2,IMPOSE);
184 maimodifie->ajouter_mg_segment(seg);
185 }
186 }
187 }
188
189 }
190 }
191 */
192 affiche((char*)" Creation des blocs de maille 1D");
193 int nbarete=geomodifie->get_nb_mg_arete();
194 nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE);
195 for (int i=0;i<nb;i++)
196 {
197 MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE,i);
198 int dim=ele->get_dimension();
199 if (dim==1)
200 {
201 multimap<unsigned long,MG_SEGMENT*,less <unsigned long> > lstseg;
202 MG_ARETE* arete=(MG_ARETE*)ele;
203 TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR ittet;
204 for (MG_ELEMENT_MAILLAGE* ele=arete->get_lien_maillage()->get_premier(ittet);ele!=NULL;ele=arete->get_lien_maillage()->get_suivant(ittet))
205 {
206 MG_SEGMENT* seg=(MG_SEGMENT*)ele;
207 if ((maiorigine->get_mg_segmentid(seg->get_id())!=NULL))
208 {
209 MG_NOEUD* n1=seg->get_noeud1();
210 MG_NOEUD* n2=seg->get_noeud2();
211 if (n1->get_nouveau_numero()==CONSERVE)
212 if (n2->get_nouveau_numero()==CONSERVE)
213 {
214 seg->change_nouveau_numero(NONTRAITE);
215 pair<unsigned long,MG_SEGMENT*> tmp(n1->get_id(),seg);
216 lstseg.insert(tmp);
217 pair<unsigned long,MG_SEGMENT*> tmp2(n2->get_id(),seg);
218 lstseg.insert(tmp2);
219 }
220 }
221 }
222
223 cree_bloc_maille(lstseg,arete);
224 }
225 }
226 // etape 6 Maillage des aretes
227 if (etape<6) return;
228 affiche((char*)" MAILLAGE 1D");
229 int nbare=geomodifie->get_nb_mg_arete();
230 for (int i=0;i<nbare;i++)
231 {
232 MG_ARETE* arete=geomodifie->get_mg_arete(i);
233 TPL_LISTE_ENTITE<double> param;
234 int typeface=arete->get_courbe()->get_type_geometrique(param);
235 multimap<unsigned long,BLOC_MAILLE_1D*,less <unsigned long> >::iterator it=lstb1d.begin();
236 while (it!=lstb1d.end())
237 {
238 BLOC_MAILLE_1D *b1d=(*it).second;
239 if ((b1d->etat==NONATTACHE) && (b1d->type==typeface))
240 {
241 MG_NOEUD* no1=b1d->lst[0]->get_noeud1();
242 MG_NOEUD* no2=b1d->lst[0]->get_noeud2();
243 double xyz1[3],xyz2[3],xyz[3];
244 transfert_coord(no1->get_coord(),xyz1);
245 transfert_coord(no2->get_coord(),xyz2);
246 double t1,t2;
247 arete->inverser(t1,xyz1);
248 arete->inverser(t2,xyz2);
249 int sensidentique=1;
250 OT_VECTEUR_3D vec(xyz1,xyz2);
251 vec.norme();
252 double dxyz1[3];
253 arete->deriver(t1,dxyz1);
254 OT_VECTEUR_3D dvec(dxyz1);
255 dvec.norme();
256 double ps=dvec*vec;
257 if (ps <0) sensidentique=0;
258 if (sensidentique==0)
259 {
260 double t=t1;
261 t1=t2;
262 t2=t;
263 double xyz[3];
264 xyz[0]=xyz1[0];
265 xyz[1]=xyz1[1];
266 xyz[2]=xyz1[2];
267 xyz1[0]=xyz2[0];
268 xyz1[1]=xyz2[1];
269 xyz1[2]=xyz2[2];
270 xyz2[0]=xyz[0];
271 xyz2[1]=xyz[1];
272 xyz2[2]=xyz[2];
273 }
274 if (arete->get_courbe()->est_periodique()==1)
275 {
276 //if (t1<arete->get_tmin()) t1=t1+arete->get_courbe()->get_periode();
277 if (t2<t1+1e-6) t2=t2+arete->get_courbe()->get_periode();
278 }
279 double xyz1tmp[3],xyz2tmp[3];
280 arete->evaluer(t1,xyz1tmp);
281 arete->evaluer(t2,xyz2tmp);
282 OT_VECTEUR_3D vec1(xyz1,xyz1tmp);
283 OT_VECTEUR_3D vec2(xyz2,xyz2tmp);
284 double metrique[9];
285 cartemod->evaluer(xyz1,metrique);
286 double eps=(1./sqrt(metrique[0])*1e-6);
287 double t=0.5*(t1+t2);
288 if ((vec1.get_longueur()<eps)&&(vec2.get_longueur()<eps))
289 {
290 if (arete->get_courbe()->est_periodique()==1)
291 if (t<arete->get_tmin()) t=t+arete->get_courbe()->get_periode();
292 if ((t>=arete->get_tmin()) && (t<=arete->get_tmax()))
293 {
294 int nbseg=b1d->lst.size();
295 b1d->etat=ATTACHE;
296 for (int i=0;i<nbseg;i++)
297 {
298 MG_SEGMENT *seg=b1d->lst[i];
299 MG_NOEUD* no1=seg->get_noeud1();
300 MG_NOEUD* no2=seg->get_noeud2();
301 double xyz1[3],xy2[3];
302 transfert_coord(no1->get_coord(),xyz1);
303 transfert_coord(no2->get_coord(),xyz2);
304 MG_NOEUD* nv1=get_liste_correspondance(no1);
305 MG_NOEUD* nv2=get_liste_correspondance(no2);
306 if (nv1==NULL) nv1=transfert_noeud(no1,arete);
307 if (nv2==NULL) nv2=transfert_noeud(no2,arete);
308 MG_SEGMENT *nvseg;
309 if (sensidentique==1)
310 nvseg=new MG_SEGMENT(arete,nv1,nv2,IMPOSE);
311 else
312 nvseg=new MG_SEGMENT(arete,nv2,nv1,IMPOSE);
313 maimodifie->ajouter_mg_segment(nvseg);
314 }
315 }
316 }
317 }
318 it++;
319 }
320 MAILLEUR1D m1d(maimodifie,geomodifie,arete,cartemod);
321 m1d.adapte();
322 }
323 //cree_quadtree(maimodifie,quad,&lstsegfront);
324 // etape 7 Destruction autour des nouveaux segments
325 if (etape<7) return;
326 affiche((char*)" Destruction autour des nouvelles aretes");
327 lstdetruire.vide();
328 /*nb=rescmp.get_nb_liste_topologie(MODIFIE_APPARUE);
329 for (int i=0;i<nb;i++)
330 {
331 MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(MODIFIE_APPARUE,i);
332 int dim=ele->get_dimension();
333 if (dim==1)
334 {
335 MG_ARETE* are=(MG_ARETE*)ele;
336 int nblien=are->get_lien_maillage()->get_nb();
337 for (int j=0;j<nblien;j++)
338 {
339 MG_SEGMENT* seg=(MG_SEGMENT*)are->get_lien_maillage()->get(j);
340 if (maimodifie->get_mg_segmentid(seg->get_id())!=NULL)
341 {
342 lstnoeuddetruire.ajouter(seg->get_noeud1());
343 lstnoeuddetruire.ajouter(seg->get_noeud2());
344 }
345 }
346 }
347 }*/
348 LISTE_MG_SEGMENT::iterator itnvseg;
349 for (MG_SEGMENT* seg=maimodifie->get_premier_segment(itnvseg);seg!=NULL;seg=maimodifie->get_suivant_segment(itnvseg))
350 {
351 if (seg->get_origine()==MAILLEUR_AUTO)
352 {
353 lstdetruire.ajouter(seg);
354 }
355 }
356 for (MG_ELEMENT_MAILLAGE* ele=lstdetruire.get_premier(it);ele!=NULL;ele=lstdetruire.get_suivant(it))
357 {
358 detruit_noeud(ele,-1,2);
359 }
360
361 // etape 8 Creation des blocs de mailles 2D
362 if (etape<8) return;
363 affiche((char*)" Creation des blocs de maille 2D");
364 nb=rescmp.get_nb_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE);
365 for (int i=0;i<nb;i++)
366 {
367 MG_ELEMENT_TOPOLOGIQUE* ele=rescmp.get_liste_topologie(VCT_COMPARAISON_RESULTAT::ORIGINE_CONSERVEE,i);
368 int dim=ele->get_dimension();
369 if (dim==2)
370 {
371 multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> > lsttri;
372 MG_FACE* face=(MG_FACE*)ele;
373 TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR ittet;
374 int compteur=0,compteur2=0;
375 for (MG_ELEMENT_MAILLAGE* ele=face->get_lien_maillage()->get_premier(ittet);ele!=NULL;ele=face->get_lien_maillage()->get_suivant(ittet))
376 {
377 MG_TRIANGLE* tri=(MG_TRIANGLE*)ele;
378 if ((maiorigine->get_mg_triangleid(tri->get_id())!=NULL))
379 {
380 MG_NOEUD* n1=tri->get_noeud1();
381 MG_NOEUD* n2=tri->get_noeud2();
382 MG_NOEUD* n3=tri->get_noeud3();
383 compteur2++;
384 if (n1->get_nouveau_numero()==CONSERVE)
385 if (n2->get_nouveau_numero()==CONSERVE)
386 if (n3->get_nouveau_numero()==CONSERVE)
387 {
388 tri->change_nouveau_numero(NONTRAITE);
389 pair<unsigned long,MG_TRIANGLE*> tmp(n1->get_id(),tri);
390 lsttri.insert(tmp);
391 pair<unsigned long,MG_TRIANGLE*> tmp2(n2->get_id(),tri);
392 lsttri.insert(tmp2);
393 pair<unsigned long,MG_TRIANGLE*> tmp3(n3->get_id(),tri);
394 lsttri.insert(tmp3);
395 compteur++;
396 }
397 }
398 }
399
400 cree_bloc_maille(lsttri,face);
401 }
402 }
403 // etape 9 MAILLAGE 2D
404 if (etape<9) return;
405 affiche((char*)" MAILLAGE 2D");
406 nbfacemod=geomodifie->get_nb_mg_face();
407 for (int i=0;i<nbfacemod;i++)
408 {
409 MG_GEOMETRIE_OUTILS ot;
410 TPL_LISTE_ENTITE<MG_TRIANGLE*> lsttricontraint;
411 MG_FACE* face=geomodifie->get_mg_face(i);
412 TPL_LISTE_ENTITE<double> param;
413 int typeface=face->get_surface()->get_type_geometrique(param);
414 multimap<unsigned long,BLOC_MAILLE_2D*,less <unsigned long> >::iterator it=lstb2d.begin();
415 while (it!=lstb2d.end())
416 {
417 BLOC_MAILLE_2D *b2d=(*it).second;
418 if ((b2d->etat==NONATTACHE) && (b2d->type==typeface))
419 {
420 MG_NOEUD* no1=b2d->lst[0]->get_noeud1();
421 MG_NOEUD* no2=b2d->lst[0]->get_noeud2();
422 MG_NOEUD* no3=b2d->lst[0]->get_noeud3();
423 double uv[2],uv1[2],uv2[2],uv3[2];
424 double xyz1[3],xyz2[3],xyz3[3];
425 transfert_coord(no1->get_coord(),xyz1);
426 transfert_coord(no2->get_coord(),xyz2);
427 transfert_coord(no3->get_coord(),xyz3);
428 face->inverser(uv1,xyz1);
429 face->inverser(uv2,xyz2);
430 face->inverser(uv3,xyz3);
431 if (face->get_surface()->est_periodique_u()==1)
432 {
433 double eps=1e-10*face->get_surface()->get_periode_u();
434 if (uv1[0]<-eps) uv1[0]=uv1[0]+face->get_surface()->get_periode_u();
435 if (uv1[0]>face->get_surface()->get_periode_u()-eps) uv1[0]=uv1[0]-face->get_surface()->get_periode_u();
436 if (uv2[0]<-eps) uv2[0]=uv2[0]+face->get_surface()->get_periode_u();
437 if (uv2[0]>face->get_surface()->get_periode_u()-eps) uv2[0]=uv2[0]-face->get_surface()->get_periode_u();
438 if (uv3[0]<-eps) uv3[0]=uv3[0]+face->get_surface()->get_periode_u();
439 if (uv3[0]>face->get_surface()->get_periode_u()-eps) uv3[0]=uv3[0]-face->get_surface()->get_periode_u();
440 }
441 if (face->get_surface()->est_periodique_v()==1)
442 {
443 double eps=1e-10*face->get_surface()->get_periode_v();
444 if (uv1[1]<-eps) uv1[1]=uv1[1]+face->get_surface()->get_periode_v();
445 if (uv1[1]>face->get_surface()->get_periode_v()-eps) uv1[1]=uv1[1]-face->get_surface()->get_periode_v();
446 if (uv2[1]<-eps) uv2[1]=uv2[1]+face->get_surface()->get_periode_v();
447 if (uv2[1]>face->get_surface()->get_periode_v()-eps) uv2[1]=uv2[1]-face->get_surface()->get_periode_v();
448 if (uv3[1]<-eps) uv3[1]=uv3[1]+face->get_surface()->get_periode_v();
449 if (uv3[1]>face->get_surface()->get_periode_v()-eps) uv3[1]=uv3[1]-face->get_surface()->get_periode_v();
450 }
451
452 double xyz1tmp[3];
453 face->evaluer(uv1,xyz1tmp);
454 double xyz2tmp[3];
455 face->evaluer(uv2,xyz2tmp);
456 double xyz3tmp[3];
457 face->evaluer(uv3,xyz3tmp);
458 OT_VECTEUR_3D vec1(xyz1,xyz1tmp);
459 OT_VECTEUR_3D vec2(xyz2,xyz2tmp);
460 OT_VECTEUR_3D vec3(xyz3,xyz3tmp);
461 double metrique[9];
462 cartemod->evaluer(xyz1,metrique);
463 double eps=(1./sqrt(metrique[0])*1e-6);
464 if (vec1.get_longueur()<eps)
465 if (vec2.get_longueur()<eps)
466 if (vec3.get_longueur()<eps)
467 {
468 OT_DECALAGE_PARAMETRE deca(face->get_surface()->get_periode_u(),face->get_surface()->get_periode_v());
469 double du=deca.calcul_decalage_parametre_u(uv1[0]);
470 double dv=deca.calcul_decalage_parametre_v(uv1[1]);
471 uv[0]=0.333333333333333333*(deca.decalage_parametre_u(uv1[0],du)+deca.decalage_parametre_u(uv2[0],du)+deca.decalage_parametre_u(uv3[0],du));
472 uv[1]=0.333333333333333333*(deca.decalage_parametre_v(uv1[1],dv)+deca.decalage_parametre_v(uv2[1],dv)+deca.decalage_parametre_v(uv3[1],dv));
473 uv[0]=deca.decalage_parametre_u(uv[0],-du);
474 uv[1]=deca.decalage_parametre_v(uv[1],-dv);
475 double dist=ot.calcule_distance_contour_face_uv(uv,face);
476 if (dist>0.)
477 {
478 int sensidentique=1;
479 OT_VECTEUR_3D vec12(xyz1,xyz2);
480 OT_VECTEUR_3D vec13(xyz1,xyz3);
481 double normal[3];
482 face->calcul_normale_unitaire(uv1,normal);
483 OT_VECTEUR_3D nor=vec13 & vec12;
484 nor.norme();
485 if ((nor*normal)<0.) sensidentique=0;
486 int nbtri=b2d->lst.size();
487 b2d->etat=ATTACHE;
488 //cout << " attachment du bloc " << b2d->face->get_id();
489 for (int i=0;i<nbtri;i++)
490 {
491 MG_TRIANGLE *tri=b2d->lst[i];
492 MG_NOEUD* no1=tri->get_noeud1();
493 MG_NOEUD* no2=tri->get_noeud2();
494 MG_NOEUD* no3=tri->get_noeud3();
495 double xyz1[3],xy2[3],xyz3[3];
496 transfert_coord(no1->get_coord(),xyz1);
497 transfert_coord(no2->get_coord(),xyz2);
498 transfert_coord(no3->get_coord(),xyz3);
499 MG_NOEUD* nvno1=new MG_NOEUD(NULL,xyz1[0],xyz1[1],xyz1[2],IMPOSE);
500 MG_NOEUD* nvno2=new MG_NOEUD(NULL,xyz2[0],xyz2[1],xyz2[2],IMPOSE);
501 MG_NOEUD* nvno3=new MG_NOEUD(NULL,xyz3[0],xyz3[1],xyz3[2],IMPOSE);
502 MG_TRIANGLE* nvtri;
503 if (sensidentique==1)
504 nvtri=new MG_TRIANGLE(NULL,nvno1,nvno2,nvno3,NULL,NULL,NULL,IMPOSE);
505 else
506 nvtri=new MG_TRIANGLE(NULL,nvno2,nvno1,nvno3,NULL,NULL,NULL,IMPOSE);
507 lsttricontraint.ajouter(nvtri);
508 nvno1->change_nouveau_numero(no1->get_id());
509 nvno2->change_nouveau_numero(no2->get_id());
510 nvno3->change_nouveau_numero(no3->get_id());
511 }
512 }
513 }
514 }
515 it++;
516 }
517 MAILLEUR2D m2d(maimodifie,geomodifie,face,cartemod);
518 char mess[50];
519 sprintf(mess," Face %i identificateur %lu",i,face->get_id());
520 affiche(mess);
521 m2d.maille(face,NULL,&(lsttricontraint));
522 int nbtricon=lsttricontraint.get_nb();
523 for (int i=0;i<nbtricon;i++)
524 {
525 MG_TRIANGLE* tri=lsttricontraint.get(0);
526 MG_NOEUD* noeud1=tri->get_noeud1();
527 MG_NOEUD* noeud2=tri->get_noeud2();
528 MG_NOEUD* noeud3=tri->get_noeud3();
529 lsttricontraint.supprimer(tri);
530 delete tri;
531 delete noeud1;
532 delete noeud2;
533 delete noeud3;
534 }
535 }
536 multimap<unsigned long,BLOC_MAILLE_2D*,less <unsigned long> >::iterator itfin2=lstb2d.begin();
537 while (itfin2!=lstb2d.end())
538 {
539 delete (*itfin2).second;
540 itfin2++;
541 }
542 cree_liste_frontiere(maimodifie);
543 // etape 10 Destruction autour des nouveaux triangles
544 if (etape<10) return;
545 affiche((char*)" Destruction autour des nouvelles faces");
546 lstdetruire.vide();
547 LISTE_MG_TRIANGLE::iterator itnvtri;
548 for (MG_TRIANGLE* tri=maimodifie->get_premier_triangle(itnvtri);tri!=NULL;tri=maimodifie->get_suivant_triangle(itnvtri))
549 {
550 if (tri->get_origine()==MAILLEUR_AUTO)
551 {
552 lstdetruire.ajouter(tri);
553 }
554 }
555 for (MG_ELEMENT_MAILLAGE* ele=lstdetruire.get_premier(it);ele!=NULL;ele=lstdetruire.get_suivant(it))
556 {
557 detruit_noeud(ele,-1,3);
558 }
559 // etape 11 Creation des blocs de mailles 3D
560 if (etape<11) return;
561 affiche((char*)" Creation des blocs de maille 3D");
562 int nbvolume=geoorigine->get_nb_mg_volume();
563 for (int i=0;i<nbvolume;i++)
564 {
565 multimap<unsigned long,MG_TETRA*,less <unsigned long> > lsttet;
566 MG_VOLUME* vol=geoorigine->get_mg_volume(i);
567 int nblien=vol->get_lien_maillage()->get_nb();
568 TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR ittet;
569 for (MG_ELEMENT_MAILLAGE* ele=vol->get_lien_maillage()->get_premier(ittet);ele!=NULL;ele=vol->get_lien_maillage()->get_suivant(ittet))
570 {
571 MG_TETRA* tet=(MG_TETRA*)ele;
572 if ((maiorigine->get_mg_tetraid(tet->get_id())!=NULL))
573 {
574 MG_NOEUD* n1=tet->get_noeud1();
575 MG_NOEUD* n2=tet->get_noeud2();
576 MG_NOEUD* n3=tet->get_noeud3();
577 MG_NOEUD* n4=tet->get_noeud4();
578 if (n1->get_nouveau_numero()==CONSERVE)
579 if (n2->get_nouveau_numero()==CONSERVE)
580 if (n3->get_nouveau_numero()==CONSERVE)
581 if (n4->get_nouveau_numero()==CONSERVE)
582 {
583 tet->change_nouveau_numero(NONTRAITE);
584 pair<unsigned long,MG_TETRA*> tmp(n1->get_id(),tet);
585 lsttet.insert(tmp);
586 pair<unsigned long,MG_TETRA*> tmp2(n2->get_id(),tet);
587 lsttet.insert(tmp2);
588 pair<unsigned long,MG_TETRA*> tmp3(n3->get_id(),tet);
589 lsttet.insert(tmp3);
590 pair<unsigned long,MG_TETRA*> tmp4(n4->get_id(),tet);
591 lsttet.insert(tmp4);
592 }
593 }
594
595 }
596 cree_bloc_maille(lsttet,vol);
597 }
598 // etape 12 Maillage 3D
599 if (etape<12) return;
600 affiche((char*)" MAILLEUR 3D");
601 nbvolume=geomodifie->get_nb_mg_volume();
602 for (int i=0;i<nbvolume;i++)
603 {
604 TPL_LISTE_ENTITE<MG_TETRA*> lsttetcontraint;
605 MG_VOLUME* vol=geomodifie->get_mg_volume(i);
606 multimap<unsigned long,BLOC_MAILLE_3D*,less <unsigned long> >::iterator it=lstb3d.begin();
607 while (it!=lstb3d.end())
608 {
609 BLOC_MAILLE_3D *b3d=(*it).second;
610 if (b3d->etat==NONATTACHE)
611 {
612 MG_NOEUD* no1=b3d->lst[0]->get_noeud1();
613 MG_NOEUD* no2=b3d->lst[0]->get_noeud2();
614 MG_NOEUD* no3=b3d->lst[0]->get_noeud3();
615 MG_NOEUD* no4=b3d->lst[0]->get_noeud4();
616 double xyz1[3],xyz2[3],xyz3[3],xyz4[3],xyz[3];
617 transfert_coord(no1->get_coord(),xyz1);
618 transfert_coord(no2->get_coord(),xyz2);
619 transfert_coord(no3->get_coord(),xyz3);
620 transfert_coord(no4->get_coord(),xyz4);
621 xyz[0]=0.25*(xyz1[0]+xyz2[0]+xyz3[0]+xyz4[0]);
622 xyz[1]=0.25*(xyz1[1]+xyz2[1]+xyz3[1]+xyz4[1]);
623 xyz[2]=0.25*(xyz1[2]+xyz2[2]+xyz3[2]+xyz4[2]);
624 int res=point_appartient_volume(xyz,i);
625 if (res==1)
626 {
627 int nbtet=b3d->lst.size();
628 b3d->etat=ATTACHE;
629 //cout << "bloc recupere " << b3d->get_id() << endl;
630 for (int j=0;j<nbtet;j++)
631 {
632 MG_TETRA *tet=b3d->lst[j];
633 MG_NOEUD* no1=tet->get_noeud1();
634 MG_NOEUD* no2=tet->get_noeud2();
635 MG_NOEUD* no3=tet->get_noeud3();
636 MG_NOEUD* no4=tet->get_noeud4();
637 double xyz1[3],xy2[3],xyz3[3],xyz4[3];
638 transfert_coord(no1->get_coord(),xyz1);
639 transfert_coord(no2->get_coord(),xyz2);
640 transfert_coord(no3->get_coord(),xyz3);
641 transfert_coord(no4->get_coord(),xyz4);
642 MG_NOEUD* nvno1=new MG_NOEUD(NULL,xyz1[0],xyz1[1],xyz1[2],IMPOSE);
643 MG_NOEUD* nvno2=new MG_NOEUD(NULL,xyz2[0],xyz2[1],xyz2[2],IMPOSE);
644 MG_NOEUD* nvno3=new MG_NOEUD(NULL,xyz3[0],xyz3[1],xyz3[2],IMPOSE);
645 MG_NOEUD* nvno4=new MG_NOEUD(NULL,xyz4[0],xyz4[1],xyz4[2],IMPOSE);
646 MG_TETRA* nvtet=new MG_TETRA(NULL,nvno1,nvno2,nvno3,nvno4,NULL,NULL,NULL,NULL,IMPOSE);
647 lsttetcontraint.ajouter(nvtet);
648 nvno1->change_nouveau_numero(no1->get_id());
649 nvno2->change_nouveau_numero(no2->get_id());
650 nvno3->change_nouveau_numero(no3->get_id());
651 nvno4->change_nouveau_numero(no4->get_id());
652 nvtet->change_nouveau_numero(tet->get_id());
653 }
654 }
655 }
656 it++;
657 }
658 MAILLEUR3D m3d(maimodifie,geomodifie,vol,cartemod);
659 char mess[50];
660 sprintf(mess," Volume %i identificateur %lu",i,vol->get_id());
661 affiche(mess);
662 m3d.active_affichage(affiche2);
663 m3d.maille(vol,NULL,&(lsttetcontraint));
664 int nbtetcon=lsttetcontraint.get_nb();
665 for (int i=0;i<nbtetcon;i++)
666 {
667 MG_TETRA* tet=lsttetcontraint.get(0);
668 MG_NOEUD* noeud1=tet->get_noeud1();
669 MG_NOEUD* noeud2=tet->get_noeud2();
670 MG_NOEUD* noeud3=tet->get_noeud3();
671 MG_NOEUD* noeud4=tet->get_noeud4();
672 lsttetcontraint.supprimer(tet);
673 delete tet;
674 delete noeud1;
675 delete noeud2;
676 delete noeud3;
677 delete noeud4;
678 }
679 }
680 multimap<unsigned long,BLOC_MAILLE_3D*,less <unsigned long> >::iterator itfin3=lstb3d.begin();
681 while (itfin3!=lstb3d.end())
682 {
683 delete (*itfin3).second;
684 itfin3++;
685 }
686 // etape 13 Recuperation des correspondances
687 if (etape<13) return;
688 affiche((char*)" Sauvegarde des correspondances");
689 LISTE_MG_NOEUD::iterator itnoeud;
690 for (MG_NOEUD* noeud=maimodifie->get_premier_noeud(itnoeud);noeud!=NULL;noeud=maimodifie->get_suivant_noeud(itnoeud))
691 {
692 if (noeud->get_origine()==IMPOSE)
693 rescmp.ajouter_correspondance(noeud->get_nouveau_numero(),noeud->get_id());
694 }
695 LISTE_MG_TETRA::iterator ittet;
696 for (MG_TETRA* tet=maimodifie->get_premier_tetra(ittet);tet!=NULL;tet=maimodifie->get_suivant_tetra(ittet))
697 {
698 if (tet->get_origine()==IMPOSE)
699 rescmp.ajouter_correspondance(tet->get_nouveau_numero(),tet->get_id());
700 }
701
702 }
703
704
705 void REMAILLEUR::initialise_octree(void)
706 {
707 FCT_GENERATEUR_3D<4> *carte2=(FCT_GENERATEUR_3D<4> *)carteori;
708 int nb_cellule=carte2->get_nb_cellule();
709 double param[32];
710 carte2->get_information(0,0,param);
711 double xmin=param[0]+0.005*(param[1]-param[0]);
712 double ymin=param[8]+0.005*(param[10]-param[8]);
713 double zmin=param[16]+0.005*(param[20]-param[16]);
714 carte2->get_information(nb_cellule-1,0,param);
715 double xmax=param[1]-0.005*(param[1]-param[0]);
716 double ymax=param[10]-0.005*(param[10]-param[8]);
717 double zmax=param[20]-0.005*(param[20]-param[16]);
718 octree=new TPL_OCTREE_FCT<MG_NOEUD*,FCT_GENERATEUR_3D<4> >;
719 octree->initialiser(*carte2,xmin,ymin,zmin,xmax,ymax,zmax);
720 LISTE_MG_NOEUD::iterator it;
721 for (MG_NOEUD* noeud=maiorigine->get_premier_noeud(it);noeud!=NULL;noeud=maiorigine->get_suivant_noeud(it))
722 {
723 noeud->change_nouveau_numero(CONSERVE);
724 octree->inserer(noeud);
725 }
726
727 }
728
729 void REMAILLEUR::detruit_noeud(MG_ELEMENT_MAILLAGE* elebase,double distance,int type)
730 {
731 if (elebase->get_origine()!=MAILLEUR_AUTO) return;
732 double xyz[3];
733 double distmin;
734 OT_VECTEUR_3D direction;
735 if (type==0)
736 {
737 MG_NOEUD* no=(MG_NOEUD*)elebase;
738 OT_VECTEUR_3D vecnoeud(no->get_coord());
739 //OT_VECTEUR_3D vecnoeud2=rescmp.change_coord_mod_ori(vecnoeud);
740 OT_VECTEUR_3D vecnoeud2=vecnoeud;
741 xyz[0]=vecnoeud2.get_x();
742 xyz[1]=vecnoeud2.get_y();
743 xyz[2]=vecnoeud2.get_z();
744 }
745 if (type==1)
746 {
747 MG_NOEUD* no=(MG_NOEUD*)elebase;
748 OT_VECTEUR_3D vecnoeud(no->get_coord());
749 OT_VECTEUR_3D vecnoeud2=rescmp.change_coord_mod_ori(vecnoeud);
750 xyz[0]=vecnoeud2.get_x();
751 xyz[1]=vecnoeud2.get_y();
752 xyz[2]=vecnoeud2.get_z();
753 }
754 if (type==2)
755 {
756 MG_SEGMENT* seg=(MG_SEGMENT*)elebase;
757 OT_VECTEUR_3D vecnoeud1(seg->get_noeud1()->get_coord());
758 OT_VECTEUR_3D vecnoeud2(seg->get_noeud2()->get_coord());
759 OT_VECTEUR_3D vecnoeud=0.5*(vecnoeud1+vecnoeud2);
760 OT_VECTEUR_3D vecnoeud3=rescmp.change_coord_mod_ori(vecnoeud);
761 xyz[0]=vecnoeud3.get_x();
762 xyz[1]=vecnoeud3.get_y();
763 xyz[2]=vecnoeud3.get_z();
764 direction=vecnoeud2-vecnoeud1;
765 distmin=direction.get_longueur()/2.;
766 direction.norme();
767 }
768 if (type==3)
769 {
770 MG_TRIANGLE* tri=(MG_TRIANGLE*)elebase;
771 OT_VECTEUR_3D vecnoeud1(tri->get_noeud1()->get_coord());
772 OT_VECTEUR_3D vecnoeud2(tri->get_noeud2()->get_coord());
773 OT_VECTEUR_3D vecnoeud3(tri->get_noeud3()->get_coord());
774 OT_VECTEUR_3D vecnoeud=0.3333333333333333*(vecnoeud1+vecnoeud2+vecnoeud3);
775 OT_VECTEUR_3D vecnoeud4=rescmp.change_coord_mod_ori(vecnoeud);
776 xyz[0]=vecnoeud4.get_x();
777 xyz[1]=vecnoeud4.get_y();
778 xyz[2]=vecnoeud4.get_z();
779 distmin=std::max((vecnoeud-vecnoeud1).get_longueur(),(vecnoeud-vecnoeud2).get_longueur());
780 distmin=std::max((vecnoeud-vecnoeud3).get_longueur(),distmin);
781 direction=(vecnoeud3-vecnoeud1)&(vecnoeud2-vecnoeud1);
782 direction.norme();
783 }
784 TPL_MAP_ENTITE<MG_NOEUD*> lst;
785 double Ea,Eb;
786 if (type!=0)
787 {
788 double m[9];
789 carteori->evaluer(xyz,m);
790 Ea=1./sqrt(m[0]);
791 }
792 octree->rechercher(xyz[0],xyz[1],xyz[2],2*Ea,lst);
793 double alpha=1.05;
794 distmin=distmin*alpha;
795 LISTE_MG_NOEUD::iterator it;
796 for (MG_NOEUD* noeud=lst.get_premier(it);noeud!=NULL;noeud=lst.get_suivant(it))
797 {
798 double *xyztmp=noeud->get_coord();
799 OT_VECTEUR_3D vec(xyz,xyztmp);
800 if (type!=0)
801 {
802 double m[9];
803 carteori->evaluer(xyztmp,m);
804 Eb=1./sqrt(m[0]);
805 distance=std::max(1.5*Ea,1.5*Eb);
806 distance=std::max(distance,2.6*fabs(Ea-Eb));
807 }
808 if (type==2)
809 {
810 double cosinus=vec*direction/vec.get_longueur();
811 double sinus=sqrt(1-cosinus*cosinus);
812 double facteur;
813 if (OPERATEUR::egal(cosinus,0.,0.000001)) facteur=distance;
814 else
815 {
816 double delta=sinus*sinus+4.*distance*distance/distmin/distmin*cosinus*cosinus;
817 facteur=(-sinus+sqrt(delta))/2./distance/cosinus/cosinus*distmin*distmin;
818 distance=facteur;
819 }
820 }
821 if (type==3)
822 {
823 double sinus=vec*direction/vec.get_longueur();
824 if (sinus>-0.0000001)
825 {
826 double cosinus=sqrt(1-sinus*sinus);
827 double facteur;
828 if (OPERATEUR::egal(cosinus,0.,0.000001)) facteur=distance;
829 else
830 {
831 double delta=sinus*sinus+4*distance*distance/distmin/distmin*cosinus*cosinus;
832 facteur=(-sinus+sqrt(delta))/2./distance/cosinus/cosinus*distmin*distmin;
833 distance=facteur;
834 }
835 }
836
837 }
838 if (vec.get_longueur()<distance)
839 {
840 noeud->change_nouveau_numero(DETRUIT);
841 }
842
843 }
844
845 }
846
847
848 /*void REMAILLEUR::cree_quadtree(MG_MAILLAGE* mai,TPL_QUADTREE<MG_SEGMENT_FRONTIERE*,MG_NOEUD*> **quad,TPL_MAP_ENTITE<MG_SEGMENT_FRONTIERE*> *lstsegfront)
849 {
850 int nbface=mai->get_mg_geometrie()->get_nb_mg_face();
851 for (int i=0;i<nbface;i++)
852 {
853 ((quad[i]))=new TPL_QUADTREE<MG_SEGMENT_FRONTIERE*,MG_NOEUD*>;
854 TPL_LISTE_ENTITE<MG_NOEUD*> lstn;
855 double umax=-1e300;
856 double vmax=-1e300;
857 double umin=1e300;
858 double vmin=1e300;
859 MG_FACE* face=mai->get_mg_geometrie()->get_mg_face(i);
860 int nbboucle=face->get_nb_mg_boucle();
861 int nbsegfravant=lstsegfront->get_nb();
862 for (int j=0;j<nbboucle;j++)
863 {
864 MG_BOUCLE* bou=face->get_mg_boucle(j);
865 int nbarete=bou->get_nb_mg_coarete();
866 for (int k=0;k<nbarete;k++)
867 {
868 MG_ARETE* are=bou->get_mg_coarete(k)->get_arete();
869 int nbseg=are->get_lien_maillage()->get_nb();
870 for (int l=0;l<nbseg;l++)
871 {
872 MG_SEGMENT* seg=(MG_SEGMENT*)are->get_lien_maillage()->get(l);
873 if (mai->get_mg_segmentid(seg->get_id())==NULL) continue;
874 double *xyz=seg->get_noeud1()->get_coord();
875 double uv[2];
876 face->inverser(uv,xyz);
877 seg->get_noeud1()->change_u(uv[0]);
878 seg->get_noeud1()->change_v(uv[1]);
879 xyz=seg->get_noeud2()->get_coord();
880 face->inverser(uv,xyz);
881 seg->get_noeud2()->change_u(uv[0]);
882 seg->get_noeud2()->change_v(uv[1]);
883 lstn.ajouter(seg->get_noeud2());
884 if (uv[0]<umin) umin=uv[0];
885 if (uv[0]>umax) umax=uv[0];
886 if (uv[1]<vmin) vmin=uv[1];
887 if (uv[1]>vmax) vmax=uv[1];
888 MG_SEGMENT_FRONTIERE* segfr=new MG_SEGMENT_FRONTIERE(seg);
889 lstsegfront->ajouter(segfr);
890 }
891 }
892 }
893 double periodeenu=face->get_surface()->get_periode_u();
894 double periodeenv=face->get_surface()->get_periode_v();
895 if (periodeenu!=0.0)
896 {
897 umin=0.;
898 umax=periodeenu;
899 }
900 else
901 {
902 double diff=umax-umin;
903 umin=umin-0.125*diff;
904 umax=umax+0.125*diff;
905 }
906 if (periodeenv!=0.0)
907 {
908 vmin=0.;
909 vmax=periodeenv;
910 }
911 else
912 {
913 double diff=vmax-vmin;
914 vmin=vmin-0.125*diff;
915 vmax=vmax+0.125*diff;
916 }
917 (quad[i])->initialiser(&lstn,1,umin,vmin,umax,vmax,face->get_surface()->get_periode_u(),face->get_surface()->get_periode_v());
918 int nbsegfr=lstsegfront->get_nb();
919 for (int j=nbsegfravant;j<nbsegfr;j++)
920 {
921 double *uv=lstsegfront->get(j)->get_uv1();
922 uv[0]=lstsegfront->get(j)->get_segment()->get_noeud1()->get_u();
923 uv[1]=lstsegfront->get(j)->get_segment()->get_noeud1()->get_v();
924 uv=lstsegfront->get(j)->get_uv2();
925 uv[0]=lstsegfront->get(j)->get_segment()->get_noeud2()->get_u();
926 uv[1]=lstsegfront->get(j)->get_segment()->get_noeud2()->get_v();
927 (quad[i])->inserer( lstsegfront->get(j));
928 }
929 }
930
931 }
932 */
933 void REMAILLEUR::cree_liste_frontiere(MG_MAILLAGE* mai)
934 {
935 int nbvolume=mai->get_mg_geometrie()->get_nb_mg_volume();
936 for (int i=0;i<nbvolume;i++)
937 {
938 MG_VOLUME* vol=mai->get_mg_geometrie()->get_mg_volume(i);
939 int nbco=vol->get_nb_mg_coquille();
940 for (int j=0;j<nbco;j++)
941 {
942 MG_COQUILLE* coq=vol->get_mg_coquille(j);
943 int nbface=coq->get_nb_mg_coface();
944 for (int k=0;k<nbface;k++)
945 {
946 MG_FACE* face=coq->get_mg_coface(k)->get_face();
947 int nbseg=face->get_lien_maillage()->get_nb();
948 for (int l=0;l<nbseg;l++)
949 {
950 MG_TRIANGLE* tri=(MG_TRIANGLE*)face->get_lien_maillage()->get(l);
951 if (mai->get_mg_triangleid(tri->get_id())==NULL) continue;
952 lsttrifront[i].ajouter(tri);
953 }
954 }
955 }
956
957
958 }
959 }
960
961
962
963 void REMAILLEUR::ajouter_liste_correspondance_noeud(MG_NOEUD* n1,MG_NOEUD* n2)
964 {
965 CORRESPONDANCENOEUD corr(n1,n2);
966 pair<unsigned long,CORRESPONDANCENOEUD> tmp(n1->get_id(),corr);
967 lstcorrnoeud.insert(tmp);
968 }
969
970 MG_NOEUD* REMAILLEUR::get_liste_correspondance(MG_NOEUD* n1)
971 {
972 MG_NOEUD *n2=NULL;
973 std::map<unsigned long,CORRESPONDANCENOEUD,std::less<unsigned long> >::iterator j=lstcorrnoeud.find(n1->get_id());
974 if (j!=lstcorrnoeud.end())
975 n2=(*j).second.noeudmod;
976 return n2;
977 }
978
979
980 MG_NOEUD* REMAILLEUR::transfert_noeud(MG_NOEUD* no,MG_ELEMENT_TOPOLOGIQUE* ele)
981 {
982 double *xyz=no->get_coord();
983 OT_VECTEUR_3D pt(xyz);
984 OT_VECTEUR_3D nvpt=rescmp.change_coord_ori_mod(pt);
985 MG_NOEUD* noeud=new MG_NOEUD(ele,nvpt.get_x(),nvpt.get_y(),nvpt.get_z(),IMPOSE);
986 maimodifie->ajouter_mg_noeud(noeud);
987 ajouter_liste_correspondance_noeud(no,noeud);
988 noeud->change_nouveau_numero(no->get_id());
989 return noeud;
990 }
991
992 void REMAILLEUR::transfert_coord(double *xyz,double *xyz2)
993 {
994 OT_VECTEUR_3D pt(xyz);
995 OT_VECTEUR_3D nvpt=rescmp.change_coord_ori_mod(pt);
996 xyz2[0]=nvpt.get_x();
997 xyz2[1]=nvpt.get_y();
998 xyz2[2]=nvpt.get_z();
999 return;
1000 }
1001
1002 void REMAILLEUR::cree_bloc_maille(multimap<unsigned long,MG_SEGMENT*,less <unsigned long> > &lst,MG_ARETE* arete)
1003 {
1004 while (lst.size()>0)
1005 {
1006 multimap<unsigned long,MG_SEGMENT*,less <unsigned long> >::iterator it=lst.begin();
1007 MG_SEGMENT* seg=(*it).second;
1008 lst.erase(it);
1009 TPL_LISTE_ENTITE<double> param;
1010 BLOC_MAILLE_1D *b1d=new BLOC_MAILLE_1D(arete,arete->get_courbe()->get_type_geometrique(param));
1011 pair<unsigned long,BLOC_MAILLE_1D*> tmp(b1d->get_id(),b1d);
1012 lstb1d.insert(tmp);
1013 b1d->lst.push_back(seg);
1014 seg->change_nouveau_numero(TRAITE);
1015 char mess[255];
1016 sprintf(mess," Arete %lu , Bloc %d",arete->get_id(),lstb1d.size());
1017 //affiche(mess);
1018 int num=0;
1019 while (num!=b1d->lst.size())
1020 {
1021 MG_SEGMENT* itseg=b1d->lst[num++];
1022 MG_SEGMENT* seg=trouve_segment(lst,itseg->get_noeud1());
1023 while (seg!=NULL)
1024 {
1025 if (seg->get_nouveau_numero()==NONTRAITE)
1026 {
1027 b1d->lst.push_back(seg);
1028 seg->change_nouveau_numero(TRAITE);
1029 }
1030 seg=trouve_segment(lst,itseg->get_noeud1());
1031 }
1032 seg=trouve_segment(lst,itseg->get_noeud2());
1033 while (seg!=NULL)
1034 {
1035 if (seg->get_nouveau_numero()==NONTRAITE)
1036 {
1037 b1d->lst.push_back(seg);
1038 seg->change_nouveau_numero(TRAITE);
1039 }
1040 seg=trouve_segment(lst,itseg->get_noeud2());
1041 }
1042 }
1043
1044 sprintf(mess," Segment %d ",b1d->lst.size());
1045 //affiche(mess);
1046 }
1047 }
1048
1049 void REMAILLEUR::cree_bloc_maille(multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> > &lst,MG_FACE* face)
1050 {
1051 while (lst.size()>0)
1052 {
1053 multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> >::iterator it=lst.begin();
1054 MG_TRIANGLE* tri=(*it).second;
1055 lst.erase(it);
1056 TPL_LISTE_ENTITE<double> param;
1057 BLOC_MAILLE_2D *b2d=new BLOC_MAILLE_2D(geoorigine->get_mg_faceid(face->get_id()),face->get_surface()->get_type_geometrique(param));
1058 pair<unsigned long,BLOC_MAILLE_2D*> tmp(b2d->get_id(),b2d);
1059 lstb2d.insert(tmp);
1060 b2d->lst.push_back(tri);
1061 tri->change_nouveau_numero(TRAITE);
1062 char mess[255];
1063 sprintf(mess," Face %lu , Bloc %d",face->get_id(),lstb2d.size());
1064 //affiche(mess);
1065 int num=0;
1066 while (num!=b2d->lst.size())
1067 {
1068 MG_TRIANGLE* ittri=b2d->lst[num++];
1069 MG_TRIANGLE* tri=trouve_triangle(lst,ittri->get_noeud1());
1070 while (tri!=NULL)
1071 {
1072 if (tri->get_nouveau_numero()==NONTRAITE)
1073 {
1074 b2d->lst.push_back(tri);
1075 tri->change_nouveau_numero(TRAITE);
1076 }
1077 tri=trouve_triangle(lst,ittri->get_noeud1());
1078 }
1079 tri=trouve_triangle(lst,ittri->get_noeud2());
1080 while (tri!=NULL)
1081 {
1082 if (tri->get_nouveau_numero()==NONTRAITE)
1083 {
1084 b2d->lst.push_back(tri);
1085 tri->change_nouveau_numero(TRAITE);
1086 }
1087 tri=trouve_triangle(lst,ittri->get_noeud2());
1088 }
1089 tri=trouve_triangle(lst,ittri->get_noeud3());
1090 while (tri!=NULL)
1091 {
1092 if (tri->get_nouveau_numero()==NONTRAITE)
1093 {
1094 b2d->lst.push_back(tri);
1095 tri->change_nouveau_numero(TRAITE);
1096 }
1097 tri=trouve_triangle(lst,ittri->get_noeud3());
1098 }
1099 }
1100
1101 sprintf(mess," Triangle %d ",b2d->lst.size());
1102 //affiche(mess);
1103 }
1104 }
1105
1106 void REMAILLEUR::cree_bloc_maille(multimap<unsigned long,MG_TETRA*,less <unsigned long> > &lst,MG_VOLUME* vol)
1107 {
1108 while (lst.size()>0)
1109 {
1110 multimap<unsigned long,MG_TETRA*,less <unsigned long> >::iterator it=lst.begin();
1111 MG_TETRA* tet=(*it).second;
1112 lst.erase(it);
1113 BLOC_MAILLE_3D *b3d=new BLOC_MAILLE_3D(geoorigine->get_mg_volumeid(vol->get_id()));
1114 pair<unsigned long,BLOC_MAILLE_3D*> tmp(b3d->get_id(),b3d);
1115 lstb3d.insert(tmp);
1116 b3d->lst.push_back(tet);
1117 tet->change_nouveau_numero(TRAITE);
1118 char mess[255];
1119 sprintf(mess," Volume %lu , Bloc %d",vol->get_id(),lstb3d.size());
1120 //affiche(mess);
1121 int num=0;
1122 while (num!=b3d->lst.size())
1123 {
1124 MG_TETRA* ittet=b3d->lst[num++];
1125 MG_TETRA* tet=trouve_tetra(lst,ittet->get_noeud1());
1126 while (tet!=NULL)
1127 {
1128 if (tet->get_nouveau_numero()==NONTRAITE)
1129 {
1130 b3d->lst.push_back(tet);
1131 tet->change_nouveau_numero(TRAITE);
1132 }
1133 tet=trouve_tetra(lst,ittet->get_noeud1());
1134 }
1135 tet=trouve_tetra(lst,ittet->get_noeud2());
1136 while (tet!=NULL)
1137 {
1138 if (tet->get_nouveau_numero()==NONTRAITE)
1139 {
1140 b3d->lst.push_back(tet);
1141 tet->change_nouveau_numero(TRAITE);
1142 }
1143 tet=trouve_tetra(lst,ittet->get_noeud2());
1144 }
1145 tet=trouve_tetra(lst,ittet->get_noeud3());
1146 while (tet!=NULL)
1147 {
1148 if (tet->get_nouveau_numero()==NONTRAITE)
1149 {
1150 b3d->lst.push_back(tet);
1151 tet->change_nouveau_numero(TRAITE);
1152 }
1153 tet=trouve_tetra(lst,ittet->get_noeud3());
1154 }
1155 tet=trouve_tetra(lst,ittet->get_noeud4());
1156 while (tet!=NULL)
1157 {
1158 if (tet->get_nouveau_numero()==NONTRAITE)
1159 {
1160 b3d->lst.push_back(tet);
1161 tet->change_nouveau_numero(TRAITE);
1162 }
1163 tet=trouve_tetra(lst,ittet->get_noeud4());
1164 }
1165 }
1166
1167 sprintf(mess," Tetra %d ",b3d->lst.size());
1168 //affiche(mess);
1169 }
1170 }
1171
1172
1173 MG_SEGMENT* REMAILLEUR::trouve_segment(multimap<unsigned long,MG_SEGMENT*,less <unsigned long> > &lst,MG_NOEUD* no)
1174 {
1175 MG_SEGMENT* seg=NULL;
1176 unsigned long id=no->get_id();
1177 multimap<unsigned long,MG_SEGMENT*,less <unsigned long> >::iterator it=lst.find(id);
1178 if (it!=lst.end())
1179 {
1180 seg=(*it).second;
1181 lst.erase(it);
1182 }
1183 return seg;
1184 }
1185
1186 MG_TRIANGLE* REMAILLEUR::trouve_triangle(multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> > &lst,MG_NOEUD* no)
1187 {
1188 MG_TRIANGLE* tri=NULL;
1189 unsigned long id=no->get_id();
1190 multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> >::iterator it=lst.find(id);
1191 if (it!=lst.end())
1192 {
1193 tri=(*it).second;
1194 lst.erase(it);
1195 }
1196 return tri;
1197 }
1198
1199 MG_TETRA* REMAILLEUR::trouve_tetra(multimap<unsigned long,MG_TETRA*,less <unsigned long> > &lst,MG_NOEUD* no)
1200 {
1201 MG_TETRA* tet=NULL;
1202 unsigned long id=no->get_id();
1203 multimap<unsigned long,MG_TETRA*,less <unsigned long> >::iterator it=lst.find(id);
1204 if (it!=lst.end())
1205 {
1206 tet=(*it).second;
1207 lst.erase(it);
1208 }
1209 return tet;
1210 }
1211
1212 int REMAILLEUR::point_appartient_volume(double *xyz,int numvol)
1213 {
1214 int ok=0;
1215 int i=0;
1216 multimap<double,double,less<double> > intersection;
1217 OT_VECTEUR_3D vecteur_dir;
1218 while (ok==0)
1219
1220 {
1221 int ok2=0;
1222 while (ok2==0)
1223 {
1224 MG_TRIANGLE* tri=lsttrifront[numvol].get(i);
1225 MG_NOEUD* no1=tri->get_noeud1();
1226 MG_NOEUD* no2=tri->get_noeud2();
1227 MG_NOEUD* no3=tri->get_noeud3();
1228 double xb=(no1->get_x()+no2->get_x()+no3->get_x())/3.;
1229 double yb=(no1->get_y()+no2->get_y()+no3->get_y())/3.;
1230 double zb=(no1->get_z()+no2->get_z()+no3->get_z())/3.;
1231 OT_VECTEUR_3D directeur(xb-xyz[0],yb-xyz[1],zb-xyz[2]);
1232 directeur.norme();
1233 vecteur_dir=directeur;
1234 OT_VECTEUR_3D n1n3(no1->get_coord(),no3->get_coord());
1235 OT_VECTEUR_3D n1n2(no1->get_coord(),no2->get_coord());
1236 OT_VECTEUR_3D n=n1n3&n1n2;
1237 n.norme();
1238 double ps1=n*directeur;
1239 if (fabs(ps1)<1e-6) i++;
1240 else ok2=1;
1241 }
1242 ok2=0;
1243 intersection.clear();
1244 for (int j=0;j<lsttrifront[numvol].get_nb();j++)
1245 {
1246 MG_TRIANGLE* tri=lsttrifront[numvol].get(j);
1247
1248 double t;
1249 int res=inter_droite_triangle(xyz,vecteur_dir,tri,&t);
1250 if (res==2) {
1251 ok2=1;
1252 break;
1253 }
1254 if (res==1)
1255 {
1256 pair<double,double> tmp(t,t);
1257 intersection.insert(tmp);
1258 }
1259 }
1260 if (ok2==0) ok=1;
1261 i++;
1262
1263 }
1264 multimap<double,double,less <double> >::iterator it=intersection.begin();
1265 int nbinterneg=0;
1266 double tavant=1e308;
1267 while (it!=intersection.end())
1268 {
1269 if ((*it).second<0.)
1270 {
1271 double t=(*it).second;
1272 if (!(OPERATEUR::egal(t,tavant,1e-6)))
1273 nbinterneg++;
1274 tavant=t;
1275 }
1276 it++;
1277 }
1278 if (nbinterneg%2==1) return 1;
1279 return 0;
1280 }
1281
1282
1283 int REMAILLEUR::inter_droite_triangle(double *xyz,double *dir,MG_TRIANGLE* tri,double *t)
1284 {
1285 int inter=0;
1286 MG_NOEUD *noeud1=tri->get_noeud1();
1287 MG_NOEUD *noeud2=tri->get_noeud2();
1288 MG_NOEUD *noeud3=tri->get_noeud3();
1289 OT_VECTEUR_3D g1(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
1290 OT_VECTEUR_3D g2(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
1291 OT_VECTEUR_3D g3(dir);
1292 OT_MATRICE_3D systeme(g1,g2,g3);
1293 double det=systeme.get_determinant();
1294 double eps=0.333333333333*(g1.diff()+g2.diff()+g3.diff());
1295 eps=eps*eps*0.0018;
1296 if (OPERATEUR::egal(det,0.0,eps)==true)
1297 {
1298 OT_VECTEUR_3D g3b(xyz[0]-noeud1->get_x(),xyz[1]-noeud1->get_y(),xyz[2]-noeud1->get_z());
1299 OT_MATRICE_3D systeme2(g1,g2,g3b);
1300 double det2=systeme2.get_determinant();
1301 double eps2=0.333333333333*(g1.diff()+g2.diff()+g3b.diff());
1302 eps2=18.*eps2*eps2*eps2*1e-6;
1303 if (OPERATEUR::egal(det2,0.0,eps)==true) // cas 2D
1304 {
1305 return 2;
1306 }
1307 else return 0;
1308 }
1309 else
1310 {
1311 double x,y,z;
1312 x=1.0/det*(g2.get_y()*g3.get_z()-g2.get_z()*g3.get_y());
1313 y=1.0/det*(g3.get_x()*g2.get_z()-g2.get_x()*g3.get_z());
1314 z=1.0/det*(g2.get_x()*g3.get_y()-g2.get_y()*g3.get_x());
1315 OT_VECTEUR_3D g1b(x,y,z);
1316 x=1.0/det*(g3.get_y()*g1.get_z()-g1.get_y()*g3.get_z());
1317 y=1.0/det*(g1.get_x()*g3.get_z()-g3.get_x()*g1.get_z());
1318 z=1.0/det*(g3.get_x()*g1.get_y()-g1.get_x()*g3.get_y());
1319 OT_VECTEUR_3D g2b(x,y,z);
1320 x=1.0/det*(g1.get_y()*g2.get_z()-g1.get_z()*g2.get_y());
1321 y=1.0/det*(g2.get_x()*g1.get_z()-g1.get_x()*g2.get_z());
1322 z=1.0/det*(g1.get_x()*g2.get_y()-g1.get_y()*g2.get_x());
1323 OT_VECTEUR_3D g3b(x,y,z);
1324 OT_VECTEUR_3D n1n4(xyz[0]-noeud1->get_x(),xyz[1]-noeud1->get_y(),xyz[2]-noeud1->get_z());
1325 double alpha1=n1n4*g1b;
1326 double alpha2=n1n4*g2b;
1327 double alpha3=-(n1n4*g3b);
1328 double alpha4=1-alpha1-alpha2;
1329 double eps=0.000001;
1330 if ((alpha1>-eps) && (alpha1<1.+eps))
1331 if ((alpha2>-eps) && (alpha2<1.+eps))
1332 if ((alpha4>-eps) && (alpha4<1.+eps)) {
1333 *t=alpha3;
1334 return 1;
1335 }
1336 return 0;
1337 }
1338
1339 }