ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/mailleur3d_outil.cpp
Revision: 532
Committed: Fri Jul 4 14:39:20 2014 UTC (10 years, 10 months ago) by francois
Original Path: magic/lib/mailleur_auto/src/mailleur3d_outil.cpp
File size: 46142 byte(s)
Log Message:
Resolution de bug avec le mailleur3D et la nouvelle formulation de carte de taille + parametrisation du pas d'integration dans le calcul des longueurs dans une metrique

File Contents

# User Rev Content
1 francois 283 //------------------------------------------------------------
2     //------------------------------------------------------------
3     // MAGiC
4     // Jean Christophe Cuilli�re et Vincent FRANCOIS
5     // D�partement de G�nie M�canique - UQTR
6     //------------------------------------------------------------
7     // Le projet MAGIC est un projet de recherche du d�partement
8     // de g�nie m�canique de l'Universit� du Qu�bec �
9     // Trois Rivi�res
10     // Les librairies ne peuvent �tre utilis�es sans l'accord
11     // des auteurs (contact : francois@uqtr.ca)
12     //------------------------------------------------------------
13     //------------------------------------------------------------
14     //
15     // mailleur3d_outil.cpp
16     //
17     //------------------------------------------------------------
18     //------------------------------------------------------------
19     // COPYRIGHT 2000
20     // Version du 02/03/2006 � 11H23
21     //------------------------------------------------------------
22     //------------------------------------------------------------
23    
24    
25     #include "gestionversion.h"
26    
27    
28     //#pragma hdrstop
29    
30     #include <math.h>
31     #include "mailleur3d.h"
32     #include "mailleur3d_front.h"
33     #include "ot_mathematique.h"
34 francois 481 #include "ot_boite_3d.h"
35 francois 283 #include "m3d_noeud.h"
36     #include "m3d_tetra.h"
37     #include "m3d_triangle.h"
38 francois 335 #include "mg_geometrie_outils.h"
39 francois 388 #include "tpl_fonctions_generiques.h"
40 francois 283 #include "const.h"
41    
42     //---------------------------------------------------------------------------
43    
44     //#pragma package(smart_init)
45    
46    
47    
48    
49    
50    
51    
52     double MAILLEUR3D::angle_front(MG_FRONT_3D* ft1,MG_FRONT_3D* ft2)
53     {
54 francois 388 MG_NOEUD* noeud1=ft1->get_noeud1();
55     MG_NOEUD* noeud2=ft1->get_noeud2();
56     MG_NOEUD* noeud3=ft1->get_noeud3();
57     MG_NOEUD* noeud4=ft2->get_noeud1();
58     MG_NOEUD* noeud5=ft2->get_noeud2();
59     MG_NOEUD* noeud6=ft2->get_noeud3();
60     return get_angle_par_noeud<MG_NOEUD*>(noeud1,noeud2,noeud3,noeud4,noeud5,noeud6);
61 francois 283 }
62    
63    
64    
65     double MAILLEUR3D::angle_front(MG_FRONT_3D* ft,MG_SEGMENT* seg)
66     {
67     MG_NOEUD* noeud1;
68     MG_NOEUD* noeud2;
69     MG_NOEUD* noeud3;
70     MG_NOEUD* noeud4;
71    
72    
73     if (seg->get_noeud1()==ft->get_noeud1())
74     {
75     noeud1=ft->get_noeud1();
76     noeud2=ft->get_noeud2();
77     noeud3=ft->get_noeud3();
78     noeud4=seg->get_noeud2();
79     }
80     if (seg->get_noeud1()==ft->get_noeud2())
81     {
82     noeud1=ft->get_noeud2();
83     noeud2=ft->get_noeud3();
84     noeud3=ft->get_noeud1();
85     noeud4=seg->get_noeud2();
86     }
87     if (seg->get_noeud1()==ft->get_noeud3())
88     {
89     noeud1=ft->get_noeud3();
90     noeud2=ft->get_noeud1();
91     noeud3=ft->get_noeud2();
92     noeud4=seg->get_noeud2();
93     }
94     if (seg->get_noeud2()==ft->get_noeud1())
95     {
96     noeud1=ft->get_noeud1();
97     noeud2=ft->get_noeud2();
98     noeud3=ft->get_noeud3();
99     noeud4=seg->get_noeud1();
100     }
101     if (seg->get_noeud2()==ft->get_noeud2())
102     {
103     noeud1=ft->get_noeud2();
104     noeud2=ft->get_noeud3();
105     noeud3=ft->get_noeud1();
106     noeud4=seg->get_noeud1();
107     }
108     if (seg->get_noeud2()==ft->get_noeud3())
109     {
110     noeud1=ft->get_noeud3();
111     noeud2=ft->get_noeud1();
112     noeud3=ft->get_noeud2();
113     noeud4=seg->get_noeud1();
114     }
115    
116     OT_VECTEUR_3D n1n2(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
117     OT_VECTEUR_3D n1n3(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
118     OT_VECTEUR_3D n1n4(noeud4->get_x()-noeud1->get_x(),noeud4->get_y()-noeud1->get_y(),noeud4->get_z()-noeud1->get_z());
119     OT_VECTEUR_3D n1n=n1n2&n1n3;
120     OT_MATRICE_3D mat(n1n2,n1n3,n1n);
121     OT_MATRICE_3D matx(n1n4,n1n3,n1n);
122     OT_MATRICE_3D maty(n1n2,n1n4,n1n);
123     OT_MATRICE_3D matz(n1n2,n1n3,n1n4);
124     double det=mat.get_determinant();
125     double xsi=matx.get_determinant()/det;
126     double eta=maty.get_determinant()/det;
127     double psi=matz.get_determinant()/det;
128    
129     double proj[3];
130     proj[0]=(1.-xsi-eta)*noeud1->get_x()+xsi*noeud2->get_x()+eta*noeud3->get_x();
131     proj[1]=(1.-xsi-eta)*noeud1->get_y()+xsi*noeud2->get_y()+eta*noeud3->get_y();
132     proj[2]=(1.-xsi-eta)*noeud1->get_z()+xsi*noeud2->get_z()+eta*noeud3->get_z();
133     OT_VECTEUR_3D n1nproj(noeud1->get_x()-proj[0],noeud1->get_y()-proj[1],noeud1->get_z()-proj[2]);
134 francois 419 //n1n4.norme();
135     //n1nproj.norme();
136 francois 283 double ps=n1n4*n1nproj;
137     if (!(OPERATEUR::egal(ps,0.0,0.0001))) ps=ps/n1n4.get_longueur()/n1nproj.get_longueur();
138     else ps=0.;
139     if (ps>1.) ps=1.;
140     if (ps<-1.) ps=-1.;
141     double angle=acos(ps);
142     if ((xsi<0.) || (eta<0.) ) angle=angle+M_PI/2.;
143     if (psi<0.) angle=2*M_PI-angle;
144     return angle;
145    
146     }
147    
148    
149    
150    
151    
152     class MG_TETRA* MAILLEUR3D::insere_tetra(class MG_VOLUME* mgvol,class MG_FRONT_3D* ft,MG_NOEUD* noeud4,int type,TPL_MAP_ENTITE<class MG_TRIANGLE*>& liste_intersection)
153     {
154     liste_intersection.vide();
155     MG_NOEUD* noeud1=ft->get_noeud1();
156     MG_NOEUD* noeud2=ft->get_noeud2();
157     MG_NOEUD* noeud3=ft->get_noeud3();
158    
159     // existence du tetra
160     MG_NOEUD* noeud;
161 francois 363 MG_TETRA* tet=mg_maillage->get_mg_tetra(noeud1,noeud2,noeud3,noeud4);
162 francois 283 if (tet!=NULL)
163     return NULL;
164     // validit� de la solution
165     TPL_MAP_ENTITE<MG_FRONT_3D*> liste;
166     double x=0.25*(noeud1->get_x()+noeud2->get_x()+noeud3->get_x()+noeud4->get_x());
167     double y=0.25*(noeud1->get_y()+noeud2->get_y()+noeud3->get_y()+noeud4->get_y());
168     double z=0.25*(noeud1->get_z()+noeud2->get_z()+noeud3->get_z()+noeud4->get_z());
169     double longueur=1.2*sqrt((x-noeud1->get_x())*(x-noeud1->get_x())+(y-noeud1->get_y())*(y-noeud1->get_y())+(z-noeud1->get_z())*(z-noeud1->get_z()));
170     octree_de_front->rechercher(x,y,z,longueur,liste);
171     int nb_liste=liste.get_nb();
172     for (int i=0;i<nb_liste;i++)
173     {
174     MG_FRONT_3D* fttmp=liste.get(i);
175     int inter=inter_tetra_triangle(noeud1,noeud2,noeud3,noeud4,fttmp->get_noeud1(),fttmp->get_noeud2(),fttmp->get_noeud3());
176 francois 419
177 francois 283 if (inter==true)
178     liste_intersection.ajouter(fttmp->get_triangle());
179     }
180     int nb_inter=liste_intersection.get_nb();
181     if (nb_inter>0)
182     return NULL;
183     // evaluation de la situation locale on ne cree pas deux tetra voisin de volume presque nul
184     MG_TRIANGLE* triangle1=mg_maillage->get_mg_triangle(noeud1->get_id(),noeud3->get_id(),noeud2->get_id());
185     MG_TRIANGLE* triangle2=mg_maillage->get_mg_triangle(noeud1->get_id(),noeud2->get_id(),noeud4->get_id());
186     MG_TRIANGLE* triangle3=mg_maillage->get_mg_triangle(noeud2->get_id(),noeud3->get_id(),noeud4->get_id());
187     MG_TRIANGLE* triangle4=mg_maillage->get_mg_triangle(noeud1->get_id(),noeud4->get_id(),noeud3->get_id());
188     double vol1=-1.;
189     double vol2=-1.;
190     double vol3=-1.;
191     double vol4=-1.;
192     if (type==GENERATION_NOEUD)
193     {
194     if (triangle1!=NULL)
195     if (triangle1->get_type_entite()==IDM3D_TRIANGLE)
196     {
197     M3D_TRIANGLE* tri=(M3D_TRIANGLE*)triangle1;
198     MG_TETRA* tet=tri->get_voisin(0);
199     if (tet->get_type_entite()==IDM3D_TETRA)
200     {
201     M3D_TETRA* mtet=(M3D_TETRA*)tet;
202     vol1=mtet->get_volume();
203     }
204     }
205     if (triangle2!=NULL)
206     if (triangle2->get_type_entite()==IDM3D_TRIANGLE)
207     {
208     M3D_TRIANGLE* tri=(M3D_TRIANGLE*)triangle2;
209     MG_TETRA* tet=tri->get_voisin(0);
210     if (tet->get_type_entite()==IDM3D_TETRA)
211     {
212     M3D_TETRA* mtet=(M3D_TETRA*)tet;
213     vol2=mtet->get_volume();
214     }
215     }
216     if (triangle3!=NULL)
217     if (triangle3->get_type_entite()==IDM3D_TRIANGLE)
218     {
219     M3D_TRIANGLE* tri=(M3D_TRIANGLE*)triangle3;
220     MG_TETRA* tet=tri->get_voisin(0);
221     if (tet->get_type_entite()==IDM3D_TETRA)
222     {
223     M3D_TETRA* mtet=(M3D_TETRA*)tet;
224     vol3=mtet->get_volume();
225     }
226     }
227     if (triangle4!=NULL)
228     if (triangle4->get_type_entite()==IDM3D_TRIANGLE)
229     {
230     M3D_TRIANGLE* tri=(M3D_TRIANGLE*)triangle4;
231     MG_TETRA* tet=tri->get_voisin(0);
232     if (tet->get_type_entite()==IDM3D_TETRA)
233     {
234     M3D_TETRA* mtet=(M3D_TETRA*)tet;
235     vol4=mtet->get_volume();
236     }
237     }
238     double vol=calcule_volume(noeud1,noeud2,noeud3,noeud4);
239     // calcul du vol du tetra ideal
240     OT_VECTEUR_3D ab(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
241     OT_VECTEUR_3D ac(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
242     OT_VECTEUR_3D da(noeud1->get_x()-noeud4->get_x(),noeud1->get_y()-noeud4->get_y(),noeud1->get_z()-noeud4->get_z());
243     OT_VECTEUR_3D vec=ab&ac;
244     vec.norme();
245     double lg_tet_ideal=fabs(vec*da*ft->get_ifail()*sqrt(3./2.));
246     double vol_ideal=sqrt(2.)/12.*lg_tet_ideal*lg_tet_ideal*lg_tet_ideal;
247     if (vol<0.1*vol_ideal)
248     {
249 francois 419 if ((vol1>0.) && (vol1<0.1*vol_ideal))
250 francois 283 return NULL;
251 francois 419 if ((vol2>0.) && (vol2<0.1*vol_ideal))
252     return NULL;
253     if ((vol3>0.) && (vol3<0.1*vol_ideal))
254 francois 283 return NULL;
255 francois 419 if ((vol4>0.) && (vol4<0.1*vol_ideal))
256 francois 283 return NULL;
257 francois 419
258 francois 283 }
259     }
260    
261    
262    
263    
264     // insertion du tetra
265     if (triangle1==NULL) triangle1=insere_triangle(mgvol,noeud1,noeud3,noeud2,MAILLEUR_AUTO);
266     if (triangle2==NULL) triangle2=insere_triangle(mgvol,noeud1,noeud2,noeud4,MAILLEUR_AUTO);
267     if (triangle3==NULL) triangle3=insere_triangle(mgvol,noeud2,noeud3,noeud4,MAILLEUR_AUTO);
268     if (triangle4==NULL) triangle4=insere_triangle(mgvol,noeud1,noeud4,noeud3,MAILLEUR_AUTO);
269     MG_TETRA* tetra=cree_tetra(mgvol,noeud1,noeud2,noeud3,noeud4,triangle1,triangle2,triangle3,triangle4,MAILLEUR_AUTO);
270     mg_maillage->ajouter_mg_tetra(tetra);
271     if (triangle1->get_type_entite()==IDM3D_TRIANGLE)
272     {
273     M3D_TRIANGLE* tri=(M3D_TRIANGLE*)triangle1;
274     tri->ajouter_voisin(tetra);
275     }
276     if (triangle2->get_type_entite()==IDM3D_TRIANGLE)
277     {
278     M3D_TRIANGLE* tri=(M3D_TRIANGLE*)triangle2;
279     tri->ajouter_voisin(tetra);
280     }
281     if (triangle3->get_type_entite()==IDM3D_TRIANGLE)
282     {
283     M3D_TRIANGLE* tri=(M3D_TRIANGLE*)triangle3;
284     tri->ajouter_voisin(tetra);
285     }
286     if (triangle4->get_type_entite()==IDM3D_TRIANGLE)
287     {
288     M3D_TRIANGLE* tri=(M3D_TRIANGLE*)triangle4;
289     tri->ajouter_voisin(tetra);
290     }
291 francois 419
292 francois 283 return tetra;
293    
294     }
295    
296    
297    
298     class MG_TRIANGLE* MAILLEUR3D::insere_triangle(class MG_VOLUME* mgvol,MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,int origine)
299     {
300 francois 419 // existence du triangle
301 francois 283 MG_TRIANGLE* triangle=mg_maillage->get_mg_triangle(noeud1->get_id(),noeud3->get_id(),noeud2->get_id());
302     if (triangle!=NULL) return NULL;
303     // insertion du triangle
304     MG_SEGMENT* segment1=mg_maillage->get_mg_segment(noeud1->get_id(),noeud2->get_id());
305     MG_SEGMENT* segment2=mg_maillage->get_mg_segment(noeud1->get_id(),noeud3->get_id());
306     MG_SEGMENT* segment3=mg_maillage->get_mg_segment(noeud2->get_id(),noeud3->get_id());
307     if (segment1==NULL) segment1=cree_segment(mgvol,noeud1,noeud2,origine);
308     if (segment2==NULL) segment2=cree_segment(mgvol,noeud1,noeud3,origine);
309     if (segment3==NULL) segment3=cree_segment(mgvol,noeud2,noeud3,origine);
310     MG_TRIANGLE* tri=cree_triangle(mgvol,noeud1,noeud2,noeud3,segment1,segment2,segment3,origine);
311     return tri;
312     }
313    
314    
315    
316     class MG_SEGMENT* MAILLEUR3D::insere_segment(class MG_VOLUME* mgvol,MG_NOEUD* noeud1,MG_NOEUD* noeud2,int origine)
317     {
318 francois 419
319 francois 283 //existence du segment
320     MG_SEGMENT* segment=mg_maillage->get_mg_segment(noeud1->get_id(),noeud2->get_id());
321     if (segment!=NULL) return NULL;
322     // insertion du segment
323     MG_SEGMENT* seg=cree_segment(mgvol,noeud1,noeud2,origine);
324     return seg;
325     }
326    
327    
328    
329    
330     MG_TETRA* MAILLEUR3D::cree_tetra(class MG_VOLUME* mgvol,MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,MG_NOEUD* noeud4,MG_TRIANGLE* triangle1,MG_TRIANGLE* triangle2,MG_TRIANGLE* triangle3,MG_TRIANGLE* triangle4,int origine)
331     {
332     M3D_TETRA* tetra=new M3D_TETRA(mgvol,noeud1,noeud2,noeud3,noeud4,triangle1,triangle2,triangle3,triangle4,origine);
333     mg_maillage->ajouter_mg_tetra(tetra);
334     return tetra;
335     }
336    
337    
338     MG_TRIANGLE* MAILLEUR3D::cree_triangle(class MG_VOLUME* mgvol,MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,MG_SEGMENT* segment1,MG_SEGMENT* segment2,MG_SEGMENT* segment3,int origine)
339     {
340     M3D_TRIANGLE* tri=new M3D_TRIANGLE(mgvol,noeud1,noeud2,noeud3,segment1,segment2,segment3,origine);
341     mg_maillage->ajouter_mg_triangle(tri);
342     return tri;
343     }
344    
345    
346    
347     MG_SEGMENT* MAILLEUR3D::cree_segment(class MG_VOLUME* mgvol,MG_NOEUD* noeud1,MG_NOEUD* noeud2,int origine)
348     {
349     MG_SEGMENT* seg=mg_maillage->ajouter_mg_segment(mgvol,noeud1,noeud2,origine);
350     return seg;
351     }
352    
353    
354    
355     double MAILLEUR3D::calcule_volume(MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,MG_NOEUD* noeud4)
356     {
357     OT_VECTEUR_3D ab(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
358     OT_VECTEUR_3D ac(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
359     OT_VECTEUR_3D da(noeud1->get_x()-noeud4->get_x(),noeud1->get_y()-noeud4->get_y(),noeud1->get_z()-noeud4->get_z());
360    
361     OT_VECTEUR_3D pvec=ab&ac;
362     double volume=pvec*da;
363     volume=fabs(volume)/3.;
364     return volume;
365     }
366    
367    
368    
369     double MAILLEUR3D::calcule_longueur_caracteristique(MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,MG_NOEUD* noeud4)
370     {
371    
372     OT_VECTEUR_3D ab(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
373     double lg=ab.get_longueur();
374     int nb=1;
375     if (noeud3!=NULL)
376     {
377     OT_VECTEUR_3D ac(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
378     OT_VECTEUR_3D bc(noeud3->get_x()-noeud2->get_x(),noeud3->get_y()-noeud2->get_y(),noeud3->get_z()-noeud2->get_z());
379     lg=lg+ac.get_longueur()+bc.get_longueur();
380     nb=3;
381     }
382     if (noeud4!=NULL)
383     {
384     OT_VECTEUR_3D da(noeud1->get_x()-noeud4->get_x(),noeud1->get_y()-noeud4->get_y(),noeud1->get_z()-noeud4->get_z());
385     OT_VECTEUR_3D bd(noeud4->get_x()-noeud2->get_x(),noeud4->get_y()-noeud2->get_y(),noeud4->get_z()-noeud2->get_z());
386     OT_VECTEUR_3D cd(noeud4->get_x()-noeud3->get_x(),noeud4->get_y()-noeud3->get_y(),noeud4->get_z()-noeud3->get_z());
387     lg=lg+da.get_longueur()+bd.get_longueur()+cd.get_longueur();
388     nb=6;
389     }
390     lg=lg*1./nb;
391     return lg;
392    
393     }
394    
395 francois 532 double MAILLEUR3D::calcul_distance_metrique_segment(MG_SEGMENT *seg,int pas)
396     {
397     double l;
398     bool res=seg->get_distance_metrique(&l);
399     if (res==false)
400     {
401     l=metrique->calcul_distance_metrique(seg,pas);
402     seg->change_distance_metrique(l);
403     }
404     return l;
405     }
406 francois 283
407    
408     int MAILLEUR3D::genere_noeud(MG_VOLUME* mgvol,MG_FRONT_3D* ft,TPL_MAP_ENTITE<MG_NOEUD*> &liste_noeud,std::vector<CAS_FRONT*> &liste_cas)
409     {
410     MG_NOEUD* noeud1=ft->get_noeud1();
411     MG_NOEUD* noeud2=ft->get_noeud2();
412     MG_NOEUD* noeud3=ft->get_noeud3();
413     OT_VECTEUR_3D n1n2(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
414     OT_VECTEUR_3D n1n3(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
415     OT_VECTEUR_3D normal=n1n2&n1n3;
416     normal.norme();
417     double xyz1[3]={noeud1->get_x(),noeud1->get_y(),noeud1->get_z()};
418     double xyz2[3]={noeud2->get_x(),noeud2->get_y(),noeud2->get_z()};
419     double xyz3[3]={noeud3->get_x(),noeud3->get_y(),noeud3->get_z()};
420     double hauteur=0.81649658*ft->get_triangle()->get_longueur();
421     double xyzg[3]={0.333333333333*(xyz1[0]+xyz2[0]+xyz3[0]),0.333333333333*(xyz1[1]+xyz2[1]+xyz3[1]),0.333333333333*(xyz1[2]+xyz2[2]+xyz3[2])};
422     double xyz4[3];
423     do
424     {
425     xyz4[0]=xyzg[0]+hauteur*normal.get_x();
426     xyz4[1]=xyzg[1]+hauteur*normal.get_y();
427     xyz4[2]=xyzg[2]+hauteur*normal.get_z();
428     hauteur=hauteur*0.9;
429     }
430     while (metrique->valide_parametre(xyz4)==0);
431 francois 532 double l1=calcul_distance_metrique_segment(ft->get_triangle()->get_segment1(),pas);
432     double l2=calcul_distance_metrique_segment(ft->get_triangle()->get_segment2(),pas);
433     double l3=calcul_distance_metrique_segment(ft->get_triangle()->get_segment3(),pas);
434 francois 283 double l=0.333333333333333*(l1+l2+l3);
435     double longueur_desiree=0.81649658*(priorite_metrique+(1-priorite_metrique)*l);
436     longueur_desiree=longueur_desiree/ft->get_ifail();
437 francois 532 metrique->ajuste_distance_metrique(xyzg,xyz4,longueur_desiree,pas);
438 francois 283 OT_VECTEUR_3D hauteurajustee(xyzg[0]-xyz4[0],xyzg[1]-xyz4[1],xyzg[2]-xyz4[2]);
439     hauteur=hauteurajustee.get_longueur();
440    
441     // recherche des voisins
442     TPL_MAP_ENTITE<MG_FRONT_3D*> liste_front;
443     TPL_MAP_ENTITE<MG_FRONT_3D*> liste_front_proche;
444     octree_de_front->rechercher(xyz4[0],xyz4[1],xyz4[2],hauteur*1.837117307,liste_front);
445     int nb_liste=liste_front.get_nb();
446     double disref=(priorite_metrique+(1-priorite_metrique)*l)/ft->get_ifail();
447     MG_NOEUD* noeudref=NULL;
448     for (int i=0;i<nb_liste;i++)
449     {
450     MG_FRONT_3D* ft2=liste_front.get(i);
451     for (int j=0;j<3;j++)
452     {
453     MG_NOEUD* noeud_tmp;
454     if (j==0) noeud_tmp=ft2->get_noeud1();
455     if (j==1) noeud_tmp=ft2->get_noeud2();
456     if (j==2) noeud_tmp=ft2->get_noeud3();
457     double *xyz_tmp=noeud_tmp->get_coord();
458 francois 532 double dis=metrique->calcul_distance_metrique(xyz_tmp,xyz4,1);
459 francois 283 if (dis<disref)
460     {
461     disref=dis;
462     noeudref=noeud_tmp;
463     }
464     }
465     if (ft2!=ft)
466     {
467     double dismin=evaluer_distance_noeud_triangle(xyz4[0],xyz4[1],xyz4[2],ft2->get_triangle());
468     if (dismin<hauteur) liste_front_proche.ajouter(ft2);
469     }
470     }
471     // il y a un noeud proche
472     if (noeudref!=NULL)
473     {
474 francois 532 double dis2=metrique->calcul_distance_metrique(noeudref->get_coord(),xyzg,2);
475 francois 283 if (dis2<longueur_desiree*1.5)
476     if (noeudref!=noeud1)
477     if (noeudref!=noeud2)
478     if (noeudref!=noeud3)
479     {
480     liste_noeud.ajouter(noeudref);
481     return EXISTE;
482     }
483     }
484     // il n y a pas de noeud proche mais il y a une face proche
485     int nb_proche=liste_front_proche.get_nb();
486     int nb_cas=liste_cas.size();
487     if (nb_proche>0)
488     {
489     for (int i=0;i<nb_proche;i++)
490     {
491     MG_FRONT_3D* fttmp=liste_front_proche.get(i);
492     int ok=0;
493     for (int j=0;j<nb_cas;j++)
494     {
495     MG_NOEUD* casnoeud=liste_cas[j]->mgnoeud;
496     if ((fttmp->get_noeud1()==casnoeud) || (fttmp->get_noeud2()==casnoeud) || (fttmp->get_noeud3()==casnoeud) )
497     {
498     ok=1;
499     break;
500     }
501     }
502     if (ok==1) continue;
503     MG_NOEUD* noeudtmp=fttmp->get_noeud1();
504     if (noeudtmp!=noeud1)
505     if (noeudtmp!=noeud2)
506     if (noeudtmp!=noeud3) liste_noeud.ajouter(noeudtmp);
507     noeudtmp=fttmp->get_noeud2();
508     if (noeudtmp!=noeud1)
509     if (noeudtmp!=noeud2)
510     if (noeudtmp!=noeud3) liste_noeud.ajouter(noeudtmp);
511     noeudtmp=fttmp->get_noeud3();
512     if (noeudtmp!=noeud1)
513     if (noeudtmp!=noeud2)
514     if (noeudtmp!=noeud3) liste_noeud.ajouter(noeudtmp);
515    
516     }
517     if (liste_noeud.get_nb()>0) return MULTIEXISTE;
518     }
519     // il faut creer un noeud c est parfait c est qu on veut.
520     M3D_NOEUD* nouveau_noeud=new M3D_NOEUD(mgvol,xyz4[0],xyz4[1],xyz4[2],MAILLEUR_AUTO);
521     mg_maillage->ajouter_mg_noeud(nouveau_noeud);
522     liste_noeud.ajouter(nouveau_noeud);
523     return CREATION;
524     }
525    
526    
527    
528     double MAILLEUR3D::evaluer_distance_noeud_triangle(double x,double y,double z,MG_TRIANGLE* tri)
529     {
530     double eps=0.0001;
531     MG_NOEUD* noeud1=tri->get_noeud1();
532     MG_NOEUD* noeud2=tri->get_noeud2();
533     MG_NOEUD* noeud3=tri->get_noeud3();
534     OT_VECTEUR_3D n1n2(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
535     OT_VECTEUR_3D n1n3(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
536     OT_VECTEUR_3D n1n4(x-noeud1->get_x(),y-noeud1->get_y(),z-noeud1->get_z());
537     OT_VECTEUR_3D n1n=n1n2&n1n3;
538     OT_MATRICE_3D matbase(n1n2,n1n3,n1n);
539     OT_MATRICE_3D matxsi(n1n4,n1n3,n1n);
540     OT_MATRICE_3D mateta(n1n2,n1n4,n1n);
541     double xsi=matxsi.get_determinant()/matbase.get_determinant();
542     double eta=mateta.get_determinant()/matbase.get_determinant();
543     double dis;
544     if ((eta>-eps) && (xsi>-eps) && (eta+xsi<1.0+eps))
545     {
546     double xproj=(1-xsi-eta)*noeud1->get_x()+xsi*noeud2->get_x()+eta*noeud3->get_x();
547     double yproj=(1-xsi-eta)*noeud1->get_y()+xsi*noeud2->get_y()+eta*noeud3->get_y();
548     double zproj=(1-xsi-eta)*noeud1->get_z()+xsi*noeud2->get_z()+eta*noeud3->get_z();
549     OT_VECTEUR_3D vec(x-xproj,y-yproj,z-zproj);
550     dis=vec.get_longueur();
551     }
552     else
553     {
554     OT_VECTEUR_3D nn1(x-noeud1->get_x(),y-noeud1->get_y(),z-noeud1->get_z());
555     OT_VECTEUR_3D nn2(x-noeud2->get_x(),y-noeud2->get_y(),z-noeud2->get_z());
556     OT_VECTEUR_3D nn3(x-noeud3->get_x(),y-noeud3->get_y(),z-noeud3->get_z());
557     dis=0.3333333333333*(nn1.get_longueur()+nn2.get_longueur()+nn3.get_longueur());
558     }
559     return dis;
560     }
561    
562    
563     void MAILLEUR3D::detruit_tetra(MG_TRIANGLE* tri)
564     {
565 francois 419
566 francois 283 M3D_TRIANGLE* mtri;
567     //if (tri->get_type_entite()==IDMG_TRIANGLE) return;
568     //else mtri=(M3D_TRIANGLE*)tri;
569     mtri=(M3D_TRIANGLE*)tri;
570    
571     M3D_TETRA* tet=(M3D_TETRA*)mtri->get_voisin(0);
572     M3D_TRIANGLE* mtriangle1=(M3D_TRIANGLE*)tet->get_triangle1();
573     M3D_TRIANGLE* mtriangle2=(M3D_TRIANGLE*)tet->get_triangle2();
574     M3D_TRIANGLE* mtriangle3=(M3D_TRIANGLE*)tet->get_triangle3();
575     M3D_TRIANGLE* mtriangle4=(M3D_TRIANGLE*)tet->get_triangle4();
576 francois 532 if (tet->get_etat()==DETRUIT)
577     std::cout << " TTTTTTTTT" << std::endl;
578 francois 283 tet->etat_detruit();
579     TPL_MAP_ENTITE<MG_FRONT_3D*> liste_de_front;
580     TPL_MAP_ENTITE<MG_FRONT_3D*> liste_de_front2;
581     MG_FRONT_3D* fttet1=NULL;
582     MG_FRONT_3D* fttet2=NULL;
583     MG_FRONT_3D* fttet3=NULL;
584     MG_FRONT_3D* fttet4=NULL;
585    
586     if (mtriangle1->get_etat_front()==FRONT_ACTIF)
587     {
588     MG_FRONT_3D* fttmp=mtriangle1->get_mgfront();
589     fttet1=fttmp;
590     liste_de_front.ajouter(fttmp->get_front_voisin(0));
591     liste_de_front.ajouter(fttmp->get_front_voisin(1));
592     liste_de_front.ajouter(fttmp->get_front_voisin(2));
593     //supprimer_front_en_reculant(fttmp);
594     mtriangle1->supprimer_voisin(tet);
595     }
596     else
597     {
598     MG_FRONT_3D* ftn1=ajouter_front_courant(PRIORITAIRE,tet->get_noeud1(),tet->get_noeud2(),tet->get_noeud3(),mtriangle1);
599     ftn1->incremente_ifail();
600     mtriangle1->supprimer_voisin(tet);
601     liste_de_front2.ajouter(ftn1);
602     }
603     if (mtriangle2->get_etat_front()==FRONT_ACTIF)
604     {
605     MG_FRONT_3D* fttmp=mtriangle2->get_mgfront();
606     fttet2=fttmp;
607     liste_de_front.ajouter(fttmp->get_front_voisin(0));
608     liste_de_front.ajouter(fttmp->get_front_voisin(1));
609     liste_de_front.ajouter(fttmp->get_front_voisin(2));
610     //supprimer_front_en_reculant(fttmp);
611     mtriangle2->supprimer_voisin(tet);
612     }
613     else
614     {
615     MG_FRONT_3D* ftn1=ajouter_front_courant(PRIORITAIRE,tet->get_noeud1(),tet->get_noeud4(),tet->get_noeud2(),mtriangle2);
616     ftn1->incremente_ifail();
617     mtriangle2->supprimer_voisin(tet);
618     liste_de_front2.ajouter(ftn1);
619     }
620     if (mtriangle3->get_etat_front()==FRONT_ACTIF)
621     {
622     MG_FRONT_3D* fttmp=mtriangle3->get_mgfront();
623     fttet3=fttmp;
624     liste_de_front.ajouter(fttmp->get_front_voisin(0));
625     liste_de_front.ajouter(fttmp->get_front_voisin(1));
626     liste_de_front.ajouter(fttmp->get_front_voisin(2));
627     //supprimer_front_en_reculant(fttmp);
628     mtriangle3->supprimer_voisin(tet);
629     }
630     else
631     {
632     MG_FRONT_3D* ftn1=ajouter_front_courant(PRIORITAIRE,tet->get_noeud2(),tet->get_noeud4(),tet->get_noeud3(),mtriangle3);
633     ftn1->incremente_ifail();
634     mtriangle3->supprimer_voisin(tet);
635     liste_de_front2.ajouter(ftn1);
636     }
637     if (mtriangle4->get_etat_front()==FRONT_ACTIF)
638     {
639     MG_FRONT_3D* fttmp=mtriangle4->get_mgfront();
640     fttet4=fttmp;
641     liste_de_front.ajouter(fttmp->get_front_voisin(0));
642     liste_de_front.ajouter(fttmp->get_front_voisin(1));
643     liste_de_front.ajouter(fttmp->get_front_voisin(2));
644     //supprimer_front_en_reculant(fttmp);
645     mtriangle4->supprimer_voisin(tet);
646     }
647     else
648     {
649     MG_FRONT_3D* ftn1=ajouter_front_courant(PRIORITAIRE,tet->get_noeud1(),tet->get_noeud3(),tet->get_noeud4(),mtriangle4);
650     ftn1->incremente_ifail();
651     mtriangle4->supprimer_voisin(tet);
652     liste_de_front2.ajouter(ftn1);
653     }
654    
655     if (fttet1!=NULL)
656     {
657     liste_de_front.supprimer(fttet1);
658     supprimer_front_en_reculant(fttet1);
659     }
660     if (fttet2!=NULL)
661     {
662     liste_de_front.supprimer(fttet2);
663     supprimer_front_en_reculant(fttet2);
664     }
665     if (fttet3!=NULL)
666     {
667     liste_de_front.supprimer(fttet3);
668     supprimer_front_en_reculant(fttet3);
669     }
670     if (fttet4!=NULL)
671     {
672     liste_de_front.supprimer(fttet4);
673     supprimer_front_en_reculant(fttet4);
674     }
675     int nb_liste1=liste_de_front.get_nb();
676     for (int i=0;i<nb_liste1;i++)
677     if (liste_de_front.get(i)!=NULL)
678     mise_a_jour_voisin_front(liste_de_front.get(i));
679     int nb_liste2=liste_de_front2.get_nb();
680     for (int i=0;i<nb_liste2;i++)
681     {
682     mise_a_jour_voisin_front(liste_de_front2.get(i));
683     mise_a_jour_voisin_front(liste_de_front2.get(i)->get_front_voisin(0));
684     mise_a_jour_voisin_front(liste_de_front2.get(i)->get_front_voisin(1));
685     mise_a_jour_voisin_front(liste_de_front2.get(i)->get_front_voisin(2));
686     }
687    
688     for (int k=0;k<4;k++)
689     {
690     MG_NOEUD* noeud;
691     if (k==0) noeud=tet->get_noeud1();
692     if (k==1) noeud=tet->get_noeud2();
693     if (k==2) noeud=tet->get_noeud3();
694     if (k==3) noeud=tet->get_noeud4();
695     int nb_tet=noeud->get_lien_tetra()->get_nb();
696    
697     if (noeud->get_type_entite()==IDM3D_NOEUD)
698     {
699     int test=false;
700     M3D_NOEUD* mnoeud=(M3D_NOEUD*)noeud;
701     for (int i=0;i<nb_tet;i++)
702     {
703     MG_TETRA* tet=noeud->get_lien_tetra()->get(i);
704     if (tet->get_type_entite()==IDM3D_TETRA)
705     {
706     M3D_TETRA* mtet=(M3D_TETRA*)tet;
707     if (mtet->get_etat()==ACTIF) test=true;
708     }
709     else test=true;
710     }
711     if (test==false) mnoeud->etat_detruit();
712     }
713     }
714 francois 419
715 francois 283 }
716    
717    
718    
719    
720 francois 532
721 francois 283 void MAILLEUR3D::detruit_element_inactif(void)
722     {
723     TPL_MAP_ENTITE<MG_TETRA*> lsttet;
724     TPL_MAP_ENTITE<MG_NOEUD*> lstnoeud;
725     LISTE_MG_TETRA::iterator it;
726    
727     for (MG_TETRA* tet=mg_maillage->get_premier_tetra(it);tet!=NULL;tet=mg_maillage->get_suivant_tetra(it))
728     {
729     M3D_TETRA* mtet=(M3D_TETRA*)tet;
730     if (mtet->get_etat()==DETRUIT)
731     lsttet.ajouter(tet);
732     }
733     LISTE_MG_NOEUD::iterator it2;
734     for (MG_NOEUD* noeud=mg_maillage->get_premier_noeud(it2);noeud!=NULL;noeud=mg_maillage->get_suivant_noeud(it2))
735     {
736     if (noeud->get_type_entite()==IDM3D_NOEUD)
737     {
738     M3D_NOEUD* mnoeud=(M3D_NOEUD*)noeud;
739     if (mnoeud->get_etat()==DETRUIT)
740     lstnoeud.ajouter(noeud);
741     }
742     }
743     TPL_MAP_ENTITE<MG_TETRA*>::ITERATEUR it3;
744     for (MG_TETRA* tet=lsttet.get_premier(it3);tet!=NULL;tet=lsttet.get_suivant(it3))
745     mg_maillage->supprimer_mg_tetraid(tet->get_id());
746     TPL_MAP_ENTITE<MG_NOEUD*>::ITERATEUR it4;
747     for (MG_NOEUD* noeud=lstnoeud.get_premier(it4);noeud!=NULL;noeud=lstnoeud.get_suivant(it4))
748     mg_maillage->supprimer_mg_noeudid(noeud->get_id());
749    
750 francois 419
751 francois 283 }
752    
753    
754    
755    
756    
757 francois 424 void MAILLEUR3D::analyse_maillage_obtenu(double &vol,int *tab,double borne1,double borne2,double borne3)
758 francois 283 {
759 francois 424 vol=0.;tab[0]=0;tab[1]=0;tab[2]=0;tab[3]=0;
760     LISTE_MG_TETRA::iterator it;
761     for (MG_TETRA* tet=mg_maillage->get_premier_tetra(it);tet!=NULL;tet=mg_maillage->get_suivant_tetra(it))
762 francois 283 {
763 francois 532 vol=vol+metrique->calcul_volume_tetra_metrique(tet);
764 francois 424 double qual=OPERATEUR::qualite_tetra(tet->get_noeud1()->get_coord(),tet->get_noeud2()->get_coord(),tet->get_noeud3()->get_coord(),tet->get_noeud4()->get_coord());
765     if (qual>=borne3) tab[3]++;
766     else if (qual>=borne2) tab[2]++;
767     else if (qual>=borne1) tab[1]++;
768     else tab[0]++;
769 francois 283 }
770     }
771    
772    
773     void MAILLEUR3D::insere_contrainte_tetra(MG_VOLUME* mgvol,TPL_LISTE_ENTITE<MG_TETRA*> *lsttet)
774     {
775     map<long,long,less<long> > lstcorrespondancenoeud;
776     int nbtet=lsttet->get_nb();
777     for (int i=0;i<nbtet;i++)
778     {
779     M3D_TETRA* tet=(M3D_TETRA*)lsttet->get(i);
780     double *xyz1=tet->get_noeud1()->get_coord();
781     double *xyz2=tet->get_noeud2()->get_coord();
782     double *xyz3=tet->get_noeud3()->get_coord();
783     double *xyz4=tet->get_noeud4()->get_coord();
784     double xyz[3];
785     xyz[0]=(xyz1[0]+xyz2[0]+xyz3[0]+xyz4[0])*0.25;
786     xyz[1]=(xyz1[1]+xyz2[1]+xyz3[1]+xyz4[1])*0.25;
787     xyz[2]=(xyz1[2]+xyz2[2]+xyz3[2]+xyz4[2])*0.25;
788     TPL_MAP_ENTITE<MG_FRONT_3D*> liste_trouvee;
789     OT_VECTEUR_3D vec1(xyz,xyz1);
790     OT_VECTEUR_3D vec2(xyz,xyz2);
791     OT_VECTEUR_3D vec3(xyz,xyz3);
792     OT_VECTEUR_3D vec4(xyz,xyz4);
793     double rayonrecherche=vec1.get_longueur();
794     if (vec2.get_longueur()>rayonrecherche) rayonrecherche=vec2.get_longueur();
795     if (vec3.get_longueur()>rayonrecherche) rayonrecherche=vec3.get_longueur();
796     if (vec4.get_longueur()>rayonrecherche) rayonrecherche=vec4.get_longueur();
797     rayonrecherche=rayonrecherche*1.1;
798     octree_de_front->rechercher(xyz[0],xyz[1],xyz[2],rayonrecherche,liste_trouvee);
799     MG_NOEUD* noeudref1=NULL;
800     if (lstcorrespondancenoeud.find(tet->get_noeud1()->get_nouveau_numero())!=lstcorrespondancenoeud.end())
801     noeudref1=mg_maillage->get_mg_noeudid(lstcorrespondancenoeud[tet->get_noeud1()->get_nouveau_numero()]);
802     MG_NOEUD* noeudref2=NULL;
803     if (lstcorrespondancenoeud.find(tet->get_noeud2()->get_nouveau_numero())!=lstcorrespondancenoeud.end())
804     noeudref2=mg_maillage->get_mg_noeudid(lstcorrespondancenoeud[tet->get_noeud2()->get_nouveau_numero()]);
805     MG_NOEUD* noeudref3=NULL;
806     if (lstcorrespondancenoeud.find(tet->get_noeud3()->get_nouveau_numero())!=lstcorrespondancenoeud.end())
807     noeudref3=mg_maillage->get_mg_noeudid(lstcorrespondancenoeud[tet->get_noeud3()->get_nouveau_numero()]);
808     MG_NOEUD* noeudref4=NULL;
809     if (lstcorrespondancenoeud.find(tet->get_noeud4()->get_nouveau_numero())!=lstcorrespondancenoeud.end())
810     noeudref4=mg_maillage->get_mg_noeudid(lstcorrespondancenoeud[tet->get_noeud4()->get_nouveau_numero()]);
811     int nbtrouve=liste_trouvee.get_nb();
812     for (int j=0;j<nbtrouve;j++)
813     {
814     if ((noeudref1!=NULL) && (noeudref2!=NULL) && (noeudref3!=NULL) && (noeudref4!=NULL) ) break;
815     MG_FRONT_3D* ft=liste_trouvee.get(j);
816     TPL_MAP_ENTITE<MG_NOEUD*> lstnoeud;
817     lstnoeud.ajouter(ft->get_noeud1());
818     lstnoeud.ajouter(ft->get_noeud2());
819     lstnoeud.ajouter(ft->get_noeud3());
820     if (noeudref1==NULL)
821     {
822     int nb=lstnoeud.get_nb();
823     for (int k=0;k<nb;k++)
824     {
825     double *xyztmp=lstnoeud.get(k)->get_coord();
826     OT_VECTEUR_3D vec(xyztmp,xyz1);
827     double longueur=vec.get_longueur();
828     if (longueur<1e-6*rayonrecherche)
829     {
830     noeudref1=lstnoeud.get(k);
831     lstnoeud.supprimer(lstnoeud.get(k));
832     break;
833     }
834     }
835     }
836     if (noeudref2==NULL)
837     {
838     int nb=lstnoeud.get_nb();
839     for (int k=0;k<nb;k++)
840     {
841     double *xyztmp=lstnoeud.get(k)->get_coord();
842     OT_VECTEUR_3D vec(xyztmp,xyz2);
843     double longueur=vec.get_longueur();
844     if (longueur<1e-6*rayonrecherche)
845     {
846     noeudref2=lstnoeud.get(k);
847     lstnoeud.supprimer(lstnoeud.get(k));
848     break;
849     }
850     }
851     }
852     if (noeudref3==NULL)
853     {
854     int nb=lstnoeud.get_nb();
855     for (int k=0;k<nb;k++)
856     {
857     double *xyztmp=lstnoeud.get(k)->get_coord();
858     OT_VECTEUR_3D vec(xyztmp,xyz3);
859     double longueur=vec.get_longueur();
860     if (longueur<1e-6*rayonrecherche)
861     {
862     noeudref3=lstnoeud.get(k);
863     lstnoeud.supprimer(lstnoeud.get(k));
864     break;
865     }
866     }
867     }
868     if (noeudref4==NULL)
869     {
870     int nb=lstnoeud.get_nb();
871     for (int k=0;k<nb;k++)
872     {
873     double *xyztmp=lstnoeud.get(k)->get_coord();
874     OT_VECTEUR_3D vec(xyztmp,xyz4);
875     double longueur=vec.get_longueur();
876     if (longueur<1e-6*rayonrecherche)
877     {
878     noeudref4=lstnoeud.get(k);
879     lstnoeud.supprimer(lstnoeud.get(k));
880     break;
881     }
882     }
883     }
884     }
885     if (noeudref1==NULL)
886     {
887     noeudref1=mg_maillage->ajouter_mg_noeud(mgvol,xyz1[0],xyz1[1],xyz1[2],IMPOSE);
888     if (tet->get_noeud1()->get_nouveau_numero()!=(-1))
889     {
890     noeudref1->change_nouveau_numero(tet->get_noeud1()->get_nouveau_numero());
891     pair<long,long> tmp(noeudref1->get_nouveau_numero(),noeudref1->get_id());
892     lstcorrespondancenoeud.insert(tmp);
893     }
894     }
895     if (noeudref2==NULL)
896     {
897     noeudref2=mg_maillage->ajouter_mg_noeud(mgvol,xyz2[0],xyz2[1],xyz2[2],IMPOSE);
898     if (tet->get_noeud2()->get_nouveau_numero()!=(-1))
899     {
900     noeudref2->change_nouveau_numero(tet->get_noeud2()->get_nouveau_numero());
901     pair<long,long> tmp(noeudref2->get_nouveau_numero(),noeudref2->get_id());
902     lstcorrespondancenoeud.insert(tmp);
903     }
904     }
905     if (noeudref3==NULL)
906     {
907     noeudref3=mg_maillage->ajouter_mg_noeud(mgvol,xyz3[0],xyz3[1],xyz3[2],IMPOSE);
908     if (tet->get_noeud3()->get_nouveau_numero()!=(-1))
909     {
910     noeudref3->change_nouveau_numero(tet->get_noeud3()->get_nouveau_numero());
911     pair<long,long> tmp(noeudref3->get_nouveau_numero(),noeudref3->get_id());
912     lstcorrespondancenoeud.insert(tmp);
913     }
914     }
915     if (noeudref4==NULL)
916     {
917     noeudref4=mg_maillage->ajouter_mg_noeud(mgvol,xyz4[0],xyz4[1],xyz4[2],IMPOSE);
918     if (tet->get_noeud4()->get_nouveau_numero()!=(-1))
919     {
920     noeudref4->change_nouveau_numero(tet->get_noeud4()->get_nouveau_numero());
921     pair<long,long> tmp(noeudref4->get_nouveau_numero(),noeudref4->get_id());
922     lstcorrespondancenoeud.insert(tmp);
923     }
924     }
925     M3D_TRIANGLE* triangle1=(M3D_TRIANGLE*)mg_maillage->get_mg_triangle(noeudref1->get_id(),noeudref3->get_id(),noeudref2->get_id());
926     M3D_TRIANGLE* triangle2=(M3D_TRIANGLE*)mg_maillage->get_mg_triangle(noeudref1->get_id(),noeudref2->get_id(),noeudref4->get_id());
927     M3D_TRIANGLE* triangle3=(M3D_TRIANGLE*)mg_maillage->get_mg_triangle(noeudref2->get_id(),noeudref3->get_id(),noeudref4->get_id());
928     M3D_TRIANGLE* triangle4=(M3D_TRIANGLE*)mg_maillage->get_mg_triangle(noeudref1->get_id(),noeudref4->get_id(),noeudref3->get_id());
929     if (triangle1==NULL) triangle1=(M3D_TRIANGLE*)insere_triangle(mgvol,noeudref1,noeudref3,noeudref2,IMPOSE);
930     if (triangle2==NULL) triangle2=(M3D_TRIANGLE*)insere_triangle(mgvol,noeudref1,noeudref2,noeudref4,IMPOSE);
931     if (triangle3==NULL) triangle3=(M3D_TRIANGLE*)insere_triangle(mgvol,noeudref2,noeudref3,noeudref4,IMPOSE);
932     if (triangle4==NULL) triangle4=(M3D_TRIANGLE*)insere_triangle(mgvol,noeudref1,noeudref4,noeudref3,IMPOSE);
933     MG_TETRA* tetra=cree_tetra(mgvol,noeudref1,noeudref2,noeudref3,noeudref4,triangle1,triangle2,triangle3,triangle4,IMPOSE);
934     tetra->change_nouveau_numero(tet->get_nouveau_numero());
935     mg_maillage->ajouter_mg_tetra(tetra);
936     triangle1->ajouter_voisin(tetra);
937     triangle2->ajouter_voisin(tetra);
938     triangle3->ajouter_voisin(tetra);
939     triangle4->ajouter_voisin(tetra);
940     MG_FRONT_3D *ftref1=NULL;
941     MG_FRONT_3D *ftref2=NULL;
942     MG_FRONT_3D *ftref3=NULL;
943     MG_FRONT_3D *ftref4=NULL;
944     if (triangle1->get_etat_front()==FRONT_ACTIF) ftref1=triangle1->get_mgfront();
945     if (triangle2->get_etat_front()==FRONT_ACTIF) ftref2=triangle2->get_mgfront();
946     if (triangle3->get_etat_front()==FRONT_ACTIF) ftref3=triangle3->get_mgfront();
947     if (triangle4->get_etat_front()==FRONT_ACTIF) ftref4=triangle4->get_mgfront();
948     TPL_MAP_ENTITE<MG_FRONT_3D*> lstfront;
949     MG_FRONT_3D *ft1=NULL;
950     MG_FRONT_3D *ft2=NULL;
951     MG_FRONT_3D *ft3=NULL;
952     MG_FRONT_3D *ft4=NULL;
953     if (ftref1==NULL) ft1=mise_a_jour_front(tetra,0,NONFORCE);
954     else
955     {
956     lstfront.ajouter(ftref1->get_front_voisin(0));
957     lstfront.ajouter(ftref1->get_front_voisin(1));
958     lstfront.ajouter(ftref1->get_front_voisin(2));
959     supprimer_front_en_avancant_sans_delete(ftref1);
960     }
961     if (ftref2==NULL) ft2=mise_a_jour_front(tetra,1,NONFORCE);
962     else
963     {
964     lstfront.ajouter(ftref2->get_front_voisin(0));
965     lstfront.ajouter(ftref2->get_front_voisin(1));
966     lstfront.ajouter(ftref2->get_front_voisin(2));
967     supprimer_front_en_avancant_sans_delete(ftref2);
968     }
969     if (ftref3==NULL) ft3=mise_a_jour_front(tetra,2,NONFORCE);
970     else
971     {
972     lstfront.ajouter(ftref3->get_front_voisin(0));
973     lstfront.ajouter(ftref3->get_front_voisin(1));
974     lstfront.ajouter(ftref3->get_front_voisin(2));
975     supprimer_front_en_avancant_sans_delete(ftref3);
976     }
977     if (ftref4==NULL) ft4=mise_a_jour_front(tetra,3,NONFORCE);
978     else
979     {
980     lstfront.ajouter(ftref4->get_front_voisin(0));
981     lstfront.ajouter(ftref4->get_front_voisin(1));
982     lstfront.ajouter(ftref4->get_front_voisin(2));
983     supprimer_front_en_avancant_sans_delete(ftref4);
984     }
985     if (ftref1!=NULL) {
986     lstfront.supprimer(ftref1);
987     delete ftref1;
988     }
989     if (ftref2!=NULL) {
990     lstfront.supprimer(ftref2);
991     delete ftref2;
992     }
993     if (ftref3!=NULL) {
994     lstfront.supprimer(ftref3);
995     delete ftref3;
996     }
997     if (ftref4!=NULL) {
998     lstfront.supprimer(ftref4);
999     delete ftref4;
1000     }
1001     if (ft1!=NULL)
1002     {
1003     mise_a_jour_voisin_front(ft1);
1004     lstfront.ajouter(ft1->get_front_voisin(0));
1005     lstfront.ajouter(ft1->get_front_voisin(1));
1006     lstfront.ajouter(ft1->get_front_voisin(2));
1007     }
1008     if (ft2!=NULL)
1009     {
1010     mise_a_jour_voisin_front(ft2);
1011     lstfront.ajouter(ft2->get_front_voisin(0));
1012     lstfront.ajouter(ft2->get_front_voisin(1));
1013     lstfront.ajouter(ft2->get_front_voisin(2));
1014     }
1015     if (ft3!=NULL)
1016     {
1017     mise_a_jour_voisin_front(ft3);
1018     lstfront.ajouter(ft3->get_front_voisin(0));
1019     lstfront.ajouter(ft3->get_front_voisin(1));
1020     lstfront.ajouter(ft3->get_front_voisin(2));
1021     }
1022     if (ft4!=NULL)
1023     {
1024     mise_a_jour_voisin_front(ft4);
1025     lstfront.ajouter(ft4->get_front_voisin(0));
1026     lstfront.ajouter(ft4->get_front_voisin(1));
1027     lstfront.ajouter(ft4->get_front_voisin(2));
1028     }
1029    
1030     for (int l=0;l<lstfront.get_nb();l++) mise_a_jour_voisin_front(lstfront.get(l));
1031    
1032     }
1033     }
1034    
1035     void MAILLEUR3D::insere_contrainte_triangle(MG_VOLUME* mgvol,TPL_LISTE_ENTITE<MG_TRIANGLE*> *lsttri)
1036     {
1037     int nbtri=lsttri->get_nb();
1038     for (int i=0;i<nbtri;i++)
1039     {
1040     M3D_TRIANGLE* tri=(M3D_TRIANGLE*)lsttri->get(i);
1041     double *xyz1=tri->get_noeud1()->get_coord();
1042     double *xyz2=tri->get_noeud2()->get_coord();
1043     double *xyz3=tri->get_noeud3()->get_coord();
1044     double xyz[3];
1045     xyz[0]=(xyz1[0]+xyz2[0]+xyz3[0])/3.;
1046     xyz[1]=(xyz1[1]+xyz2[1]+xyz3[1])/3.;
1047     xyz[2]=(xyz1[2]+xyz2[2]+xyz3[2])/3.;
1048     TPL_MAP_ENTITE<MG_FRONT_3D*> liste_trouvee;
1049     OT_VECTEUR_3D vec1(xyz,xyz1);
1050     OT_VECTEUR_3D vec2(xyz,xyz2);
1051     OT_VECTEUR_3D vec3(xyz,xyz3);
1052     double rayonrecherche=vec1.get_longueur();
1053     if (vec2.get_longueur()>rayonrecherche) rayonrecherche=vec2.get_longueur();
1054     if (vec3.get_longueur()>rayonrecherche) rayonrecherche=vec3.get_longueur();
1055     rayonrecherche=rayonrecherche*1.1;
1056     octree_de_front->rechercher(xyz[0],xyz[2],xyz[2],rayonrecherche,liste_trouvee);
1057     MG_NOEUD* noeudref1=NULL;
1058     MG_NOEUD* noeudref2=NULL;
1059     MG_NOEUD* noeudref3=NULL;
1060     int nbtrouve=liste_trouvee.get_nb();
1061     for (int j=0;j<nbtrouve;j++)
1062     {
1063     MG_FRONT_3D* ft=liste_trouvee.get(j);
1064     TPL_MAP_ENTITE<MG_NOEUD*> lstnoeud;
1065     lstnoeud.ajouter(ft->get_noeud1());
1066     lstnoeud.ajouter(ft->get_noeud2());
1067     lstnoeud.ajouter(ft->get_noeud3());
1068     int nb=lstnoeud.get_nb();
1069     for (int k=0;k<nb;k++)
1070     {
1071     double *xyztmp=lstnoeud.get(k)->get_coord();
1072     OT_VECTEUR_3D vec(xyztmp,xyz1);
1073     double longueur=vec.get_longueur();
1074     if (longueur<1e-10*rayonrecherche)
1075     {
1076     noeudref1=lstnoeud.get(k);
1077     lstnoeud.supprimer(lstnoeud.get(k));
1078     break;
1079     }
1080     }
1081     nb=lstnoeud.get_nb();
1082     for (int k=0;k<nb;k++)
1083     {
1084     double *xyztmp=lstnoeud.get(k)->get_coord();
1085     OT_VECTEUR_3D vec(xyztmp,xyz2);
1086     double longueur=vec.get_longueur();
1087     if (longueur<1e-10*rayonrecherche)
1088     {
1089     noeudref2=lstnoeud.get(k);
1090     lstnoeud.supprimer(lstnoeud.get(k));
1091     break;
1092     }
1093     }
1094     nb=lstnoeud.get_nb();
1095     for (int k=0;k<nb;k++)
1096     {
1097     double *xyztmp=lstnoeud.get(k)->get_coord();
1098     OT_VECTEUR_3D vec(xyztmp,xyz3);
1099     double longueur=vec.get_longueur();
1100     if (longueur<1e-10*rayonrecherche)
1101     {
1102     noeudref3=lstnoeud.get(k);
1103     lstnoeud.supprimer(lstnoeud.get(k));
1104     break;
1105     }
1106     }
1107    
1108     }
1109     if (noeudref1==NULL)
1110     {
1111     noeudref1=mg_maillage->ajouter_mg_noeud(mgvol,xyz1[0],xyz1[1],xyz1[2],IMPOSE);
1112     noeudref1->change_nouveau_numero(tri->get_noeud1()->get_nouveau_numero());
1113     }
1114     if (noeudref2==NULL)
1115     {
1116     noeudref2=mg_maillage->ajouter_mg_noeud(mgvol,xyz2[0],xyz2[1],xyz2[2],IMPOSE);
1117     noeudref2->change_nouveau_numero(tri->get_noeud2()->get_nouveau_numero());
1118     }
1119     if (noeudref3==NULL)
1120     {
1121     noeudref3=mg_maillage->ajouter_mg_noeud(mgvol,xyz3[0],xyz3[1],xyz3[2],IMPOSE);
1122     noeudref3->change_nouveau_numero(tri->get_noeud3()->get_nouveau_numero());
1123     }
1124    
1125     MG_SEGMENT* segment1=mg_maillage->get_mg_segment(noeudref1->get_id(),noeudref2->get_id());
1126     MG_SEGMENT* segment2=mg_maillage->get_mg_segment(noeudref1->get_id(),noeudref3->get_id());
1127     MG_SEGMENT* segment3=mg_maillage->get_mg_segment(noeudref2->get_id(),noeudref3->get_id());
1128     if (segment1==NULL) segment1=cree_segment(mgvol,noeudref1,noeudref2,IMPOSE);
1129     if (segment2==NULL) segment2=cree_segment(mgvol,noeudref1,noeudref3,IMPOSE);
1130     if (segment3==NULL) segment3=cree_segment(mgvol,noeudref2,noeudref3,IMPOSE);
1131     MG_TRIANGLE* nvtri=cree_triangle(mgvol,noeudref1,noeudref2,noeudref3,segment1,segment2,segment3,IMPOSE);
1132     MG_TRIANGLE* nvtri2=cree_triangle(mgvol,noeudref1,noeudref3,noeudref2,segment1,segment2,segment3,IMPOSE);
1133     //MG_TRIANGLE* nvtri=insere_triangle(mgvol,noeudref1,noeudref2,noeudref3,IMPOSE);
1134     //MG_TRIANGLE* nvtri2=insere_triangle(mgvol,noeudref1,noeudref3,noeudref2,IMPOSE);
1135     MG_FRONT_3D *ft=ajouter_front_courant(NONFORCE,noeudref1,noeudref2,noeudref3,nvtri);
1136     MG_FRONT_3D *ft2=ajouter_front_courant(NONFORCE,noeudref1,noeudref3,noeudref2,nvtri2);
1137     mise_a_jour_voisin_front(ft);
1138     mise_a_jour_voisin_front(ft2);
1139     mise_a_jour_voisin_front(ft->get_front_voisin(0));
1140     mise_a_jour_voisin_front(ft->get_front_voisin(1));
1141     mise_a_jour_voisin_front(ft->get_front_voisin(2));
1142     mise_a_jour_voisin_front(ft2->get_front_voisin(0));
1143     mise_a_jour_voisin_front(ft2->get_front_voisin(1));
1144     mise_a_jour_voisin_front(ft2->get_front_voisin(2));
1145    
1146     }
1147    
1148    
1149    
1150    
1151     }
1152