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

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