ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/remailleur.cpp
Revision: 1189
Committed: Tue Feb 4 17:26:49 2025 UTC (3 months ago) by francois
File size: 55009 byte(s)
Log Message:
Version 5.0 de MAGIC. Integration de ALGLIB pour faire de l'optimisation. ALGLIB se download automatiquement en executant un script dans le repertoire config update_magic.bash


File Contents

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