ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur_auto/src/remailleur.cpp
Revision: 980
Committed: Fri Oct 19 00:13:38 2018 UTC (6 years, 6 months ago) by francois
File size: 54883 byte(s)
Log Message:
main de polycristaux dans microstructure avait disparu :)

File Contents

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