ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/remailleur.cpp
Revision: 258
Committed: Thu Aug 12 19:10:34 2010 UTC (15 years ago) by francois
File size: 44892 byte(s)
Log Message:
Mise a jour toxfem + parametrisation compilation toxfem + bug 
comparaison

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