ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur_auto/src/remailleur.cpp
Revision: 210
Committed: Fri Jul 31 16:19:04 2009 UTC (15 years, 9 months ago) by francois
Original Path: magic/lib/mailleur/mailleur/src/remailleur.cpp
File size: 43589 byte(s)
Log Message:
Nouvelle facon de representer la deforme + groupe topologique pour le maillage + bug idmax

File Contents

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