ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/remailleur.cpp
Revision: 174
Committed: Wed Apr 22 21:46:09 2009 UTC (16 years ago) by francois
Original Path: magic/lib/mailleur/mailleur/src/remailleur.cpp
File size: 43264 byte(s)
Log Message:
Plusieurs modif importantes : 
1) la fonction un point appartient a une face est generique donc le mailleurbloc et le remailleur utilise la methode generic
2) dans une boucle il y une methode qui permet de connaitre la coarete suivante et precedente
3) Des solutions sont possibles aux elements. Attention le format de fichier magic est change pour SOLUTION. Il y a un parametre different. Si il y a des solutions dans le fichier il n y a pas compatibilite avec le format d'avant

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     void REMAILLEUR::maille(void)
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 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 147 }
592     }
593     }
594     it++;
595     }
596     MAILLEUR3D m3d(maimodifie,geomodifie,vol,cartemod);
597     char mess[50];
598     sprintf(mess," Volume %i identificateur %lu",i,vol->get_id());
599     affiche(mess);
600     m3d.active_affichage(affiche);
601     m3d.maille(vol,NULL,&(lsttetcontraint));
602     int nbtetcon=lsttetcontraint.get_nb();
603     for (int i=0;i<nbtetcon;i++)
604     {
605     MG_TETRA* tet=lsttetcontraint.get(0);
606     MG_NOEUD* noeud1=tet->get_noeud1();
607     MG_NOEUD* noeud2=tet->get_noeud2();
608     MG_NOEUD* noeud3=tet->get_noeud3();
609     MG_NOEUD* noeud4=tet->get_noeud4();
610     lsttetcontraint.supprimer(tet);
611     delete tet;
612     delete noeud1;
613     delete noeud2;
614     delete noeud3;
615     delete noeud4;
616     }
617     }
618     multimap<unsigned long,BLOC_MAILLE_3D*,less <unsigned long> >::iterator itfin3=lstb3d.begin();
619     while (itfin3!=lstb3d.end())
620     {
621     delete (*itfin3).second;
622     itfin3++;
623     }
624 francois 155 // etape 13 Recuperation des correspondances
625     if (etape<13) return;
626     affiche(" Sauvegarde des correspondances");
627     LISTE_MG_NOEUD::iterator itnoeud;
628     for (MG_NOEUD* noeud=maimodifie->get_premier_noeud(itnoeud);noeud!=NULL;noeud=maimodifie->get_suivant_noeud(itnoeud))
629     {
630     if (noeud->get_origine()==IMPOSE)
631     rescmp.ajouter_correspondance_noeud(noeud->get_nouveau_numero(),noeud->get_id());
632     }
633 francois 147
634    
635     }
636    
637    
638     void REMAILLEUR::initialise_octree(void)
639     {
640     FCT_GENERATEUR_3D<4> *carte2=(FCT_GENERATEUR_3D<4> *)carteori;
641     int nb_cellule=carte2->get_nb_cellule();
642     double param[32];
643     carte2->get_information(0,0,param);
644     double xmin=param[0]+0.005*(param[1]-param[0]);
645     double ymin=param[8]+0.005*(param[10]-param[8]);
646     double zmin=param[16]+0.005*(param[20]-param[16]);
647     carte2->get_information(nb_cellule-1,0,param);
648     double xmax=param[1]-0.005*(param[1]-param[0]);
649     double ymax=param[10]-0.005*(param[10]-param[8]);
650     double zmax=param[20]-0.005*(param[20]-param[16]);
651     octree=new TPL_OCTREE_FCT<MG_NOEUD*,FCT_GENERATEUR_3D<4> >;
652     octree->initialiser(*carte2,xmin,ymin,zmin,xmax,ymax,zmax);
653     LISTE_MG_NOEUD::iterator it;
654     for (MG_NOEUD* noeud=maiorigine->get_premier_noeud(it);noeud!=NULL;noeud=maiorigine->get_suivant_noeud(it))
655     {
656     noeud->change_nouveau_numero(CONSERVE);
657     octree->inserer(noeud);
658     }
659    
660     }
661    
662 francois 168 void REMAILLEUR::detruit_noeud(MG_ELEMENT_MAILLAGE* elebase,double distance,int type)
663 francois 147 {
664 francois 168 if (elebase->get_origine()!=MAILLEUR_AUTO) return;
665 francois 155 double xyz[3];
666 francois 168 double distmin;
667     OT_VECTEUR_3D direction;
668     if (type==0)
669     {
670     MG_NOEUD* no=(MG_NOEUD*)elebase;
671     OT_VECTEUR_3D vecnoeud(no->get_coord());
672     OT_VECTEUR_3D vecnoeud2=rescmp.get_coord_inverse(vecnoeud);
673     xyz[0]=vecnoeud2.get_x();
674     xyz[1]=vecnoeud2.get_y();
675     xyz[2]=vecnoeud2.get_z();
676     }
677     if (type==1)
678     {
679     MG_NOEUD* no=(MG_NOEUD*)elebase;
680     OT_VECTEUR_3D vecnoeud(no->get_coord());
681     OT_VECTEUR_3D vecnoeud2=rescmp.get_coord_inverse(vecnoeud);
682     xyz[0]=vecnoeud2.get_x();
683     xyz[1]=vecnoeud2.get_y();
684     xyz[2]=vecnoeud2.get_z();
685     }
686     if (type==2)
687     {
688     MG_SEGMENT* seg=(MG_SEGMENT*)elebase;
689     OT_VECTEUR_3D vecnoeud1(seg->get_noeud1()->get_coord());
690     OT_VECTEUR_3D vecnoeud2(seg->get_noeud2()->get_coord());
691     OT_VECTEUR_3D vecnoeud=0.5*(vecnoeud1+vecnoeud2);
692     OT_VECTEUR_3D vecnoeud3=rescmp.get_coord_inverse(vecnoeud);
693     xyz[0]=vecnoeud3.get_x();
694     xyz[1]=vecnoeud3.get_y();
695     xyz[2]=vecnoeud3.get_z();
696     direction=vecnoeud2-vecnoeud1;
697     distmin=direction.get_longueur()/2.;
698     direction.norme();
699     }
700     if (type==3)
701     {
702     MG_TRIANGLE* tri=(MG_TRIANGLE*)elebase;
703     OT_VECTEUR_3D vecnoeud1(tri->get_noeud1()->get_coord());
704     OT_VECTEUR_3D vecnoeud2(tri->get_noeud2()->get_coord());
705     OT_VECTEUR_3D vecnoeud3(tri->get_noeud3()->get_coord());
706     OT_VECTEUR_3D vecnoeud=0.3333333333333333*(vecnoeud1+vecnoeud2+vecnoeud3);
707     OT_VECTEUR_3D vecnoeud4=rescmp.get_coord_inverse(vecnoeud);
708     xyz[0]=vecnoeud4.get_x();
709     xyz[1]=vecnoeud4.get_y();
710     xyz[2]=vecnoeud4.get_z();
711     distmin=std::max((vecnoeud-vecnoeud1).get_longueur(),(vecnoeud-vecnoeud2).get_longueur());
712     distmin=std::max((vecnoeud-vecnoeud3).get_longueur(),distmin);
713     direction=(vecnoeud3-vecnoeud1)&(vecnoeud2-vecnoeud1);
714     direction.norme();
715     }
716 francois 147 TPL_MAP_ENTITE<MG_NOEUD*> lst;
717     double Ea,Eb;
718     if (type!=0)
719     {
720     double m[9];
721     carteori->evaluer(xyz,m);
722 francois 168 Ea=1./sqrt(m[0]);
723 francois 147 }
724 francois 168 octree->rechercher(xyz[0],xyz[1],xyz[2],2*Ea,lst);
725     double alpha=1.05;
726     distmin=distmin*alpha;
727 francois 147 LISTE_MG_NOEUD::iterator it;
728     for (MG_NOEUD* noeud=lst.get_premier(it);noeud!=NULL;noeud=lst.get_suivant(it))
729     {
730     double *xyztmp=noeud->get_coord();
731     OT_VECTEUR_3D vec(xyz,xyztmp);
732     if (type!=0)
733     {
734     double m[9];
735     carteori->evaluer(xyztmp,m);
736 francois 168 Eb=1./sqrt(m[0]);
737     distance=std::max(1.5*Ea,1.5*Eb);
738     distance=std::max(distance,2.6*fabs(Ea-Eb));
739 francois 147 }
740 francois 155 if (type==2)
741 francois 147 {
742 francois 168 double cosinus=vec*direction/vec.get_longueur();
743     double sinus=sqrt(1-cosinus*cosinus);
744     double facteur;
745     if (OPERATEUR::egal(cosinus,0.,0.000001)) facteur=distance;
746     else
747 francois 155 {
748 francois 168 double delta=sinus*sinus+4.*distance*distance/distmin/distmin*cosinus*cosinus;
749     facteur=(-sinus+sqrt(delta))/2./distance/cosinus/cosinus*distmin*distmin;
750     distance=facteur;
751 francois 155 }
752 francois 168 }
753     if (type==3)
754     {
755     double sinus=vec*direction/vec.get_longueur();
756     if (sinus>-0.0000001)
757 francois 155 {
758 francois 168 double cosinus=sqrt(1-sinus*sinus);
759     double facteur;
760     if (OPERATEUR::egal(cosinus,0.,0.000001)) facteur=distance;
761     else
762     {
763     double delta=sinus*sinus+4*distance*distance/distmin/distmin*cosinus*cosinus;
764     facteur=(-sinus+sqrt(delta))/2./distance/cosinus/cosinus*distmin*distmin;
765     distance=facteur;
766     }
767 francois 155 }
768 francois 168
769 francois 155 }
770 francois 168 if (vec.get_longueur()<distance)
771 francois 155 {
772 francois 168 noeud->change_nouveau_numero(DETRUIT);
773 francois 147 }
774 francois 168
775 francois 147 }
776    
777     }
778    
779    
780 francois 174 /*void REMAILLEUR::cree_quadtree(MG_MAILLAGE* mai,TPL_QUADTREE<MG_SEGMENT_FRONTIERE*,MG_NOEUD*> **quad,TPL_MAP_ENTITE<MG_SEGMENT_FRONTIERE*> *lstsegfront)
781 francois 147 {
782     int nbface=mai->get_mg_geometrie()->get_nb_mg_face();
783     for (int i=0;i<nbface;i++)
784     {
785     ((quad[i]))=new TPL_QUADTREE<MG_SEGMENT_FRONTIERE*,MG_NOEUD*>;
786     TPL_LISTE_ENTITE<MG_NOEUD*> lstn;
787     double umax=-1e300;
788     double vmax=-1e300;
789     double umin=1e300;
790     double vmin=1e300;
791     MG_FACE* face=mai->get_mg_geometrie()->get_mg_face(i);
792     int nbboucle=face->get_nb_mg_boucle();
793     int nbsegfravant=lstsegfront->get_nb();
794     for (int j=0;j<nbboucle;j++)
795     {
796     MG_BOUCLE* bou=face->get_mg_boucle(j);
797     int nbarete=bou->get_nb_mg_coarete();
798     for (int k=0;k<nbarete;k++)
799     {
800     MG_ARETE* are=bou->get_mg_coarete(k)->get_arete();
801     int nbseg=are->get_lien_maillage()->get_nb();
802     for (int l=0;l<nbseg;l++)
803     {
804     MG_SEGMENT* seg=(MG_SEGMENT*)are->get_lien_maillage()->get(l);
805     if (mai->get_mg_segmentid(seg->get_id())==NULL) continue;
806     double *xyz=seg->get_noeud1()->get_coord();
807     double uv[2];
808     face->inverser(uv,xyz);
809     seg->get_noeud1()->change_u(uv[0]);
810     seg->get_noeud1()->change_v(uv[1]);
811     xyz=seg->get_noeud2()->get_coord();
812     face->inverser(uv,xyz);
813     seg->get_noeud2()->change_u(uv[0]);
814     seg->get_noeud2()->change_v(uv[1]);
815     lstn.ajouter(seg->get_noeud2());
816     if (uv[0]<umin) umin=uv[0];
817     if (uv[0]>umax) umax=uv[0];
818     if (uv[1]<vmin) vmin=uv[1];
819     if (uv[1]>vmax) vmax=uv[1];
820     MG_SEGMENT_FRONTIERE* segfr=new MG_SEGMENT_FRONTIERE(seg);
821     lstsegfront->ajouter(segfr);
822     }
823     }
824     }
825     double periodeenu=face->get_surface()->get_periode_u();
826     double periodeenv=face->get_surface()->get_periode_v();
827     if (periodeenu!=0.0)
828     {
829     umin=0.;
830     umax=periodeenu;
831     }
832     else
833     {
834     double diff=umax-umin;
835     umin=umin-0.125*diff;
836     umax=umax+0.125*diff;
837     }
838     if (periodeenv!=0.0)
839     {
840     vmin=0.;
841     vmax=periodeenv;
842     }
843     else
844     {
845     double diff=vmax-vmin;
846     vmin=vmin-0.125*diff;
847     vmax=vmax+0.125*diff;
848     }
849     (quad[i])->initialiser(&lstn,1,umin,vmin,umax,vmax,face->get_surface()->get_periode_u(),face->get_surface()->get_periode_v());
850     int nbsegfr=lstsegfront->get_nb();
851     for (int j=nbsegfravant;j<nbsegfr;j++)
852     {
853     double *uv=lstsegfront->get(j)->get_uv1();
854     uv[0]=lstsegfront->get(j)->get_segment()->get_noeud1()->get_u();
855     uv[1]=lstsegfront->get(j)->get_segment()->get_noeud1()->get_v();
856     uv=lstsegfront->get(j)->get_uv2();
857     uv[0]=lstsegfront->get(j)->get_segment()->get_noeud2()->get_u();
858     uv[1]=lstsegfront->get(j)->get_segment()->get_noeud2()->get_v();
859     (quad[i])->inserer( lstsegfront->get(j));
860     }
861     }
862    
863     }
864 francois 174 */
865 francois 147 void REMAILLEUR::cree_liste_frontiere(MG_MAILLAGE* mai)
866     {
867     int nbvolume=mai->get_mg_geometrie()->get_nb_mg_volume();
868     for (int i=0;i<nbvolume;i++)
869     {
870     MG_VOLUME* vol=mai->get_mg_geometrie()->get_mg_volume(i);
871     int nbco=vol->get_nb_mg_coquille();
872     for (int j=0;j<nbco;j++)
873     {
874     MG_COQUILLE* coq=vol->get_mg_coquille(j);
875     int nbface=coq->get_nb_mg_coface();
876     for (int k=0;k<nbface;k++)
877     {
878     MG_FACE* face=coq->get_mg_coface(k)->get_face();
879     int nbseg=face->get_lien_maillage()->get_nb();
880     for (int l=0;l<nbseg;l++)
881     {
882     MG_TRIANGLE* tri=(MG_TRIANGLE*)face->get_lien_maillage()->get(l);
883     if (mai->get_mg_triangleid(tri->get_id())==NULL) continue;
884     lsttrifront[i].ajouter(tri);
885     }
886     }
887     }
888    
889    
890     }
891     }
892    
893    
894    
895     void REMAILLEUR::ajouter_liste_correspondance_noeud(MG_NOEUD* n1,MG_NOEUD* n2)
896     {
897     CORRESPONDANCENOEUD corr(n1,n2);
898     pair<unsigned long,CORRESPONDANCENOEUD> tmp(n1->get_id(),corr);
899     lstcorrnoeud.insert(tmp);
900     }
901    
902     MG_NOEUD* REMAILLEUR::get_liste_correspondance(MG_NOEUD* n1)
903     {
904     MG_NOEUD *n2=NULL;
905     std::map<unsigned long,CORRESPONDANCENOEUD,std::less<unsigned long> >::iterator j=lstcorrnoeud.find(n1->get_id());
906     if (j!=lstcorrnoeud.end())
907     n2=(*j).second.noeudmod;
908     return n2;
909     }
910    
911    
912     MG_NOEUD* REMAILLEUR::transfert_noeud(MG_NOEUD* no,MG_ELEMENT_TOPOLOGIQUE* ele)
913     {
914     double *xyz=no->get_coord();
915     OT_VECTEUR_3D pt(xyz);
916     OT_VECTEUR_3D nvpt=rescmp.get_coord(pt);
917     MG_NOEUD* noeud=new MG_NOEUD(ele,nvpt.get_x(),nvpt.get_y(),nvpt.get_z(),IMPOSE);
918     maimodifie->ajouter_mg_noeud(noeud);
919     ajouter_liste_correspondance_noeud(no,noeud);
920 francois 155 noeud->change_nouveau_numero(no->get_id());
921 francois 147 return noeud;
922     }
923    
924     void REMAILLEUR::transfert_coord(double *xyz,double *xyz2)
925     {
926     OT_VECTEUR_3D pt(xyz);
927     OT_VECTEUR_3D nvpt=rescmp.get_coord(pt);
928     xyz2[0]=nvpt.get_x();
929     xyz2[1]=nvpt.get_y();
930     xyz2[2]=nvpt.get_z();
931     return;
932     }
933    
934     void REMAILLEUR::cree_bloc_maille(multimap<unsigned long,MG_SEGMENT*,less <unsigned long> > &lst,MG_ARETE* arete)
935     {
936     while (lst.size()>0)
937     {
938     multimap<unsigned long,MG_SEGMENT*,less <unsigned long> >::iterator it=lst.begin();
939     MG_SEGMENT* seg=(*it).second;
940     lst.erase(it);
941     TPL_LISTE_ENTITE<double> param;
942     BLOC_MAILLE_1D *b1d=new BLOC_MAILLE_1D(arete,arete->get_courbe()->get_type_geometrique(param));
943     pair<unsigned long,BLOC_MAILLE_1D*> tmp(b1d->get_id(),b1d);
944     lstb1d.insert(tmp);
945     b1d->lst.push_back(seg);
946     seg->change_nouveau_numero(TRAITE);
947 francois 164 char mess[255];
948     sprintf(mess," Arete %lu , Bloc %d",arete->get_id(),lstb1d.size());
949 francois 147 //affiche(mess);
950     int num=0;
951     while (num!=b1d->lst.size())
952     {
953     MG_SEGMENT* itseg=b1d->lst[num++];
954     MG_SEGMENT* seg=trouve_segment(lst,itseg->get_noeud1());
955     while (seg!=NULL)
956     {
957     if (seg->get_nouveau_numero()==NONTRAITE)
958     {
959     b1d->lst.push_back(seg);
960     seg->change_nouveau_numero(TRAITE);
961     }
962     seg=trouve_segment(lst,itseg->get_noeud1());
963     }
964     seg=trouve_segment(lst,itseg->get_noeud2());
965     while (seg!=NULL)
966     {
967     if (seg->get_nouveau_numero()==NONTRAITE)
968     {
969     b1d->lst.push_back(seg);
970     seg->change_nouveau_numero(TRAITE);
971     }
972     seg=trouve_segment(lst,itseg->get_noeud2());
973     }
974     }
975    
976 francois 164 sprintf(mess," Segment %d ",b1d->lst.size());
977 francois 147 //affiche(mess);
978     }
979     }
980    
981     void REMAILLEUR::cree_bloc_maille(multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> > &lst,MG_FACE* face)
982     {
983     while (lst.size()>0)
984     {
985     multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> >::iterator it=lst.begin();
986     MG_TRIANGLE* tri=(*it).second;
987     lst.erase(it);
988     TPL_LISTE_ENTITE<double> param;
989     BLOC_MAILLE_2D *b2d=new BLOC_MAILLE_2D(geoorigine->get_mg_faceid(face->get_id()),face->get_surface()->get_type_geometrique(param));
990     pair<unsigned long,BLOC_MAILLE_2D*> tmp(b2d->get_id(),b2d);
991     lstb2d.insert(tmp);
992     b2d->lst.push_back(tri);
993     tri->change_nouveau_numero(TRAITE);
994 francois 164 char mess[255];
995     sprintf(mess," Face %lu , Bloc %d",face->get_id(),lstb2d.size());
996 francois 147 //affiche(mess);
997     int num=0;
998     while (num!=b2d->lst.size())
999     {
1000     MG_TRIANGLE* ittri=b2d->lst[num++];
1001     MG_TRIANGLE* tri=trouve_triangle(lst,ittri->get_noeud1());
1002     while (tri!=NULL)
1003     {
1004     if (tri->get_nouveau_numero()==NONTRAITE)
1005     {
1006     b2d->lst.push_back(tri);
1007     tri->change_nouveau_numero(TRAITE);
1008     }
1009     tri=trouve_triangle(lst,ittri->get_noeud1());
1010     }
1011     tri=trouve_triangle(lst,ittri->get_noeud2());
1012     while (tri!=NULL)
1013     {
1014     if (tri->get_nouveau_numero()==NONTRAITE)
1015     {
1016     b2d->lst.push_back(tri);
1017     tri->change_nouveau_numero(TRAITE);
1018     }
1019     tri=trouve_triangle(lst,ittri->get_noeud2());
1020     }
1021     tri=trouve_triangle(lst,ittri->get_noeud3());
1022     while (tri!=NULL)
1023     {
1024     if (tri->get_nouveau_numero()==NONTRAITE)
1025     {
1026     b2d->lst.push_back(tri);
1027     tri->change_nouveau_numero(TRAITE);
1028     }
1029     tri=trouve_triangle(lst,ittri->get_noeud3());
1030     }
1031     }
1032    
1033 francois 164 sprintf(mess," Triangle %d ",b2d->lst.size());
1034 francois 147 //affiche(mess);
1035     }
1036     }
1037    
1038     void REMAILLEUR::cree_bloc_maille(multimap<unsigned long,MG_TETRA*,less <unsigned long> > &lst,MG_VOLUME* vol)
1039     {
1040     while (lst.size()>0)
1041     {
1042     multimap<unsigned long,MG_TETRA*,less <unsigned long> >::iterator it=lst.begin();
1043     MG_TETRA* tet=(*it).second;
1044     lst.erase(it);
1045     BLOC_MAILLE_3D *b3d=new BLOC_MAILLE_3D(geoorigine->get_mg_volumeid(vol->get_id()));
1046     pair<unsigned long,BLOC_MAILLE_3D*> tmp(b3d->get_id(),b3d);
1047     lstb3d.insert(tmp);
1048     b3d->lst.push_back(tet);
1049     tet->change_nouveau_numero(TRAITE);
1050 francois 164 char mess[255];
1051     sprintf(mess," Volume %lu , Bloc %d",vol->get_id(),lstb3d.size());
1052 francois 147 //affiche(mess);
1053     int num=0;
1054     while (num!=b3d->lst.size())
1055     {
1056     MG_TETRA* ittet=b3d->lst[num++];
1057     MG_TETRA* tet=trouve_tetra(lst,ittet->get_noeud1());
1058     while (tet!=NULL)
1059     {
1060     if (tet->get_nouveau_numero()==NONTRAITE)
1061     {
1062     b3d->lst.push_back(tet);
1063     tet->change_nouveau_numero(TRAITE);
1064     }
1065     tet=trouve_tetra(lst,ittet->get_noeud1());
1066     }
1067     tet=trouve_tetra(lst,ittet->get_noeud2());
1068     while (tet!=NULL)
1069     {
1070     if (tet->get_nouveau_numero()==NONTRAITE)
1071     {
1072     b3d->lst.push_back(tet);
1073     tet->change_nouveau_numero(TRAITE);
1074     }
1075     tet=trouve_tetra(lst,ittet->get_noeud2());
1076     }
1077     tet=trouve_tetra(lst,ittet->get_noeud3());
1078     while (tet!=NULL)
1079     {
1080     if (tet->get_nouveau_numero()==NONTRAITE)
1081     {
1082     b3d->lst.push_back(tet);
1083     tet->change_nouveau_numero(TRAITE);
1084     }
1085     tet=trouve_tetra(lst,ittet->get_noeud3());
1086     }
1087     tet=trouve_tetra(lst,ittet->get_noeud4());
1088     while (tet!=NULL)
1089     {
1090     if (tet->get_nouveau_numero()==NONTRAITE)
1091     {
1092     b3d->lst.push_back(tet);
1093     tet->change_nouveau_numero(TRAITE);
1094     }
1095     tet=trouve_tetra(lst,ittet->get_noeud4());
1096     }
1097     }
1098    
1099 francois 164 sprintf(mess," Tetra %d ",b3d->lst.size());
1100 francois 147 //affiche(mess);
1101     }
1102     }
1103    
1104    
1105     MG_SEGMENT* REMAILLEUR::trouve_segment(multimap<unsigned long,MG_SEGMENT*,less <unsigned long> > &lst,MG_NOEUD* no)
1106     {
1107     MG_SEGMENT* seg=NULL;
1108     unsigned long id=no->get_id();
1109     multimap<unsigned long,MG_SEGMENT*,less <unsigned long> >::iterator it=lst.find(id);
1110     if (it!=lst.end())
1111     {
1112     seg=(*it).second;
1113     lst.erase(it);
1114     }
1115     return seg;
1116     }
1117    
1118     MG_TRIANGLE* REMAILLEUR::trouve_triangle(multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> > &lst,MG_NOEUD* no)
1119     {
1120     MG_TRIANGLE* tri=NULL;
1121     unsigned long id=no->get_id();
1122     multimap<unsigned long,MG_TRIANGLE*,less <unsigned long> >::iterator it=lst.find(id);
1123     if (it!=lst.end())
1124     {
1125     tri=(*it).second;
1126     lst.erase(it);
1127     }
1128     return tri;
1129     }
1130    
1131     MG_TETRA* REMAILLEUR::trouve_tetra(multimap<unsigned long,MG_TETRA*,less <unsigned long> > &lst,MG_NOEUD* no)
1132     {
1133     MG_TETRA* tet=NULL;
1134     unsigned long id=no->get_id();
1135     multimap<unsigned long,MG_TETRA*,less <unsigned long> >::iterator it=lst.find(id);
1136     if (it!=lst.end())
1137     {
1138     tet=(*it).second;
1139     lst.erase(it);
1140     }
1141     return tet;
1142     }
1143    
1144     int REMAILLEUR::point_appartient_volume(double *xyz,int numvol)
1145     {
1146     int ok=0;
1147     int i=0;
1148     multimap<double,double,less<double> > intersection;
1149     OT_VECTEUR_3D vecteur_dir;
1150     while (ok==0)
1151    
1152     {
1153     int ok2=0;
1154     while (ok2==0)
1155     {
1156     MG_TRIANGLE* tri=lsttrifront[numvol].get(i);
1157     MG_NOEUD* no1=tri->get_noeud1();
1158     MG_NOEUD* no2=tri->get_noeud2();
1159     MG_NOEUD* no3=tri->get_noeud3();
1160     double xb=(no1->get_x()+no2->get_x()+no3->get_x())/3.;
1161     double yb=(no1->get_y()+no2->get_y()+no3->get_y())/3.;
1162     double zb=(no1->get_z()+no2->get_z()+no3->get_z())/3.;
1163     OT_VECTEUR_3D directeur(xb-xyz[0],yb-xyz[1],zb-xyz[2]);
1164     directeur.norme();
1165     vecteur_dir=directeur;
1166     OT_VECTEUR_3D n1n3(no1->get_coord(),no3->get_coord());
1167     OT_VECTEUR_3D n1n2(no1->get_coord(),no2->get_coord());
1168     OT_VECTEUR_3D n=n1n3&n1n2;
1169     n.norme();
1170     double ps1=n*directeur;
1171     if (fabs(ps1)<1e-6) i++;
1172     else ok2=1;
1173     }
1174     ok2=0;
1175     intersection.clear();
1176     for (int j=0;j<lsttrifront[numvol].get_nb();j++)
1177     {
1178     MG_TRIANGLE* tri=lsttrifront[numvol].get(j);
1179    
1180     double t;
1181     int res=inter_droite_triangle(xyz,vecteur_dir,tri,&t);
1182     if (res==2) {ok2=1;break;}
1183     if (res==1)
1184     {
1185     pair<double,double> tmp(t,t);
1186     intersection.insert(tmp);
1187     }
1188     }
1189     if (ok2==0) ok=1;
1190     i++;
1191    
1192     }
1193     multimap<double,double,less <double> >::iterator it=intersection.begin();
1194     int nbinterneg=0;
1195 francois 167 double tavant=1e308;
1196 francois 147 while (it!=intersection.end())
1197     {
1198 francois 167 if ((*it).second<0.)
1199     {
1200     double t=(*it).second;
1201     if (!(OPERATEUR::egal(t,tavant,1e-6)))
1202     nbinterneg++;
1203     tavant=t;
1204     }
1205 francois 147 it++;
1206     }
1207     if (nbinterneg%2==1) return 1;
1208     return 0;
1209     }
1210    
1211    
1212     int REMAILLEUR::inter_droite_triangle(double *xyz,double *dir,MG_TRIANGLE* tri,double *t)
1213     {
1214     int inter=0;
1215     MG_NOEUD *noeud1=tri->get_noeud1();
1216     MG_NOEUD *noeud2=tri->get_noeud2();
1217     MG_NOEUD *noeud3=tri->get_noeud3();
1218     OT_VECTEUR_3D g1(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
1219     OT_VECTEUR_3D g2(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
1220     OT_VECTEUR_3D g3(dir);
1221     OT_MATRICE_3D systeme(g1,g2,g3);
1222     double det=systeme.get_determinant();
1223     double eps=0.333333333333*(g1.diff()+g2.diff()+g3.diff());
1224     eps=eps*eps*0.0018;
1225     if (OPERATEUR::egal(det,0.0,eps)==true)
1226     {
1227     OT_VECTEUR_3D g3b(xyz[0]-noeud1->get_x(),xyz[1]-noeud1->get_y(),xyz[2]-noeud1->get_z());
1228     OT_MATRICE_3D systeme2(g1,g2,g3b);
1229     double det2=systeme2.get_determinant();
1230     double eps2=0.333333333333*(g1.diff()+g2.diff()+g3b.diff());
1231     eps2=18.*eps2*eps2*eps2*1e-6;
1232     if (OPERATEUR::egal(det2,0.0,eps)==true) // cas 2D
1233     {
1234     return 2;
1235     }
1236     else return 0;
1237     }
1238     else
1239     {
1240     double x,y,z;
1241     x=1.0/det*(g2.get_y()*g3.get_z()-g2.get_z()*g3.get_y());
1242     y=1.0/det*(g3.get_x()*g2.get_z()-g2.get_x()*g3.get_z());
1243     z=1.0/det*(g2.get_x()*g3.get_y()-g2.get_y()*g3.get_x());
1244     OT_VECTEUR_3D g1b(x,y,z);
1245     x=1.0/det*(g3.get_y()*g1.get_z()-g1.get_y()*g3.get_z());
1246     y=1.0/det*(g1.get_x()*g3.get_z()-g3.get_x()*g1.get_z());
1247     z=1.0/det*(g3.get_x()*g1.get_y()-g1.get_x()*g3.get_y());
1248     OT_VECTEUR_3D g2b(x,y,z);
1249     x=1.0/det*(g1.get_y()*g2.get_z()-g1.get_z()*g2.get_y());
1250     y=1.0/det*(g2.get_x()*g1.get_z()-g1.get_x()*g2.get_z());
1251     z=1.0/det*(g1.get_x()*g2.get_y()-g1.get_y()*g2.get_x());
1252     OT_VECTEUR_3D g3b(x,y,z);
1253     OT_VECTEUR_3D n1n4(xyz[0]-noeud1->get_x(),xyz[1]-noeud1->get_y(),xyz[2]-noeud1->get_z());
1254     double alpha1=n1n4*g1b;
1255     double alpha2=n1n4*g2b;
1256     double alpha3=-(n1n4*g3b);
1257     double alpha4=1-alpha1-alpha2;
1258     double eps=0.000001;
1259     if ((alpha1>-eps) && (alpha1<1.+eps))
1260     if ((alpha2>-eps) && (alpha2<1.+eps))
1261     if ((alpha4>-eps) && (alpha4<1.+eps)) {*t=alpha3;return 1;}
1262     return 0;
1263     }
1264    
1265     }