ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/mailleur2d.cpp
Revision: 210
Committed: Fri Jul 31 16:19:04 2009 UTC (15 years, 9 months ago) by francois
Original Path: magic/lib/mailleur/mailleur/src/mailleur2d.cpp
File size: 65200 byte(s)
Log Message: 
Nouvelle facon de representer la deforme + groupe topologique pour le maillage + bug idmax

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