ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur_auto/src/mailleur2d.cpp
Revision: 632
Committed: Thu Jan 15 18:40:00 2015 UTC (10 years, 3 months ago) by francois
File size: 63029 byte(s)
Log Message:
Changement de l'espace de voisinage dans le mailleur 2D et 3D. On utilise un ntree (octree "anisotrope" qui ne se divise pas necessairement en 8). En 2D l'espace de voisinage est maintenant sur l'espace reelle au lieu d'être sur l'espace parametrique.
+Mise sous forme de parametres de certains choix stratégiques du mailleur

File Contents

# User Rev Content
1 5 //------------------------------------------------------------
2     //------------------------------------------------------------
3     // MAGiC
4 francois 155 // Jean Christophe Cuilli�re et Vincent FRANCOIS
5     // D�partement de G�nie M�canique - UQTR
6 5 //------------------------------------------------------------
7 francois 155 // 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 5 // des auteurs (contact : francois@uqtr.ca)
12     //------------------------------------------------------------
13     //------------------------------------------------------------
14     //
15     // mailleur2d.cpp
16     //
17     //------------------------------------------------------------
18     //------------------------------------------------------------
19     // COPYRIGHT 2000
20 francois 155 // Version du 02/03/2006 � 11H23
21 5 //------------------------------------------------------------
22     //------------------------------------------------------------
23    
24    
25     #include "gestionversion.h"
26     #include "mailleur2d.h"
27     #include "m3d_triangle.h"
28     #include "ot_mathematique.h"
29     #include "mg_gestionnaire.h"
30     //#include "message.h"
31     //#include "affiche.h"
32     #include <fstream>
33     #include <math.h>
34 francois 447 #include "mailleur2d_optimisation.h"
35 5
36    
37 francois 315
38 francois 494 MAILLEUR2D::MAILLEUR2D(MG_MAILLAGE* mgmai,MG_GEOMETRIE *mggeo,FCT_TAILLE* fct_taille,MG_FACE* mgface):MAILLEUR(false),mg_maillage(mgmai),mg_geometrie(mggeo),mg_face(mgface),metrique(fct_taille),niveau_optimisation(2)
39 5 {
40     }
41    
42    
43    
44     MAILLEUR2D::~MAILLEUR2D()
45     {
46     }
47    
48    
49    
50    
51 francois 210 void MAILLEUR2D::maille(MG_GROUPE_TOPOLOGIQUE* mggt)
52 5 {
53     if (mg_face!=NULL) maille(mg_face);
54     else
55     {
56 francois 210 TPL_MAP_ENTITE<MG_ELEMENT_TOPOLOGIQUE*> lst;
57     if (mggt!=NULL)
58     {
59     int nb=mggt->get_nb();
60     for (int i=0;i<nb;i++)
61     {
62     lst.ajouter(mggt->get(i));
63     mggt->get(i)->get_topologie_sousjacente(&lst);
64     }
65     }
66 5 int nb_face=mg_geometrie->get_nb_mg_face();
67     for (int i=0;i<nb_face;i++)
68     {
69     MG_FACE* mgface=mg_geometrie->get_mg_face(i);
70 francois 210 if (mggt!=NULL)
71     if (lst.existe(mgface)==0) continue;
72 5 maille(mgface);
73     }
74     }
75     }
76    
77    
78     void MAILLEUR2D::maille(MG_FACE* mgface,TPL_LISTE_ENTITE<MG_SEGMENT*> *lstseg,TPL_LISTE_ENTITE<MG_TRIANGLE*> *lsttri)
79    
80     {
81     periode_u=mgface->get_surface()->get_periode_u();
82     periode_v=mgface->get_surface()->get_periode_v();
83     decalage=new OT_DECALAGE_PARAMETRE(periode_u,periode_v);
84     //afficheur << MAILLAGEFACE << mgface->get_id() << endaff;
85     initialise_frontiere(mgface);
86 francois 632 cree_ntree(mgface);
87 5 initialise_front(mgface);
88     if (lstseg!=NULL) insere_contrainte_segment(mgface,lstseg);
89     if (lsttri!=NULL) insere_contrainte_triangle(mgface,lsttri);
90     progresse_front(mgface);
91 francois 447 MAILLEUR2D_OPTIMISATION opt(mg_maillage,niveau_optimisation);
92     opt.optimise(mgface);
93 francois 446
94 francois 632 delete ntree_de_segment;
95     //delete quadtree_de_frontiere;
96     delete ntree_de_front;
97 5 delete decalage;
98     segment_frontiere.vide();
99     noeud_frontiere.vide();
100     //ofstream o3("test_2D.mai",ios::out|ios::trunc);
101     //o3.precision(16);
102     //o3.setf(ios::showpoint);
103     //mg_maillage->enregistrer_sous_mesh_2D(o3);
104     }
105    
106    
107    
108    
109    
110     void MAILLEUR2D::initialise_frontiere(MG_FACE* mgface)
111     {
112 francois 632 umin=1e308;vmin=1e-308;
113     umax=-1e308;vmax=-1e-308;
114 5 int nb_boucle=mgface->get_nb_mg_boucle();
115     for (int i=0;i<nb_boucle;i++)
116     {
117     MG_BOUCLE* mgboucle=mgface->get_mg_boucle(i);
118     int nb_coarete=mgboucle->get_nb_mg_coarete();
119     for (int j=0;j<nb_coarete;j++)
120     {
121     MG_COARETE* coarete=mgboucle->get_mg_coarete(j);
122     MG_ARETE* mgarete=coarete->get_arete();
123     int nb_segment=mgarete->get_lien_maillage()->get_nb();
124     for (int k=0;k<nb_segment;k++)
125     {
126     refresh();
127     MG_SEGMENT* mgsegment=(MG_SEGMENT*)mgarete->get_lien_maillage()->get(k);
128     MG_SEGMENT* mgsegtemp=mg_maillage->get_mg_segmentid(mgsegment->get_id());
129     if (mgsegtemp==NULL) continue;
130     MG_NOEUD* mgnoeud=mgsegment->get_noeud1();
131     segment_frontiere.ajouter(mgsegment);
132 francois 61 if (!noeud_frontiere.est_dans_la_liste(mgnoeud))
133     {
134     double uv[2];
135     double *coo=mgnoeud->get_coord();
136     mgface->inverser(uv,coo,1e-4);
137     if (mgface->get_surface()->est_periodique_u()==1)
138     {
139 francois 315 double eps=1e-10*mgface->get_surface()->get_periode_u();
140     if (uv[0]<-eps) uv[0]=uv[0]+mgface->get_surface()->get_periode_u();
141     if (uv[0]>=mgface->get_surface()->get_periode_u()-eps) uv[0]=uv[0]-mgface->get_surface()->get_periode_u();
142 francois 61 }
143     if (mgface->get_surface()->est_periodique_v()==1)
144     {
145 francois 315 double eps=1e-10*mgface->get_surface()->get_periode_v();
146     if (uv[1]<-eps) uv[1]=uv[1]+mgface->get_surface()->get_periode_v();
147     if (uv[1]>=mgface->get_surface()->get_periode_v()-eps) uv[1]=uv[1]-mgface->get_surface()->get_periode_v();
148 francois 61 }
149 francois 632 if (uv[0]<umin) umin=uv[0];
150     if (uv[1]<vmin) vmin=uv[1];
151     if (uv[0]>umax) umax=uv[0];
152     if (uv[1]>vmax) vmax=uv[1];
153 francois 61 mgnoeud->change_u(uv[0]);
154     mgnoeud->change_v(uv[1]);
155     noeud_frontiere.ajouter(mgnoeud);
156     }
157 5 mgnoeud=mgsegment->get_noeud2();
158 francois 61 if (!noeud_frontiere.est_dans_la_liste(mgnoeud))
159 5 {
160     double uv[2];
161     double *coo=mgnoeud->get_coord();
162     mgface->inverser(uv,coo,1e-4);
163     if (mgface->get_surface()->est_periodique_u()==1)
164     {
165 francois 315 double eps=1e-10*mgface->get_surface()->get_periode_u();
166     if (uv[0]<-eps) uv[0]=uv[0]+mgface->get_surface()->get_periode_u();
167     if (uv[0]>=mgface->get_surface()->get_periode_u()-eps) uv[0]=uv[0]-mgface->get_surface()->get_periode_u();
168 5 }
169     if (mgface->get_surface()->est_periodique_v()==1)
170     {
171 francois 315 double eps=1e-10*mgface->get_surface()->get_periode_v();
172     if (uv[1]<-eps) uv[1]=uv[1]+mgface->get_surface()->get_periode_v();
173     if (uv[1]>=mgface->get_surface()->get_periode_v()-eps) uv[1]=uv[1]-mgface->get_surface()->get_periode_v();
174 5 }
175 francois 632 if (uv[0]<umin) umin=uv[0];
176     if (uv[1]<vmin) vmin=uv[1];
177     if (uv[0]>umax) umax=uv[0];
178     if (uv[1]>vmax) vmax=uv[1];
179 5 mgnoeud->change_u(uv[0]);
180     mgnoeud->change_v(uv[1]);
181     noeud_frontiere.ajouter(mgnoeud);
182     }
183     }
184    
185     }
186     }
187    
188     }
189    
190 francois 632 void MAILLEUR2D::cree_ntree(MG_FACE* mgface)
191 5 {
192     int nb_noeud=noeud_frontiere.get_nb();
193 francois 632 int nb_echan=(int)param.get_valeur((char*)"Echantillon_face");
194     BOITE_3D boiteface=mgface->get_boite_3D(nb_echan);
195     boiteface.change_grosseur(1.1);
196     BOITE_3D boitemetrique=metrique->get_boite_3D();
197     BOITE_3D boite=boiteface&boitemetrique;
198     boite=boitemetrique;
199     double xmin=boite.get_xmin();
200     double ymin=boite.get_ymin();
201     double zmin=boite.get_zmin();
202     double xmax=boite.get_xmax();
203     double ymax=boite.get_ymax();
204     double zmax=boite.get_zmax();
205     ntree_de_segment=new TPL_NTREE_FCT<MG_SEGMENT*,FCT_TAILLE>;
206     ntree_de_front=new TPL_OCTREE<MG_FRONT_2D*,MG_NOEUD*>;
207     ntree_de_segment->initialiser(*metrique,xmin,ymin,zmin,xmax,ymax,zmax);
208     ntree_de_front->initialiser(ntree_de_segment);
209 5 int nb_segment=segment_frontiere.get_nb();
210     for (int j=0;j<nb_segment;j++)
211     {
212 francois 632 MG_SEGMENT* seg=segment_frontiere.get(j);
213     ntree_de_segment->inserer(segment_frontiere.get(j));
214 5 }
215     }
216    
217    
218     void MAILLEUR2D::initialise_front(MG_FACE* mgface)
219     {
220     MG_FRONT_2D::initialise_compteur_id();
221     int nb_boucle=mgface->get_nb_mg_boucle();
222     for (int iboucle=0;iboucle<nb_boucle;iboucle++)
223     {
224     unsigned int nb_front_avant_cette_boucle=get_nb_front(front_courant);
225     MG_BOUCLE* mgboucle=mgface->get_mg_boucle(iboucle);
226     unsigned int nbcoaretetotale=mgboucle->get_nb_mg_coarete();
227     MG_COARETE* mgcoarete=mgboucle->get_mg_coarete(0);
228     MG_NOEUD* noeud_de_depart;
229     int bon_noeud=0;
230     int numnoeud=0;
231     while (bon_noeud==0)
232     {
233     if (mgcoarete->get_orientation()==MEME_SENS) noeud_de_depart=(MG_NOEUD*)mgcoarete->get_arete()->get_cosommet1()->get_sommet()->get_lien_maillage()->get(numnoeud);
234     else noeud_de_depart=(MG_NOEUD*)mgcoarete->get_arete()->get_cosommet2()->get_sommet()->get_lien_maillage()->get(numnoeud);
235     MG_NOEUD* noeudtemp=mg_maillage->get_mg_noeudid(noeud_de_depart->get_id());
236     if (noeudtemp==NULL) numnoeud++; else bon_noeud=1;
237     }
238     MG_FRONT_2D* premier_front;
239     MG_NOEUD* noeud_courant=noeud_de_depart;
240     MG_SEGMENT* segment_courant=NULL;
241     MG_FRONT_2D* front_precedent=NULL;
242     unsigned int nbcoaretetraite=0;
243     do
244     {
245     nbcoaretetraite++;
246     MG_NOEUD* noeud_d_arrivee;
247     int bon_noeud=0;
248     int numnoeud=0;
249     while (bon_noeud==0)
250     {
251     if (mgcoarete->get_orientation()==MEME_SENS) noeud_d_arrivee=(MG_NOEUD*)mgcoarete->get_arete()->get_cosommet2()->get_sommet()->get_lien_maillage()->get(numnoeud);
252     else noeud_d_arrivee=(MG_NOEUD*)mgcoarete->get_arete()->get_cosommet1()->get_sommet()->get_lien_maillage()->get(numnoeud);
253     MG_NOEUD* noeudtemp=mg_maillage->get_mg_noeudid(noeud_d_arrivee->get_id());
254     if (noeudtemp==NULL) numnoeud++; else bon_noeud=1;
255     }
256     MG_ARETE* mgarete=mgcoarete->get_arete();
257     int passe_aretefermee=0;
258     do
259     {
260     int trouve=0;
261     int i;
262     if ((mgcoarete->get_arete()->get_cosommet1()->get_sommet()==mgcoarete->get_arete()->get_cosommet2()->get_sommet())&&(passe_aretefermee==0))
263     {
264     passe_aretefermee=1;
265     MG_SEGMENT* mgsegment_depart1=NULL;
266     MG_SEGMENT* mgsegment_depart2=NULL;
267     TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR it;
268     MG_SEGMENT* mgsegment=(MG_SEGMENT*)mgarete->get_lien_maillage()->get_premier(it);
269     do
270     {
271     if ((mgsegment->get_noeud1()==noeud_courant) || (mgsegment->get_noeud2()==noeud_courant)) mgsegment_depart1=mgsegment;
272     mgsegment=(MG_SEGMENT*)mgarete->get_lien_maillage()->get_suivant(it);
273     }
274     while (mgsegment_depart1==NULL);
275     do
276     {
277     if ((mgsegment->get_noeud1()==noeud_courant) || (mgsegment->get_noeud2()==noeud_courant)) mgsegment_depart2=mgsegment;
278     mgsegment=(MG_SEGMENT*)mgarete->get_lien_maillage()->get_suivant(it);
279     }
280     while (mgsegment_depart2==NULL);
281     MG_NOEUD* noeud1_2;
282     MG_NOEUD* noeud2_2;
283     if (mgsegment_depart1->get_noeud1()==noeud_courant) noeud1_2=mgsegment_depart1->get_noeud2(); else noeud1_2=mgsegment_depart1->get_noeud1();
284     if (mgsegment_depart2->get_noeud1()==noeud_courant) noeud2_2=mgsegment_depart2->get_noeud2(); else noeud2_2=mgsegment_depart2->get_noeud1();
285     OT_VECTEUR_3D vecteur1(noeud1_2->get_x()-noeud_courant->get_x(),noeud1_2->get_y()-noeud_courant->get_y(),noeud1_2->get_z()-noeud_courant->get_z());
286     double coo[3];
287     mgarete->deriver(mgarete->get_tmin(),coo);
288     OT_VECTEUR_3D tangeante(coo[0]*mgcoarete->get_orientation(),coo[1]*mgcoarete->get_orientation(),coo[2]*mgcoarete->get_orientation());
289     MG_NOEUD* noeud_suivant;
290     if (vecteur1*tangeante>0.)
291     {
292     mgsegment=mgsegment_depart1;
293     noeud_suivant=noeud1_2;
294     }
295     else
296     {
297     mgsegment=mgsegment_depart2;
298     noeud_suivant=noeud2_2;
299     }
300     MG_FRONT_2D* front=ajouter_front(front_courant,noeud_courant,noeud_suivant,mgsegment);
301 francois 632 ntree_de_segment->inserer(mgsegment);
302 5 if (get_nb_front(front_courant)==nb_front_avant_cette_boucle+1) premier_front=front;
303     front->changer_front_precedent(front_precedent);
304     front_precedent=front;
305     noeud_courant=noeud_suivant;
306     segment_courant=mgsegment;
307     }
308     else
309     {
310     trouve=0;
311     TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR it;
312     MG_SEGMENT* mgsegment=(MG_SEGMENT*)mgarete->get_lien_maillage()->get_premier(it);
313     do
314     {
315     if ((mgsegment->get_noeud1()==noeud_courant) || (mgsegment->get_noeud2()==noeud_courant))
316     if (mgsegment!=segment_courant)
317     {
318     trouve=1;
319     MG_NOEUD* noeud_suivant;
320     if (mgsegment->get_noeud1()==noeud_courant) noeud_suivant=mgsegment->get_noeud2();
321     else if (mgsegment->get_noeud2()==noeud_courant) noeud_suivant=mgsegment->get_noeud1();
322     MG_FRONT_2D* front=ajouter_front(front_courant,noeud_courant,noeud_suivant,mgsegment);
323 francois 632 ntree_de_segment->inserer(mgsegment);
324 5 if (get_nb_front(front_courant)==nb_front_avant_cette_boucle+1) premier_front=front;
325     else front_precedent->changer_front_suivant(front);
326     front->changer_front_precedent(front_precedent);
327     front_precedent=front;
328     noeud_courant=noeud_suivant;
329     segment_courant=mgsegment;
330     }
331     mgsegment=(MG_SEGMENT*)mgarete->get_lien_maillage()->get_suivant(it);
332     }
333     while (trouve==0);
334     }
335     }
336     while (noeud_courant!=noeud_d_arrivee);
337     int trouve=0;
338     MG_SOMMET* mgsommet;
339     if (mgcoarete->get_orientation()==MEME_SENS) mgsommet=mgcoarete->get_arete()->get_cosommet2()->get_sommet();
340     else mgsommet=mgcoarete->get_arete()->get_cosommet1()->get_sommet();
341     double xyz[3];
342     mgsommet->get_point()->evaluer(xyz);
343     int nbcoaretesuivante=mgboucle->get_nb_mg_coarete();
344     MG_COARETE* mgcoaretetmp=NULL;
345     double angleref=0.;
346     double angleref2=0.;
347     for (int i=0;i<nbcoaretesuivante;i++)
348     {
349     MG_COARETE* mgcoarete_suivante=mgboucle->get_mg_coarete(i);
350     MG_SOMMET* mgsommet_suivant;
351     if (mgcoarete_suivante->get_orientation()==MEME_SENS) mgsommet_suivant=mgcoarete_suivante->get_arete()->get_cosommet1()->get_sommet();
352     else mgsommet_suivant=mgcoarete_suivante->get_arete()->get_cosommet2()->get_sommet();
353     if (mgsommet==mgsommet_suivant)
354     {
355     if (trouve==0)
356     {
357     trouve=1;
358     mgcoaretetmp=mgcoarete_suivante;
359     }
360     else
361     {
362     //angleref
363     double uv[2];
364     mgface->inverser(uv,xyz);
365     double normal[3];
366     mgface->calcul_normale_unitaire(uv,normal);
367     double dir1[3];
368     double dir2[3];
369     MG_COSOMMET* mgcosommet;
370     if (mgcoarete->get_orientation()==MEME_SENS) mgcosommet=mgcoarete->get_arete()->get_cosommet2();
371     else mgcosommet=mgcoarete->get_arete()->get_cosommet1();
372     mgcoarete->get_arete()->deriver(mgcosommet->get_t(),dir1);
373     if (mgcoaretetmp->get_orientation()==MEME_SENS) mgcosommet=mgcoaretetmp->get_arete()->get_cosommet1();
374     else mgcosommet=mgcoaretetmp->get_arete()->get_cosommet2();
375     mgcoaretetmp->get_arete()->deriver(mgcosommet->get_t(),dir2);
376     int sens1=mgcoarete->get_orientation();
377     int sens2=mgcoaretetmp->get_orientation();
378     OT_VECTEUR_3D vec1(-dir1[0]*sens1,-dir1[1]*sens1,-dir1[2]*sens1);
379     OT_VECTEUR_3D vec2(dir2[0]*sens2,dir2[1]*sens2,dir2[2]*sens2);
380     vec1.norme();
381     vec2.norme();
382     OT_VECTEUR_3D nor=normal&vec1;
383     double cs=vec1*vec2;
384     double sn=-1*(vec2*nor);
385     angleref=acos(cs);
386     if (sn<0) angleref=-angleref;
387     if (angleref<0) angleref=angleref+2.*M_PI;
388     if (mgcoarete->get_arete()==mgcoaretetmp->get_arete()) angleref=2*M_PI;
389     if (angleref<1e-10)
390     {
391     double dt=1e-3*(mgcoaretetmp->get_arete()->get_tmax()-mgcoaretetmp->get_arete()->get_tmax());
392     mgcoaretetmp->get_arete()->deriver(mgcosommet->get_t()+dt,dir2);
393     OT_VECTEUR_3D vec2(dir2[0]*sens2,dir2[1]*sens2,dir2[2]*sens2);
394     vec2.norme();
395     double cs=vec1*vec2;
396     double sn=-1*(vec2*nor);
397     angleref=acos(cs);
398     if (sn<0) angleref2=-angleref2;
399     if (angleref<0) angleref2=angleref2+2.*M_PI;
400     }
401     //angle
402     if (mgcoarete->get_orientation()==MEME_SENS) mgcosommet=mgcoarete->get_arete()->get_cosommet2();
403     else mgcosommet=mgcoarete->get_arete()->get_cosommet1();
404     mgcoarete->get_arete()->deriver(mgcosommet->get_t(),dir1);
405     if (mgcoarete_suivante->get_orientation()==MEME_SENS) mgcosommet=mgcoarete_suivante->get_arete()->get_cosommet1();
406     else mgcosommet=mgcoarete_suivante->get_arete()->get_cosommet2();
407     mgcoarete_suivante->get_arete()->deriver(mgcosommet->get_t(),dir2);
408     sens1=mgcoarete->get_orientation();
409     sens2=mgcoarete_suivante->get_orientation();
410     vec1.change_x(-dir1[0]*sens1);
411     vec1.change_y(-dir1[1]*sens1);
412     vec1.change_z(-dir1[2]*sens1);
413     vec2.change_x(dir2[0]*sens2);
414     vec2.change_y(dir2[1]*sens2);
415     vec2.change_z(dir2[2]*sens2);
416     vec1.norme();
417     vec2.norme();
418     nor=normal&vec1;
419     cs=vec1*vec2;
420     sn=-1*(vec2*nor);
421     double angle=acos(cs);
422     if (sn<0) angle=-angle;
423     if (angle<0) angle=angle+2.*M_PI;
424     if (mgcoarete->get_arete()==mgcoarete_suivante->get_arete()) angle=2*M_PI;
425     if (angle<angleref)
426     {
427     angleref=angle;
428     mgcoaretetmp=mgcoarete_suivante;
429     }
430     else if ((angle<1e-10) && (angleref<1e-10))
431     {
432     double dt=1e-3*(mgcoarete_suivante->get_arete()->get_tmax()-mgcoarete_suivante->get_arete()->get_tmax());
433     mgcoarete_suivante->get_arete()->deriver(mgcosommet->get_t()+dt,dir2);
434     OT_VECTEUR_3D vec2(dir2[0]*sens2,dir2[1]*sens2,dir2[2]*sens2);
435     vec2.norme();
436     double cs=vec1*vec2;
437     double sn=-1*(vec2*nor);
438     angle=acos(cs);
439     if (sn<0) angle=-angle;
440     if (angle<0) angle=angle+2.*M_PI;
441     if (angle<angleref2)
442     mgcoaretetmp=mgcoarete_suivante;
443    
444     }
445    
446     }
447     }
448     }
449     mgcoarete=mgcoaretetmp;
450     }
451     while ((noeud_courant!=noeud_de_depart) || (nbcoaretetraite!=nbcoaretetotale));
452     front_precedent->changer_front_suivant(premier_front);
453     premier_front->changer_front_precedent(front_precedent);
454     }
455    
456     }
457    
458     void MAILLEUR2D::progresse_front(MG_FACE* mgface)
459     {
460     int compteur=0;
461     while (get_nb_front(front_courant)!=0)
462     {
463 francois 33 //ofstream o3("c:\\void_1D.mai",ios::out|ios::trunc);
464     //o3.precision(16);
465     //o3.setf(ios::showpoint);
466     //mg_maillage->enregistrer_sous_mesh_1D(o3);
467     //ofstream o4("c:\\void_2D.mai",ios::out|ios::trunc);
468     //o4.precision(16);
469     //o4.setf(ios::showpoint);
470     //mg_maillage->enregistrer_sous_mesh_2D(o4);
471 francois 35 //mg_maillage->get_gestionnaire()->enregistrer("c:\\void.magic");
472 5 compteur++;
473 francois 315 /*if (mg_maillage->get_nb_mg_triangle()+558>=1160)
474     {int i=0;
475     LISTE_MG_NOEUD::iterator it;for (MG_NOEUD* noeud=mg_maillage->get_premier_noeud(it);noeud;noeud=mg_maillage->get_suivant_noeud(it))
476     {
477     noeud->change_nouveau_numero(i+1);
478     std::cout << i+1 << " " << noeud->get_u() << " " << noeud->get_v() << std::endl; ;
479     i++;
480     }
481     mg_maillage->get_gestionnaire()->enregistrer("toto.magic");
482     }*/
483 francois 632 //mg_maillage->get_gestionnaire()->enregistrer("toto.magic");
484 5 refresh();
485     MG_FRONT_2D* ft=get_front(front_courant,0);
486     MG_FRONT_2D* ftp=ft->get_front_precedent();
487     MG_FRONT_2D* fts=ft->get_front_suivant();
488     double eps=0.0001*ft->get_segment()->get_longueur();
489     MG_NOEUD* nop=ftp->get_noeud1();
490     MG_NOEUD* no1=ft->get_noeud1();
491     MG_NOEUD* no2=ft->get_noeud2();
492     MG_NOEUD* nos=fts->get_noeud2();
493     /* angle precedent et angle suivant */
494     OT_VECTEUR_3D n1n2(no1->get_coord(),no2->get_coord());
495     OT_VECTEUR_3D n2s(no2->get_coord(),nos->get_coord());
496     OT_VECTEUR_3D pn1(nop->get_coord(),no1->get_coord());
497     n1n2.norme();
498     n2s.norme();
499     pn1.norme();
500     double coo[3];
501     double uv[2];
502     uv[0]=no1->get_u();
503     uv[1]=no1->get_v();
504     mgface->calcul_normale_unitaire(uv,coo);
505     OT_VECTEUR_3D nplan(coo);
506     OT_VECTEUR_3D n=nplan&n1n2;
507     n.norme();
508     double pcp=(-1)*(n1n2*pn1);
509     double psp=(-1)*(n*pn1);
510     int anglep;
511     if ((pcp>=0.1786481777)&&(psp>eps)) anglep=1; else anglep=0;
512     uv[0]=no2->get_u();
513     uv[1]=no2->get_v();
514     mgface->calcul_normale_unitaire(uv,coo);
515     OT_VECTEUR_3D nplan2(coo);
516     n=nplan2&n1n2;
517     n.norme();
518     double pcs=(-1.)*(n1n2*n2s);
519     double pss=n*n2s;
520     int angles;
521     if ((pcs>=0.1786481777)&&(pss>eps)) angles=1; else angles=0;
522     int type_cas_front[7];
523     if ( (ftp==fts->get_front_suivant()) && (no1!=nos) && (no2!=nop) ) type_cas_front[CAS_FRONT_3]=1; else type_cas_front[CAS_FRONT_3]=0;
524     if (ftp->get_front_precedent()==fts->get_front_suivant()) type_cas_front[CAS_FRONT_4]=1; else type_cas_front[CAS_FRONT_4]=0;
525     if ( (anglep==1) && (angles==1) ) type_cas_front[CAS_FERME_CAVITE]=1; else type_cas_front[CAS_FERME_CAVITE]=0;
526     if ( anglep==1 ) type_cas_front[CAS_FERME_CAVITE_P]=1; else type_cas_front[CAS_FERME_CAVITE_P]=0;
527     if ( angles==1 ) type_cas_front[CAS_FERME_CAVITE_S]=1; else type_cas_front[CAS_FERME_CAVITE_S]=0;
528     type_cas_front[CAS_GENERAL]=1;
529     int solution=PASTROUVE;
530     if (type_cas_front[CAS_FRONT_3]) solution=traite_front(CAS_FRONT_3,mgface,ftp,ft,fts);
531     if ((solution==PASTROUVE) && (type_cas_front[CAS_FRONT_4])) solution=traite_front(CAS_FRONT_4,mgface,ftp,ft,fts);
532     if ((solution==PASTROUVE) && (type_cas_front[CAS_FERME_CAVITE])) solution=traite_front(CAS_FERME_CAVITE,mgface,ftp,ft,fts);
533     if ((solution==PASTROUVE) && (type_cas_front[CAS_FERME_CAVITE_P])) solution=traite_front(CAS_FERME_CAVITE_P,mgface,ftp,ft,fts);
534     if ((solution==PASTROUVE) && (type_cas_front[CAS_FERME_CAVITE_S])) solution=traite_front(CAS_FERME_CAVITE_S,mgface,ftp,ft,fts);
535     if ((solution==PASTROUVE) && (type_cas_front[CAS_GENERAL]))
536     {
537     solution=traite_front(CAS_GENERAL,mgface,ftp,ft,fts);
538     if ((solution!=PASTROUVE)&&(solution!=TROUVE))
539     if (!type_cas_front[solution]) solution=traite_front(solution,mgface,ftp,ft,fts); else solution=PASTROUVE;
540     if (solution==PASTROUVE)
541     {
542     echange_de_front(front_courant,front_attente,ft);
543     ft->incremente_ifail();
544     solution=TROUVE;
545     }
546     }
547     if (get_nb_front(front_courant)==0) front_courant.swap(front_attente);
548     }
549     }
550    
551    
552    
553 francois 447
554 5 void MAILLEUR2D::insere_contrainte_segment(MG_FACE* mgface,TPL_LISTE_ENTITE<MG_SEGMENT*> *lstseg)
555     {
556     int nbseg=lstseg->get_nb();
557     for (int i=0;i<nbseg;i++)
558     {
559     MG_SEGMENT* segori=lstseg->get(i);
560     MG_NOEUD* noeudori1=segori->get_noeud1();
561     MG_NOEUD* noeudori2=segori->get_noeud2();
562     double uv1[2];
563     double uv2[2];
564     double *coo=noeudori1->get_coord();
565     mgface->inverser(uv1,coo,1e-10);
566     if (mgface->get_surface()->est_periodique_u()==1)
567     {
568     if (uv1[0]<0.0) uv1[0]=uv1[0]+mgface->get_surface()->get_periode_u();
569     if (uv1[0]>=mgface->get_surface()->get_periode_u()) uv1[0]=uv1[0]-mgface->get_surface()->get_periode_u();
570     }
571     if (mgface->get_surface()->est_periodique_v()==1)
572     {
573     if (uv1[1]<0.0) uv1[1]=uv1[1]+mgface->get_surface()->get_periode_v();
574     if (uv1[1]>=mgface->get_surface()->get_periode_v()) uv1[1]=uv1[1]-mgface->get_surface()->get_periode_v();
575     }
576     coo=noeudori2->get_coord();
577     mgface->inverser(uv2,coo,1e-10);
578     if (mgface->get_surface()->est_periodique_u()==1)
579     {
580     if (uv2[0]<0.0) uv2[0]=uv2[0]+mgface->get_surface()->get_periode_u();
581     if (uv2[0]>=mgface->get_surface()->get_periode_u()) uv2[0]=uv2[0]-mgface->get_surface()->get_periode_u();
582     }
583     if (mgface->get_surface()->est_periodique_v()==1)
584     {
585     if (uv2[1]<0.0) uv2[1]=uv2[1]+mgface->get_surface()->get_periode_v();
586     if (uv2[1]>=mgface->get_surface()->get_periode_v()) uv2[1]=uv2[1]-mgface->get_surface()->get_periode_v();
587     }
588    
589     TPL_MAP_ENTITE<MG_FRONT_2D*> liste_trouvee;
590 francois 632 double *xyz1=noeudori1->get_coord();
591     double *xyz2=noeudori2->get_coord();
592     double xyz[3];
593     xyz[0]=0.5*(xyz1[0]+xyz2[0]);
594     xyz[1]=0.5*(xyz1[1]+xyz2[1]);
595     xyz[2]=0.5*(xyz1[2]+xyz2[2]);
596     OT_VECTEUR_3D vec(xyz1,xyz);
597     //double rayonrecherche=sqrt((u-uv1[0])*(u-uv1[0])+(v-uv1[1])*(v-uv1[1]))*1.1;
598     double rayonrecherche=vec.get_longueur()*1.1;
599     ntree_de_front->rechercher(xyz[0],xyz[1],xyz[2],rayonrecherche,liste_trouvee);
600    
601     double u=0.5*(uv1[0]+uv2[0]);
602 5 double v=0.5*(uv1[1]+uv2[1]);
603     int nb_entite=liste_trouvee.get_nb();
604     double du=decalage->calcul_decalage_parametre_u(u);
605     double dv=decalage->calcul_decalage_parametre_v(v);
606     double u1=decalage->decalage_parametre_u(uv1[0],du);
607     double v1=decalage->decalage_parametre_v(uv1[1],dv);
608     double u2=decalage->decalage_parametre_u(uv2[0],du);
609     double v2=decalage->decalage_parametre_v(uv2[1],dv);
610     MG_FRONT_2D *ref1=NULL,*ref2=NULL;
611     double angle1,angle2;
612     for(int j=0;j<nb_entite;j++)
613     {
614     MG_FRONT_2D *ft=liste_trouvee.get(j);
615     MG_NOEUD* noeudfront=ft->get_noeud2();
616     double uvfront[2];
617     uvfront[0]=noeudfront->get_u();
618     uvfront[1]=noeudfront->get_v();
619     double ufront=decalage->decalage_parametre_u(uvfront[0],du);
620     double vfront=decalage->decalage_parametre_v(uvfront[1],dv);
621     //noeud 1 du segment a inserer
622 francois 632 double *xyzfront=noeudfront->get_coord();
623     //mgface->evaluer(uvfront,xyzfront);
624     OT_VECTEUR_3D vec(xyz1,xyzfront);
625     double dist=vec.get_longueur();
626     //double dist=sqrt((ufront-u1)*(ufront-u1)+(vfront-v1)*(vfront-v1));
627 francois 61 if (dist<1e-5*rayonrecherche)
628 5 {
629     if (ref1==NULL)
630     {
631     ref1=ft;
632     OT_VECTEUR_3D vecteur_baseu(u2-u1,v2-v1,0.);
633     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
634     double uvfrontautre[2];
635     uvfrontautre[0]=noeudfrontautre->get_u();
636     uvfrontautre[1]=noeudfrontautre->get_v();
637     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
638     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
639     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
640     vecteur_baseu.norme();
641     vecteur_front.norme();
642     double cosangle=vecteur_baseu*vecteur_front;
643     if (cosangle>1.) cosangle=1.;
644     if (cosangle<-1.) cosangle=(-1.);
645     angle1=acos(cosangle);
646     }
647     else
648     {
649     OT_VECTEUR_3D vecteur_baseu(u2-u1,v2-v1,0.);
650     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
651     double uvfrontautre[2];
652     uvfrontautre[0]=noeudfrontautre->get_u();
653     uvfrontautre[1]=noeudfrontautre->get_v();
654     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
655     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
656     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
657     vecteur_baseu.norme();
658     vecteur_front.norme();
659     double cosangle=vecteur_baseu*vecteur_front;
660     if (cosangle>1.) cosangle=1.;
661     if (cosangle<-1.) cosangle=(-1.);
662     double angle=acos(cosangle);
663     if (angle<angle1) {angle1=angle;ref1=ft;}
664     }
665     }
666     //noeud 2 du segment a inserer
667 francois 632 OT_VECTEUR_3D vec2(xyz2,xyzfront);
668     dist=vec2.get_longueur();
669     //dist=sqrt((ufront-u2)*(ufront-u2)+(vfront-v2)*(vfront-v2));
670 francois 61 if (dist<1e-5*rayonrecherche)
671 5 {
672     if (ref2==NULL)
673     {
674     ref2=ft;
675     OT_VECTEUR_3D vecteur_baseu(u1-u2,v1-v2,0.);
676     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
677     double uvfrontautre[2];
678     uvfrontautre[0]=noeudfrontautre->get_u();
679     uvfrontautre[1]=noeudfrontautre->get_v();
680     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
681     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
682     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
683     vecteur_baseu.norme();
684     vecteur_front.norme();
685     double cosangle=vecteur_baseu*vecteur_front;
686     if (cosangle>1.) cosangle=1.;
687     if (cosangle<-1.) cosangle=(-1.);
688     angle2=acos(cosangle);
689     }
690     else
691     {
692     OT_VECTEUR_3D vecteur_baseu(u1-u2,v1-v2,0.);
693     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
694     double uvfrontautre[2];
695     uvfrontautre[0]=noeudfrontautre->get_u();
696     uvfrontautre[1]=noeudfrontautre->get_v();
697     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
698     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
699     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
700     vecteur_baseu.norme();
701     vecteur_front.norme();
702     double cosangle=vecteur_baseu*vecteur_front;
703     if (cosangle>1.) cosangle=1.;
704     if (cosangle<-1.) cosangle=(-1.);
705     double angle=acos(cosangle);
706     if (angle<angle2) {angle1=angle;ref2=ft;}
707     }
708     }
709    
710     }
711     // creation du segemnt
712     MG_NOEUD* noeud1;
713     MG_NOEUD* noeud2;
714     if (ref1!=NULL) noeud1=ref1->get_noeud2();
715     else
716     {
717     double *coord=noeudori1->get_coord();
718 francois 35 noeud1=mg_maillage->ajouter_mg_noeud(mgface,coord[0],coord[1],coord[2],IMPOSE);
719 5 noeud1->change_u(uv1[0]);
720     noeud1->change_v(uv1[1]);
721 francois 447 noeud1->change_nouveau_numero(noeudori1->get_nouveau_numero());
722 francois 155 }
723 5 if (ref2!=NULL) noeud2=ref2->get_noeud2();
724     else
725     {
726     double *coord=noeudori2->get_coord();
727 francois 35 noeud2=mg_maillage->ajouter_mg_noeud(mgface,coord[0],coord[1],coord[2],IMPOSE);
728 5 noeud2->change_u(uv2[0]);
729     noeud2->change_v(uv2[1]);
730 francois 155 noeud2->change_nouveau_numero(noeudori2->get_nouveau_numero());
731     }
732 francois 35 MG_SEGMENT* mgsegment=mg_maillage->ajouter_mg_segment(mgface,noeud1,noeud2,IMPOSE);
733 5 MG_FRONT_2D *nv_ft=ajouter_front(front_courant,noeud1,noeud2,mgsegment);
734     MG_FRONT_2D *nv_ft2=ajouter_front(front_courant,noeud2,noeud1,mgsegment);
735     if (ref1==NULL)
736     {
737     nv_ft->changer_front_precedent(nv_ft2);
738     nv_ft2->changer_front_suivant(nv_ft);
739     }
740     else
741     {
742     MG_FRONT_2D *tmp=ref1->get_front_suivant();
743     ref1->changer_front_suivant(nv_ft);
744     nv_ft->changer_front_precedent(ref1);
745     nv_ft2->changer_front_suivant(tmp);
746     tmp->changer_front_precedent(nv_ft2);
747     }
748     if (ref2==NULL)
749     {
750     nv_ft->changer_front_suivant(nv_ft2);
751     nv_ft2->changer_front_precedent(nv_ft);
752     }
753     else
754     {
755     MG_FRONT_2D *tmp=ref2->get_front_suivant();
756     ref2->changer_front_suivant(nv_ft2);
757     nv_ft2->changer_front_precedent(ref2);
758     nv_ft->changer_front_suivant(tmp);
759     tmp->changer_front_precedent(nv_ft);
760     }
761    
762     }
763     }
764    
765     void MAILLEUR2D::insere_contrainte_triangle(MG_FACE* mgface,TPL_LISTE_ENTITE<MG_TRIANGLE*> *lsttri)
766     {
767     int nbtri=lsttri->get_nb();
768     for (int i=0;i<nbtri;i++)
769     {
770     MG_TRIANGLE* triori=lsttri->get(i);
771     MG_NOEUD* noeudori1=triori->get_noeud1();
772     MG_NOEUD* noeudori2=triori->get_noeud2();
773     MG_NOEUD* noeudori3=triori->get_noeud3();
774     double uv1[2];
775     double uv2[2];
776     double uv3[2];
777     double *coo=noeudori1->get_coord();
778     mgface->inverser(uv1,coo,1e-10);
779     if (mgface->get_surface()->est_periodique_u()==1)
780     {
781     if (uv1[0]<0.0) uv1[0]=uv1[0]+mgface->get_surface()->get_periode_u();
782     if (uv1[0]>=mgface->get_surface()->get_periode_u()) uv1[0]=uv1[0]-mgface->get_surface()->get_periode_u();
783     }
784     if (mgface->get_surface()->est_periodique_v()==1)
785     {
786     if (uv1[1]<0.0) uv1[1]=uv1[1]+mgface->get_surface()->get_periode_v();
787     if (uv1[1]>=mgface->get_surface()->get_periode_v()) uv1[1]=uv1[1]-mgface->get_surface()->get_periode_v();
788     }
789     coo=noeudori2->get_coord();
790     mgface->inverser(uv2,coo,1e-10);
791     if (mgface->get_surface()->est_periodique_u()==1)
792     {
793     if (uv2[0]<0.0) uv2[0]=uv2[0]+mgface->get_surface()->get_periode_u();
794     if (uv2[0]>=mgface->get_surface()->get_periode_u()) uv2[0]=uv2[0]-mgface->get_surface()->get_periode_u();
795     }
796     if (mgface->get_surface()->est_periodique_v()==1)
797     {
798     if (uv2[1]<0.0) uv2[1]=uv2[1]+mgface->get_surface()->get_periode_v();
799     if (uv2[1]>=mgface->get_surface()->get_periode_v()) uv2[1]=uv2[1]-mgface->get_surface()->get_periode_v();
800     }
801     coo=noeudori3->get_coord();
802     mgface->inverser(uv3,coo,1e-10);
803     if (mgface->get_surface()->est_periodique_u()==1)
804     {
805 francois 33 if (uv3[0]<0.0) uv3[0]=uv3[0]+mgface->get_surface()->get_periode_u();
806     if (uv3[0]>=mgface->get_surface()->get_periode_u()) uv3[0]=uv3[0]-mgface->get_surface()->get_periode_u();
807 5 }
808     if (mgface->get_surface()->est_periodique_v()==1)
809     {
810     if (uv3[1]<0.0) uv3[1]=uv3[1]+mgface->get_surface()->get_periode_v();
811     if (uv3[1]>=mgface->get_surface()->get_periode_v()) uv3[1]=uv3[1]-mgface->get_surface()->get_periode_v();
812     }
813     TPL_MAP_ENTITE<MG_FRONT_2D*> liste_trouvee;
814     double u=0.3333333333333333*(uv1[0]+uv2[0]+uv3[0]);
815     double v=0.3333333333333333*(uv1[1]+uv2[1]+uv3[1]);
816 francois 33 double rayonrecherche1=sqrt((u-uv1[0])*(u-uv1[0])+(v-uv1[1])*(v-uv1[1]))*2.;
817     double rayonrecherche2=sqrt((u-uv2[0])*(u-uv2[0])+(v-uv2[1])*(v-uv2[1]))*2.;
818     double rayonrecherche3=sqrt((u-uv3[0])*(u-uv3[0])+(v-uv3[1])*(v-uv3[1]))*2.;
819     double rayonrecherche=rayonrecherche1;
820     if (rayonrecherche2>rayonrecherche) rayonrecherche=rayonrecherche2;
821     if (rayonrecherche3>rayonrecherche) rayonrecherche=rayonrecherche3;
822 francois 632 double *xyz1=noeudori1->get_coord();
823     double *xyz2=noeudori2->get_coord();
824     double *xyz3=noeudori3->get_coord();
825     double xyz[3];
826     xyz[0]=0.333333333333333333333*(xyz1[0]+xyz2[0]+xyz3[0]);
827     xyz[1]=0.333333333333333333333*(xyz1[1]+xyz2[1]+xyz3[1]);
828     xyz[2]=0.333333333333333333333*(xyz1[2]+xyz2[2]+xyz3[2]);
829     OT_VECTEUR_3D vec1(xyz1,xyz);
830     OT_VECTEUR_3D vec2(xyz2,xyz);
831     OT_VECTEUR_3D vec3(xyz3,xyz);
832     rayonrecherche=std::max(vec1.get_longueur(),vec2.get_longueur());
833     rayonrecherche=std::max(rayonrecherche,vec3.get_longueur())*1.1;
834     //double rayonrecherche=sqrt((u-uv1[0])*(u-uv1[0])+(v-uv1[1])*(v-uv1[1]))*1.1;
835     ntree_de_front->rechercher(xyz[0],xyz[1],xyz[2],rayonrecherche,liste_trouvee);
836 5 int nb_entite=liste_trouvee.get_nb();
837     double du=decalage->calcul_decalage_parametre_u(u);
838     double dv=decalage->calcul_decalage_parametre_v(v);
839     double u1=decalage->decalage_parametre_u(uv1[0],du);
840     double v1=decalage->decalage_parametre_v(uv1[1],dv);
841     double u2=decalage->decalage_parametre_u(uv2[0],du);
842     double v2=decalage->decalage_parametre_v(uv2[1],dv);
843     double u3=decalage->decalage_parametre_u(uv3[0],du);
844     double v3=decalage->decalage_parametre_v(uv3[1],dv);
845     MG_FRONT_2D *ref1=NULL,*ref2=NULL,*ref3=NULL;
846     double angle1,angle2,angle3;
847     for(int j=0;j<nb_entite;j++)
848     {
849     MG_FRONT_2D *ft=liste_trouvee.get(j);
850     MG_NOEUD* noeudfront=ft->get_noeud2();
851     double uvfront[2];
852     uvfront[0]=noeudfront->get_u();
853     uvfront[1]=noeudfront->get_v();
854     double ufront=decalage->decalage_parametre_u(uvfront[0],du);
855     double vfront=decalage->decalage_parametre_v(uvfront[1],dv);
856     //noeud 1 du triangle a inserer
857 francois 632 //double dist=sqrt((ufront-u1)*(ufront-u1)+(vfront-v1)*(vfront-v1));
858     double *xyzfront=noeudfront->get_coord();
859     //mgface->evaluer(uvfront,xyzfront);
860     OT_VECTEUR_3D vec(xyz1,xyzfront);
861     double dist=vec.get_longueur();
862 francois 54 if (dist<1e-5*rayonrecherche)
863 5 {
864     if (ref1==NULL)
865     {
866     ref1=ft;
867     OT_VECTEUR_3D vecteur_baseu(u2-u1,v2-v1,0.);
868     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
869     double uvfrontautre[2];
870     uvfrontautre[0]=noeudfrontautre->get_u();
871     uvfrontautre[1]=noeudfrontautre->get_v();
872     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
873     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
874     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
875     vecteur_baseu.norme();
876     vecteur_front.norme();
877     double cosangle=vecteur_baseu*vecteur_front;
878     if (cosangle>1.) cosangle=1.;
879     if (cosangle<-1.) cosangle=(-1.);
880     angle1=acos(cosangle);
881     OT_VECTEUR_3D w(0.,0.,1.);
882     OT_VECTEUR_3D vecteur_basev=w&vecteur_baseu;
883     double sangle=vecteur_basev*vecteur_front;
884     if (sangle<0.) angle1=2*M_PI-angle1;
885     }
886     else
887     {
888     OT_VECTEUR_3D vecteur_baseu(u2-u1,v2-v1,0.);
889     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
890     double uvfrontautre[2];
891     uvfrontautre[0]=noeudfrontautre->get_u();
892     uvfrontautre[1]=noeudfrontautre->get_v();
893     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
894     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
895     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
896     vecteur_baseu.norme();
897     vecteur_front.norme();
898     double cosangle=vecteur_baseu*vecteur_front;
899     if (cosangle>1.) cosangle=1.;
900     if (cosangle<-1.) cosangle=(-1.);
901     double angle=acos(cosangle);
902     OT_VECTEUR_3D w(0.,0.,1.);
903     OT_VECTEUR_3D vecteur_basev=w&vecteur_baseu;
904     double sangle=vecteur_basev*vecteur_front;
905     if (sangle<0.) angle=2*M_PI-angle;
906     if (angle<angle1) {angle1=angle;ref1=ft;}
907     }
908     }
909     //noeud 2 du triangle a inserer
910 francois 632 //mgface->evaluer(uvfront,xyzfront);
911     OT_VECTEUR_3D vec2(xyz2,xyzfront);
912     dist=vec2.get_longueur();
913     //dist=sqrt((ufront-u2)*(ufront-u2)+(vfront-v2)*(vfront-v2));
914 francois 54 if (dist<1e-5*rayonrecherche)
915 5 {
916     if (ref2==NULL)
917     {
918     ref2=ft;
919     OT_VECTEUR_3D vecteur_baseu(u3-u2,v3-v2,0.);
920     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
921     double uvfrontautre[2];
922     uvfrontautre[0]=noeudfrontautre->get_u();
923     uvfrontautre[1]=noeudfrontautre->get_v();
924     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
925     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
926     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
927     vecteur_baseu.norme();
928     vecteur_front.norme();
929     double cosangle=vecteur_baseu*vecteur_front;
930     if (cosangle>1.) cosangle=1.;
931     if (cosangle<-1.) cosangle=(-1.);
932     angle2=acos(cosangle);
933     OT_VECTEUR_3D w(0.,0.,1.);
934     OT_VECTEUR_3D vecteur_basev=w&vecteur_baseu;
935     double sangle=vecteur_basev*vecteur_front;
936     if (sangle<0.) angle2=2*M_PI-angle2;
937     }
938     else
939     {
940     OT_VECTEUR_3D vecteur_baseu(u3-u2,v3-v2,0.);
941     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
942     double uvfrontautre[2];
943     uvfrontautre[0]=noeudfrontautre->get_u();
944     uvfrontautre[1]=noeudfrontautre->get_v();
945     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
946     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
947     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
948     vecteur_baseu.norme();
949     vecteur_front.norme();
950     double cosangle=vecteur_baseu*vecteur_front;
951     if (cosangle>1.) cosangle=1.;
952     if (cosangle<-1.) cosangle=(-1.);
953     double angle=acos(cosangle);
954     OT_VECTEUR_3D w(0.,0.,1.);
955     OT_VECTEUR_3D vecteur_basev=w&vecteur_baseu;
956     double sangle=vecteur_basev*vecteur_front;
957     if (sangle<0.) angle=2*M_PI-angle;
958     if (angle<angle2) {angle2=angle;ref2=ft;}
959     }
960     }
961     //noeud 3 du triangle a inserer
962 francois 632 //dist=sqrt((ufront-u3)*(ufront-u3)+(vfront-v3)*(vfront-v3));
963     //mgface->evaluer(uvfront,xyzfront);
964     OT_VECTEUR_3D vec3(xyz3,xyzfront);
965     dist=vec3.get_longueur();
966 francois 54 if (dist<1e-5*rayonrecherche)
967 5 {
968     if (ref3==NULL)
969     {
970     ref3=ft;
971     OT_VECTEUR_3D vecteur_baseu(u1-u3,v1-v3,0.);
972     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
973     double uvfrontautre[2];
974     uvfrontautre[0]=noeudfrontautre->get_u();
975     uvfrontautre[1]=noeudfrontautre->get_v();
976     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
977     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
978     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
979     vecteur_baseu.norme();
980     vecteur_front.norme();
981     double cosangle=vecteur_baseu*vecteur_front;
982     if (cosangle>1.) cosangle=1.;
983     if (cosangle<-1.) cosangle=(-1.);
984     angle3=acos(cosangle);
985     OT_VECTEUR_3D w(0.,0.,1.);
986     OT_VECTEUR_3D vecteur_basev=w&vecteur_baseu;
987     double sangle=vecteur_basev*vecteur_front;
988     if (sangle<0.) angle3=2*M_PI-angle3;
989     }
990     else
991     {
992     OT_VECTEUR_3D vecteur_baseu(u1-u3,v1-v3,0.);
993     MG_NOEUD* noeudfrontautre=ft->get_noeud1();
994     double uvfrontautre[2];
995     uvfrontautre[0]=noeudfrontautre->get_u();
996     uvfrontautre[1]=noeudfrontautre->get_v();
997     double ufrontautre=decalage->decalage_parametre_u(uvfrontautre[0],du);
998     double vfrontautre=decalage->decalage_parametre_v(uvfrontautre[1],dv);
999     OT_VECTEUR_3D vecteur_front(ufrontautre-ufront,vfrontautre-vfront,0.);
1000     vecteur_baseu.norme();
1001     vecteur_front.norme();
1002     double cosangle=vecteur_baseu*vecteur_front;
1003     if (cosangle>1.) cosangle=1.;
1004     if (cosangle<-1.) cosangle=(-1.);
1005     double angle=acos(cosangle);
1006     OT_VECTEUR_3D w(0.,0.,1.);
1007     OT_VECTEUR_3D vecteur_basev=w&vecteur_baseu;
1008     double sangle=vecteur_basev*vecteur_front;
1009     if (sangle<0.) angle=2*M_PI-angle;
1010     if (angle<angle3) {angle3=angle;ref3=ft;}
1011    
1012     }
1013     }
1014     }
1015     // creation du triangle
1016     MG_NOEUD* noeud1;
1017     MG_NOEUD* noeud2;
1018     MG_NOEUD* noeud3;
1019     if (ref1!=NULL) noeud1=ref1->get_noeud2();
1020     else
1021     {
1022     double *coord=noeudori1->get_coord();
1023 francois 35 noeud1=mg_maillage->ajouter_mg_noeud(mgface,coord[0],coord[1],coord[2],IMPOSE);
1024 5 noeud1->change_u(uv1[0]);
1025     noeud1->change_v(uv1[1]);
1026 francois 447 noeud1->change_nouveau_numero(noeudori1->get_nouveau_numero());
1027 5 }
1028     if (ref2!=NULL) noeud2=ref2->get_noeud2();
1029     else
1030     {
1031     double *coord=noeudori2->get_coord();
1032 francois 35 noeud2=mg_maillage->ajouter_mg_noeud(mgface,coord[0],coord[1],coord[2],IMPOSE);
1033 5 noeud2->change_u(uv2[0]);
1034     noeud2->change_v(uv2[1]);
1035 francois 155 noeud2->change_nouveau_numero(noeudori2->get_nouveau_numero());
1036 5 }
1037     if (ref3!=NULL) noeud3=ref3->get_noeud2();
1038     else
1039     {
1040     double *coord=noeudori3->get_coord();
1041 francois 35 noeud3=mg_maillage->ajouter_mg_noeud(mgface,coord[0],coord[1],coord[2],IMPOSE);
1042 5 noeud3->change_u(uv3[0]);
1043     noeud3->change_v(uv3[1]);
1044 francois 447 noeud3->change_nouveau_numero(noeudori3->get_nouveau_numero());
1045 5 }
1046     MG_SEGMENT* mgsegment1=mg_maillage->get_mg_segment(noeud1->get_id(),noeud2->get_id());
1047     MG_SEGMENT* mgsegment2=mg_maillage->get_mg_segment(noeud2->get_id(),noeud3->get_id());
1048     MG_SEGMENT* mgsegment3=mg_maillage->get_mg_segment(noeud3->get_id(),noeud1->get_id());
1049     int seg1=1,seg2=1,seg3=1;
1050 francois 35 if (mgsegment1==NULL) mgsegment1=mg_maillage->ajouter_mg_segment(mgface,noeud1,noeud2,IMPOSE); else seg1=0;
1051     if (mgsegment2==NULL) mgsegment2=mg_maillage->ajouter_mg_segment(mgface,noeud2,noeud3,IMPOSE); else seg2=0;
1052     if (mgsegment3==NULL) mgsegment3=mg_maillage->ajouter_mg_segment(mgface,noeud3,noeud1,IMPOSE); else seg3=0;
1053 francois 54 M3D_TRIANGLE* mtriangle=new M3D_TRIANGLE(mgface,noeud1,noeud2,noeud3,mgsegment1,mgsegment2,mgsegment3,IMPOSE);
1054 5 mg_maillage->ajouter_mg_triangle(mtriangle);
1055 francois 632 int segtotal=seg1+seg2+seg3;
1056 5 if (segtotal==3)
1057     {
1058     MG_FRONT_2D* ft1=ajouter_front(front_courant,noeud1,noeud2,mgsegment1);
1059     MG_FRONT_2D* ft2=ajouter_front(front_courant,noeud2,noeud3,mgsegment2);
1060     MG_FRONT_2D* ft3=ajouter_front(front_courant,noeud3,noeud1,mgsegment3);
1061     if (ref1==NULL)
1062     {
1063     ft1->changer_front_precedent(ft3);
1064     ft3->changer_front_suivant(ft1);
1065     }
1066     else
1067     {
1068     MG_FRONT_2D* tmp=ref1->get_front_suivant();
1069     ft1->changer_front_precedent(ref1);
1070     ref1->changer_front_suivant(ft1);
1071     ft3->changer_front_suivant(tmp);
1072     tmp->changer_front_precedent(ft3);
1073     }
1074     if (ref2==NULL)
1075     {
1076     ft1->changer_front_suivant(ft2);
1077     ft2->changer_front_precedent(ft1);
1078     }
1079     else
1080     {
1081     MG_FRONT_2D* tmp=ref2->get_front_suivant();
1082     ft2->changer_front_precedent(ref2);
1083     ref2->changer_front_suivant(ft2);
1084     ft1->changer_front_suivant(tmp);
1085     tmp->changer_front_precedent(ft1);
1086     }
1087     if (ref3==NULL)
1088     {
1089     ft2->changer_front_suivant(ft3);
1090     ft3->changer_front_precedent(ft2);
1091     }
1092     else
1093     {
1094     MG_FRONT_2D* tmp=ref3->get_front_suivant();
1095     ft3->changer_front_precedent(ref3);
1096     ref3->changer_front_suivant(ft3);
1097     ft2->changer_front_suivant(tmp);
1098     tmp->changer_front_precedent(ft2);
1099     }
1100     }
1101     if (segtotal==2)
1102     {
1103     MG_FRONT_2D *front1,*front2,*front3;
1104     MG_SEGMENT *segcre1,*segcre2;
1105     if (seg3==0)
1106     {
1107     front1=ref1;front2=ref3;front3=ref2;
1108     segcre1=mgsegment1;segcre2=mgsegment2;
1109     }
1110     if (seg1==0)
1111     {
1112     front1=ref2;front2=ref1;front3=ref3;
1113     segcre1=mgsegment2;segcre2=mgsegment3;
1114     }
1115     if (seg2==0)
1116     {
1117     front1=ref3;front2=ref2;front3=ref1;
1118     segcre1=mgsegment3;segcre2=mgsegment1;
1119     }
1120     MG_FRONT_2D* ft1=ajouter_front(front_courant,segcre1->get_noeud1(),segcre1->get_noeud2(),segcre1);
1121     MG_FRONT_2D* ft2=ajouter_front(front_courant,segcre2->get_noeud1(),segcre2->get_noeud2(),segcre2);
1122     front1->changer_front_suivant(ft1);
1123     ft1->changer_front_precedent(front1);
1124     if (front3==NULL)
1125     {
1126     ft1->changer_front_suivant(ft2);
1127     ft2->changer_front_precedent(ft1);
1128     }
1129     else
1130     {
1131     MG_FRONT_2D* tmp=front3->get_front_suivant();
1132     ft1->changer_front_suivant(tmp);
1133     tmp->changer_front_precedent(ft1);
1134     ft2->changer_front_precedent(front3);
1135     front3->changer_front_suivant(ft2);
1136     }
1137     ft2->changer_front_suivant(front2->get_front_suivant());
1138     front2->get_front_suivant()->changer_front_precedent(ft2);
1139     supprimer_front(front2);
1140     }
1141     if (segtotal==1)
1142     {
1143     MG_FRONT_2D *front1,*front2,*front3;
1144     MG_SEGMENT *segcre1;
1145     if (seg1==1)
1146     {
1147     front1=ref1;front2=ref3;front3=ref2;
1148     segcre1=mgsegment1;
1149     }
1150     if (seg2==1)
1151     {
1152     front1=ref2;front2=ref1;front3=ref3;
1153     segcre1=mgsegment2;
1154     }
1155     if (seg3==1)
1156     {
1157     front1=ref3;front2=ref2;front3=ref1;
1158     segcre1=mgsegment3;
1159     }
1160     MG_FRONT_2D* ft1=ajouter_front(front_courant,segcre1->get_noeud1(),segcre1->get_noeud2(),segcre1);
1161     front1->changer_front_suivant(ft1);
1162     ft1->changer_front_precedent(front1);
1163     ft1->changer_front_suivant(front3->get_front_suivant());
1164     front3->get_front_suivant()->changer_front_precedent(ft1);
1165     supprimer_front(front2);
1166     supprimer_front(front3);
1167     }
1168     if (segtotal==0)
1169     {
1170     supprimer_front(ref1);
1171     supprimer_front(ref2);
1172     supprimer_front(ref3);
1173     }
1174    
1175    
1176     }
1177     }
1178    
1179 francois 447 void MAILLEUR2D::change_niveau_optimisation(int num)
1180     {
1181     niveau_optimisation=num;
1182     }
1183 5
1184 francois 447 int MAILLEUR2D::get_niveau_optimisation(void)
1185     {
1186     return niveau_optimisation;
1187     }
1188    
1189    
1190 5 // FONCTIONS GERANT L ENTITE FRONT (ajout suppression et tout le bordel)
1191     MG_FRONT_2D* MAILLEUR2D::ajouter_front(FRONT& front,MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_SEGMENT* segment)
1192     {
1193     MG_FRONT_2D* ft=new MG_FRONT_2D(noeud1,noeud2,segment);
1194     ajouter_front(front,ft);
1195 francois 632 ntree_de_front->inserer(ft);
1196 5 return ft;
1197     }
1198    
1199     void MAILLEUR2D::ajouter_front(FRONT& front,MG_FRONT_2D *ft)
1200     {
1201     std::pair<const double,MG_FRONT_2D*> tmp(ft->get_segment()->get_longueur(),ft);
1202     front.insert(tmp);
1203     }
1204    
1205     MG_FRONT_2D* MAILLEUR2D::get_front(FRONT& front,unsigned int num)
1206     {
1207     FRONT::iterator i=front.begin();
1208     for (unsigned long j=0;j<num;j++) i++;
1209     return ((*i).second);
1210     }
1211    
1212    
1213     unsigned int MAILLEUR2D::get_nb_front(FRONT& front)
1214     {
1215     return front.size();
1216     }
1217    
1218    
1219     void MAILLEUR2D::supprimer_front(MG_FRONT_2D* ft)
1220     {
1221 francois 632 ntree_de_front->supprimer(ft);
1222 5
1223     FRONT::iterator j=front_courant.lower_bound(ft->get_segment()->get_longueur());
1224     int ok=0;
1225     if (j==front_courant.end()) ok=2;
1226     while (ok==0)
1227     {
1228     MG_FRONT_2D* fttmp=(*j).second;
1229     if (ft==fttmp) {ok=1;front_courant.erase(j);}
1230     if (fttmp->get_segment()->get_longueur()>ft->get_segment()->get_longueur()) ok=2;
1231     j++;
1232     }
1233     if (ok!=1)
1234     {
1235     j=front_attente.lower_bound(ft->get_segment()->get_longueur());
1236     while (ft!=(*j).second) j++;
1237     front_attente.erase(j);
1238     }
1239    
1240     delete ft;
1241     }
1242    
1243     void MAILLEUR2D::echange_de_front(FRONT& front_original,FRONT& front_destination,MG_FRONT_2D* ft)
1244     {
1245     FRONT::iterator j=front_original.lower_bound(ft->get_segment()->get_longueur());
1246     while (ft!=(*j).second) j++;
1247     front_original.erase(j);
1248     ajouter_front(front_destination,ft);
1249     }
1250 francois 632
1251