ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/mailleur_fem.cpp
Revision: 883
Committed: Thu Apr 20 13:38:18 2017 UTC (8 years ago) by francois
File size: 36055 byte(s)
Log Message:
Creation d'un mailleur FEM pour pouvoir avoir des stratégies paramétrées de maillage. L'ancienne méthode construit disparait et est remplacée par l'utilisation d'un MAILLEUR_FEM.
Stratégie de quadratisation mise en place : déplacer les noeuds pour s'assurer que les tetras quadratiques présentent une distortion au dessu d'une borne inférieure. Mais ces noeuds quittent la géométrie.

Les paramètres dans ~/.magic.

File Contents

# User Rev Content
1 francois 883 //------------------------------------------------------------
2     //------------------------------------------------------------
3     // MAGiC
4     // Jean Christophe Cuilli�re et Vincent FRANCOIS
5     // D�partement de G�nie M�canique - UQTR
6     //------------------------------------------------------------
7     // Le projet MAGIC est un projet de recherche du d�partement
8     // de g�nie m�canique de l'Universit� du Qu�bec �
9     // Trois Rivi�res
10     // Les librairies ne peuvent �tre utilis�es sans l'accord
11     // des auteurs (contact : francois@uqtr.ca)
12     //------------------------------------------------------------
13     //------------------------------------------------------------
14     //
15     // mailleur.cpp
16     //
17     //------------------------------------------------------------
18     //------------------------------------------------------------
19     // COPYRIGHT 2000
20     // Version du 02/03/2006 � 11H23
21     //------------------------------------------------------------
22     //------------------------------------------------------------
23    
24    
25     #include "gestionversion.h"
26     #include "mailleur_fem.h"
27     #include "mg_gestionnaire.h"
28     #include "fem_segment2.h"
29     #include "fem_mini_segment2.h"
30     #include "fem_segment3.h"
31     #include "fem_triangle3.h"
32     #include "fem_triangle6.h"
33     #include "fem_quadrangle4.h"
34     #include "fem_quadrangle8.h"
35     #include "fem_tetra4.h"
36     #include "fem_tetra10.h"
37     #include "fem_hexa8.h"
38     #include "fem_hexa20.h"
39     #include "fem_penta6.h"
40     #include "fem_penta15.h"
41     #include "fem_graph_noeud.h"
42     #include "fem_maillage_quadratique_outils.h"
43     #include "ot_decalage_parametre.h"
44    
45    
46    
47     MAILLEUR_FEM::MAILLEUR_FEM(OT_CPU* comp):MAILLEUR(false,comp)
48     {
49    
50     }
51    
52     MAILLEUR_FEM::MAILLEUR_FEM(MAILLEUR_FEM &mdd):MAILLEUR(mdd)
53     {
54     ini_param();
55     }
56    
57    
58     MAILLEUR_FEM::~MAILLEUR_FEM()
59     {
60     }
61    
62    
63     void MAILLEUR_FEM::maille(FEM_MAILLAGE* fem,int num)
64     {
65     int degre=fem->get_degre();
66     MG_MAILLAGE* maillage=fem->get_mg_maillage();
67     MG_GEOMETRIE* mg_geometrie=maillage->get_mg_geometrie();
68     TPL_MAP_ENTITE<MG_SEGMENT*> mini;
69     if (mg_geometrie!=NULL)
70     if (mg_geometrie->get_nb_mg_poutre()!=0)
71     if (mg_geometrie->get_nb_mg_volume()!=0)
72     recherche_connexion_multidimension(mini,maillage,mg_geometrie);
73     if (degre==1) construire_lineaire(fem,mini,maillage,mg_geometrie);
74     if (degre==2) construire_quadratique(fem,mini,maillage,mg_geometrie);
75     if (num>0) optimise_numerotation(fem);
76     if (num>1) copie_numerotation_opt(fem);
77     BOITE_3D b;
78     int nx,ny,nz;
79     bool structure=maillage->get_param_structure(b,nx,ny,nz);
80     if (structure==true)
81     fem->change_param_structure(b,nx,ny,nz);
82     }
83    
84    
85     void MAILLEUR_FEM::recherche_connexion_multidimension(TPL_MAP_ENTITE<MG_SEGMENT*> &mini,MG_MAILLAGE* maillage,MG_GEOMETRIE* mg_geometrie)
86     {
87     LISTE_MG_NOEUD::iterator it;
88     for (MG_NOEUD* no=maillage->get_premier_noeud(it);no!=NULL;no=maillage->get_suivant_noeud(it))
89     {
90     if ((no->get_lien_hexa()->get_nb()>0) || (no->get_lien_tetra()->get_nb()))
91     if (no->get_lien_segment()->get_nb()>0)
92     {
93     bool estjoint=false;
94     for (int i=0;i<no->get_lien_segment()->get_nb();i++)
95     {
96     LISTE_MG_POUTRE::iterator itp;
97     for (MG_POUTRE* pou=mg_geometrie->get_premier_poutre(itp);pou!=NULL;pou=mg_geometrie->get_suivant_poutre(itp))
98     if (pou->contient_element(no->get_lien_segment()->get(i))) {estjoint=true;}
99     }
100     if (estjoint)
101     {
102     for (int i=0;i<no->get_lien_segment()->get_nb();i++)
103     {
104     bool bord=false;
105     LISTE_MG_POUTRE::iterator itp;
106     int nbele2=no->get_lien_segment()->get(i)->get_lien_triangle()->get_nb();
107     for (int j=0;j<nbele2;j++)
108     if (no->get_lien_segment()->get(i)->get_lien_triangle()->get(j)->get_lien_tetra()->get_nb()>0)
109     if (no->get_lien_segment()->get(i)->get_lien_triangle()->get(j)->get_lien_tetra()->get_nb()<2) bord=true;
110     nbele2=no->get_lien_segment()->get(i)->get_lien_quadrangle()->get_nb();
111     for (int j=0;j<nbele2;j++)
112     if (no->get_lien_segment()->get(i)->get_lien_quadrangle()->get(j)->get_lien_hexa()->get_nb()>0)
113     if (no->get_lien_segment()->get(i)->get_lien_quadrangle()->get(j)->get_lien_hexa()->get_nb()<2) bord=true;
114     if (bord==true) mini.ajouter(no->get_lien_segment()->get(i));
115     }
116     //std::cout << no->get_id() << " est un joint 1D-3D : " << mini.get_nb() << " mini " << std::endl;
117    
118     }
119     }
120    
121     }
122     LISTE_MG_TRIANGLE::iterator itt;
123     for (MG_TRIANGLE* tri=maillage->get_premier_triangle(itt);tri!=NULL;tri=maillage->get_suivant_triangle(itt))
124     if (tri->get_origine()==MAGIC::ORIGINE::SECTION)
125     {
126     mini.ajouter(tri->get_segment1());
127     mini.ajouter(tri->get_segment2());
128     mini.ajouter(tri->get_segment3());
129     }
130    
131     LISTE_MG_QUADRANGLE::iterator itq;
132     for (MG_QUADRANGLE* qua=maillage->get_premier_quadrangle(itq);qua!=NULL;qua=maillage->get_suivant_quadrangle(itq))
133     if (qua->get_origine()==MAGIC::ORIGINE::SECTION)
134     {
135     mini.ajouter(qua->get_segment1());
136     mini.ajouter(qua->get_segment2());
137     mini.ajouter(qua->get_segment3());
138     mini.ajouter(qua->get_segment4());
139     }
140    
141    
142     }
143    
144    
145    
146    
147     void MAILLEUR_FEM::construire_lineaire(FEM_MAILLAGE* fem,TPL_MAP_ENTITE<MG_SEGMENT*> &mini,MG_MAILLAGE* maillage,MG_GEOMETRIE* mg_geometrie)
148     {
149     LISTE_MG_NOEUD::iterator it_noeud;
150     MG_NOEUD* mgnoeud=maillage->get_premier_noeud(it_noeud);
151     std::vector<FEM_NOEUD*> lstnoeuddirect;
152     unsigned int i=0;
153     while (mgnoeud)
154     {
155     mgnoeud->change_nouveau_numero(i);
156     FEM_NOEUD* femnoeud=new FEM_NOEUD(mgnoeud);
157     femnoeud->change_numero(i+1);
158     i++;
159     fem->ajouter_fem_noeud(femnoeud);
160     lstnoeuddirect.insert(lstnoeuddirect.end(),femnoeud);
161     if (mgnoeud->get_lien_topologie()!=NULL)
162     if (mgnoeud->get_lien_topologie()->get_dimension()==0)
163     {
164     FEM_ELEMENT0* femele=new FEM_ELEMENT0(mgnoeud->get_lien_topologie(),mgnoeud,&femnoeud);
165     fem-> ajouter_fem_element0(femele);
166     }
167     mgnoeud=maillage->get_suivant_noeud(it_noeud);
168    
169     }
170     int dimsansgeo=0;
171     if (mg_geometrie==NULL)
172     {
173     if (maillage->get_nb_mg_tetra()+maillage->get_nb_mg_hexa()!=0) dimsansgeo=3;
174     else if (maillage->get_nb_mg_triangle()+maillage->get_nb_mg_quadrangle()!=0) dimsansgeo=2;
175     else if (maillage->get_nb_mg_segment()!=0) dimsansgeo=1;
176     }
177     if (mg_geometrie!=NULL)
178     if (strcmp(mg_geometrie->get_type_geometrie(),"VIRTUEL")==0)
179     {
180     if (maillage->get_nb_mg_tetra()+maillage->get_nb_mg_hexa()!=0) dimsansgeo=3;
181     else if (maillage->get_nb_mg_triangle()+maillage->get_nb_mg_quadrangle()!=0) dimsansgeo=2;
182     else if (maillage->get_nb_mg_segment()!=0) dimsansgeo=1;
183     }
184    
185     if ((mg_geometrie!=NULL) || ( (mg_geometrie==NULL) && (dimsansgeo==1)))
186     {
187     LISTE_MG_SEGMENT::iterator it_seg;
188     MG_SEGMENT * mgseg = maillage->get_premier_segment(it_seg);
189     while (mgseg)
190     {
191     if (mgseg->get_lien_topologie()!=NULL)
192     if (mgseg->get_lien_topologie()->get_dimension()!=1)
193     {
194     mgseg=maillage->get_suivant_segment(it_seg);
195     continue;
196     }
197     if (mg_geometrie!=NULL)
198     if (mgseg->get_lien_topologie()==NULL)
199     if (dimsansgeo!=1)
200     {
201     mgseg=maillage->get_suivant_segment(it_seg);
202     continue;
203     }
204     FEM_NOEUD *tabnoeud[2];
205     tabnoeud[0]=lstnoeuddirect[mgseg->get_noeud1()->get_nouveau_numero()];
206     tabnoeud[1]=lstnoeuddirect[mgseg->get_noeud2()->get_nouveau_numero()];
207     FEM_SEGMENT2* seg=new FEM_SEGMENT2(mgseg,tabnoeud);
208     fem->ajouter_fem_element1(seg);
209     mgseg=maillage->get_suivant_segment(it_seg);
210     }
211     TPL_MAP_ENTITE<MG_SEGMENT*>::ITERATEUR itm;
212     for (MG_SEGMENT* mgseg=mini.get_premier(itm);mgseg!=NULL;mgseg=mini.get_suivant(itm))
213     {
214     FEM_NOEUD *tabnoeud[2];
215     tabnoeud[0]=lstnoeuddirect[mgseg->get_noeud1()->get_nouveau_numero()];
216     tabnoeud[1]=lstnoeuddirect[mgseg->get_noeud2()->get_nouveau_numero()];
217     FEM_SEGMENT2* seg=new FEM_MINI_SEGMENT2(mgseg,tabnoeud);
218     fem->ajouter_fem_element1(seg);
219     }
220     }
221     if ((mg_geometrie!=NULL) || ( (mg_geometrie==NULL) && (dimsansgeo==2)))
222     {
223     LISTE_MG_TRIANGLE::iterator it_tri;
224     MG_TRIANGLE* mgtri=maillage->get_premier_triangle(it_tri);
225     while (mgtri)
226     {
227     if (mgtri->get_lien_topologie()!=NULL) if (mgtri->get_lien_topologie()->get_dimension()!=2)
228     {
229     mgtri=maillage->get_suivant_triangle(it_tri);
230     continue;
231     }
232     if (mg_geometrie!=NULL)
233     if (mgtri->get_lien_topologie()==NULL)
234     if (dimsansgeo!=2)
235     {
236     mgtri=maillage->get_suivant_triangle(it_tri);
237     continue;
238     }
239     FEM_NOEUD *tabnoeud[3];
240     tabnoeud[0]=lstnoeuddirect[mgtri->get_noeud1()->get_nouveau_numero()];
241     tabnoeud[1]=lstnoeuddirect[mgtri->get_noeud2()->get_nouveau_numero()];
242     tabnoeud[2]=lstnoeuddirect[mgtri->get_noeud3()->get_nouveau_numero()];
243     FEM_TRIANGLE3* tri=new FEM_TRIANGLE3(mgtri,tabnoeud);
244     fem->ajouter_fem_element2(tri);
245     mgtri=maillage->get_suivant_triangle(it_tri);
246     }
247     LISTE_MG_QUADRANGLE::iterator it_quad;
248     MG_QUADRANGLE* mgquad=maillage->get_premier_quadrangle(it_quad);
249     while (mgquad)
250     {
251     if (mgquad->get_lien_topologie()!=NULL) if (mgquad->get_lien_topologie()->get_dimension()!=2) {
252     mgquad=maillage->get_suivant_quadrangle(it_quad);
253     continue;
254     }
255     if (mg_geometrie!=NULL)
256     if (mgquad->get_lien_topologie()==NULL)
257     if (dimsansgeo!=2)
258     {
259     mgquad=maillage->get_suivant_quadrangle(it_quad);
260     continue;
261     }
262     FEM_NOEUD *tabnoeud[4];
263     tabnoeud[0]=lstnoeuddirect[mgquad->get_noeud1()->get_nouveau_numero()];
264     tabnoeud[1]=lstnoeuddirect[mgquad->get_noeud2()->get_nouveau_numero()];
265     tabnoeud[2]=lstnoeuddirect[mgquad->get_noeud3()->get_nouveau_numero()];
266     tabnoeud[3]=lstnoeuddirect[mgquad->get_noeud4()->get_nouveau_numero()];
267     FEM_QUADRANGLE4* quad=new FEM_QUADRANGLE4(mgquad,tabnoeud);
268     fem->ajouter_fem_element2(quad);
269     mgquad=maillage->get_suivant_quadrangle(it_quad);
270     }
271     }
272     if ((mg_geometrie!=NULL) || ( (mg_geometrie==NULL) && (dimsansgeo==3)))
273     {
274     LISTE_MG_TETRA::iterator it_tetra;
275     MG_TETRA* mgtetra=maillage->get_premier_tetra(it_tetra);
276     while (mgtetra)
277     {
278     FEM_NOEUD *tabnoeud[4];
279     tabnoeud[0]=lstnoeuddirect[mgtetra->get_noeud1()->get_nouveau_numero()];
280     tabnoeud[1]=lstnoeuddirect[mgtetra->get_noeud2()->get_nouveau_numero()];
281     tabnoeud[2]=lstnoeuddirect[mgtetra->get_noeud3()->get_nouveau_numero()];
282     tabnoeud[3]=lstnoeuddirect[mgtetra->get_noeud4()->get_nouveau_numero()];
283     FEM_TETRA4* tet=new FEM_TETRA4(mgtetra,tabnoeud);
284     fem->ajouter_fem_element3(tet);
285     mgtetra=maillage->get_suivant_tetra(it_tetra);
286     }
287     LISTE_MG_HEXA::iterator it_hex;
288     MG_HEXA* mghex=maillage->get_premier_hexa(it_hex);
289     while (mghex)
290     {
291     FEM_NOEUD *tabnoeud[8];
292     tabnoeud[0]=lstnoeuddirect[mghex->get_noeud1()->get_nouveau_numero()];
293     tabnoeud[1]=lstnoeuddirect[mghex->get_noeud2()->get_nouveau_numero()];
294     tabnoeud[2]=lstnoeuddirect[mghex->get_noeud3()->get_nouveau_numero()];
295     tabnoeud[3]=lstnoeuddirect[mghex->get_noeud4()->get_nouveau_numero()];
296     tabnoeud[4]=lstnoeuddirect[mghex->get_noeud5()->get_nouveau_numero()];
297     tabnoeud[5]=lstnoeuddirect[mghex->get_noeud6()->get_nouveau_numero()];
298     tabnoeud[6]=lstnoeuddirect[mghex->get_noeud7()->get_nouveau_numero()];
299     tabnoeud[7]=lstnoeuddirect[mghex->get_noeud8()->get_nouveau_numero()];
300     FEM_HEXA8* hex=new FEM_HEXA8(mghex,tabnoeud);
301     fem->ajouter_fem_element3(hex);
302     mghex=maillage->get_suivant_hexa(it_hex);
303     }
304     LISTE_MG_PENTA::iterator it_pen;
305     MG_PENTA* mgpen=maillage->get_premier_penta(it_pen);
306     while (mgpen)
307     {
308     FEM_NOEUD *tabnoeud[8];
309     tabnoeud[0]=lstnoeuddirect[mgpen->get_noeud1()->get_nouveau_numero()];
310     tabnoeud[1]=lstnoeuddirect[mgpen->get_noeud2()->get_nouveau_numero()];
311     tabnoeud[2]=lstnoeuddirect[mgpen->get_noeud3()->get_nouveau_numero()];
312     tabnoeud[3]=lstnoeuddirect[mgpen->get_noeud4()->get_nouveau_numero()];
313     tabnoeud[4]=lstnoeuddirect[mgpen->get_noeud5()->get_nouveau_numero()];
314     tabnoeud[5]=lstnoeuddirect[mgpen->get_noeud6()->get_nouveau_numero()];
315     FEM_PENTA6* pen=new FEM_PENTA6(mgpen,tabnoeud);
316     fem->ajouter_fem_element3(pen);
317     mgpen=maillage->get_suivant_penta(it_pen);
318     }
319     }
320     }
321    
322    
323    
324    
325     void MAILLEUR_FEM::construire_quadratique(FEM_MAILLAGE* fem,TPL_MAP_ENTITE<MG_SEGMENT*> &mini,MG_MAILLAGE* maillage,MG_GEOMETRIE* mg_geometrie)
326     {
327     LISTE_MG_NOEUD::iterator it_nd;
328     MG_NOEUD* mgnoeud=maillage->get_premier_noeud(it_nd);
329     std::vector<FEM_NOEUD*> lstnoeuddirect;
330     unsigned int i=0;
331     while (mgnoeud)
332     {
333     mgnoeud->change_nouveau_numero(i);
334     FEM_NOEUD* femnoeud=new FEM_NOEUD(mgnoeud);
335     femnoeud->change_numero(i+1);
336     i++ ;
337     fem->ajouter_fem_noeud(femnoeud);
338     //mgnoeud=maillage->get_suivant_noeud(it_nd);
339     lstnoeuddirect.insert(lstnoeuddirect.end(),femnoeud);
340     if (mgnoeud->get_lien_topologie()!=NULL)
341     if (mgnoeud->get_lien_topologie()->get_dimension()==0)
342     {
343     FEM_ELEMENT0* femele=new FEM_ELEMENT0(mgnoeud->get_lien_topologie(),mgnoeud,&femnoeud);
344     fem->ajouter_fem_element0(femele);
345     }
346     mgnoeud=maillage->get_suivant_noeud(it_nd);
347    
348     }
349     int dimsansgeo=0;
350     if (mg_geometrie==NULL)
351     {
352     if (maillage->get_nb_mg_tetra()!=0) dimsansgeo=3;
353     else if (maillage->get_nb_mg_triangle()!=0) dimsansgeo=2;
354     else if (maillage->get_nb_mg_segment()!=0) dimsansgeo=1;
355     }
356     if (mg_geometrie!=NULL)
357     if (strcmp(mg_geometrie->get_type_geometrie(),"VIRTUEL")==0)
358     {
359     if (maillage->get_nb_mg_tetra()+maillage->get_nb_mg_hexa()!=0) dimsansgeo=3;
360     else if (maillage->get_nb_mg_triangle()+maillage->get_nb_mg_quadrangle()!=0) dimsansgeo=2;
361     else if (maillage->get_nb_mg_segment()!=0) dimsansgeo=1;
362     }
363     int nbmgnoeud = maillage->get_nb_mg_noeud();
364    
365     LISTE_MG_SEGMENT::iterator it_seg;
366     MG_SEGMENT* mgseg=maillage->get_premier_segment(it_seg);
367     FEM_MAILLAGE_QUADRATIQUE_OUTILS ou ;
368     i=0;
369     while (mgseg)
370     {
371     mgseg->change_nouveau_numero(i++);
372     MG_NOEUD *noeud1=mgseg->get_noeud1();
373     MG_NOEUD *noeud2=mgseg->get_noeud2();
374     double *xyz1=noeud1->get_coord();
375     double *xyz2=noeud2->get_coord();
376     double xori=0.5*(xyz1[0]+xyz2[0]);
377     double yori=0.5*(xyz1[1]+xyz2[1]);
378     double zori=0.5*(xyz1[2]+xyz2[2]);
379     FEM_NOEUD* femnoeud;
380     if (mgseg->get_lien_topologie()!=NULL)
381     if (mgseg->get_lien_topologie()->get_dimension()==1)
382     {
383     MG_ARETE* arete=(MG_ARETE*)(mgseg->get_lien_topologie());
384     double t1,t2,tm;
385     arete->inverser(t1,xyz1);
386     arete->inverser(t2,xyz2);
387     if (arete->get_courbe()->est_periodique())
388     if (t1>t2) t2=t2+arete->get_courbe()->get_periode();
389     double t=0.5*(t1+t2);
390     double xyz[2],xyz3[2];
391     arete->evaluer(t,xyz);
392     double L3=0,L4=0,j=0.0;
393     double L1= ou.get_distance_curviligne(0,0.5,xyz1,xyz,xyz2);
394     double L2= ou.get_distance_curviligne(0.5,1,xyz1,xyz,xyz2);
395    
396     if (L1/L2 >1.01)
397     {
398     do
399     {
400     j+=0.0001;
401     tm=(-t+t1)*j+t;
402     arete->evaluer(tm,xyz3);
403     L3= ou.get_distance_curviligne(0,0.5,xyz1,xyz3,xyz2);
404     L4= ou.get_distance_curviligne(0.5,1,xyz1,xyz3,xyz2);
405    
406     }while (0.5-(L4/(L3+L4))>=0.0001);
407    
408     femnoeud= new FEM_NOEUD(mgseg->get_lien_topologie(),xyz3[0],xyz3[1],xyz3[2],xori,yori,zori);
409     }
410    
411     if (L2/L1 >1.01)
412     {
413     do
414     {
415     j+=0.0001;
416     tm=(-t+t2)*j+t;
417     arete->evaluer(tm,xyz3);
418     L3= ou.get_distance_curviligne(0,0.5,xyz1,xyz3,xyz2);
419     L4= ou.get_distance_curviligne(0.5,1,xyz1,xyz3,xyz2);
420    
421     }while (0.5-(L3/(L3+L4))>=0.0001);
422    
423     femnoeud= new FEM_NOEUD(mgseg->get_lien_topologie(),xyz3[0],xyz3[1],xyz3[2],xori,yori,zori);
424     }
425    
426     if(j == 0)
427     {
428     femnoeud= new FEM_NOEUD(mgseg->get_lien_topologie(),xyz[0],xyz[1],xyz[2],xori,yori,zori);
429     }
430    
431     }
432     if (mgseg->get_lien_topologie()!=NULL)
433     if (mgseg->get_lien_topologie()->get_dimension()==2)
434     {
435     MG_FACE* face=(MG_FACE*)(mgseg->get_lien_topologie());
436     double uv1[2],uv2[2];
437     OT_DECALAGE_PARAMETRE decalage(face->get_surface()->get_periode_u(),face->get_surface()->get_periode_v());
438     face->inverser(uv1,xyz1);
439     face->inverser(uv2,xyz2);
440     if (face->get_surface()->est_periodique_u()==1)
441     {
442     double eps=1e-10*face->get_surface()->get_periode_u();
443     if (uv1[0]<-eps) uv1[0]=uv1[0]+face->get_surface()->get_periode_u();
444     if (uv1[0]>face->get_surface()->get_periode_u()-eps) uv1[0]=uv1[0]-face->get_surface()->get_periode_u();
445     if (uv2[0]<-eps) uv2[0]=uv2[0]+face->get_surface()->get_periode_u();
446     if (uv2[0]>face->get_surface()->get_periode_u()-eps) uv2[0]=uv2[0]-face->get_surface()->get_periode_u();
447     }
448     if (face->get_surface()->est_periodique_v()==1)
449     {
450     double eps=1e-10*face->get_surface()->get_periode_v();
451     if (uv1[1]<-eps) uv1[1]=uv1[1]+face->get_surface()->get_periode_v();
452     if (uv1[1]>face->get_surface()->get_periode_v()-eps) uv1[1]=uv1[1]-face->get_surface()->get_periode_v();
453     if (uv2[1]<-eps) uv2[1]=uv2[1]+face->get_surface()->get_periode_v();
454     if (uv2[1]>face->get_surface()->get_periode_v()-eps) uv2[1]=uv2[1]-face->get_surface()->get_periode_v();
455     }
456     double du=decalage.calcul_decalage_parametre_u(uv1[0]);
457     double dv=decalage.calcul_decalage_parametre_v(uv1[1]);
458     double u1=decalage.decalage_parametre_u(uv1[0],du);
459     double v1=decalage.decalage_parametre_v(uv1[1],dv);
460     double u2=decalage.decalage_parametre_u(uv2[0],du);
461     double v2=decalage.decalage_parametre_v(uv2[1],dv);
462     double uv[2],uv3[2],xyz[3],xyz3[3];
463     double milieu[3];
464     milieu[0]=0.5*(xyz1[0]+xyz2[0]);
465     milieu[1]=0.5*(xyz1[1]+xyz2[1]);
466     milieu[2]=0.5*(xyz1[2]+xyz2[2]);
467     face->inverser(uv,milieu);
468    
469     face->evaluer(uv,xyz);
470     double L3=0,L4=0,j=0;
471     double L1= ou.get_distance_curviligne(0,0.5,xyz1,xyz,xyz2);
472     double L2= ou.get_distance_curviligne(0.5,1,xyz1,xyz,xyz2);
473     if (L1/L2 >1.01)
474     {
475     do
476     {
477     j+=0.0001;
478     if ((uv1[0]<uv[0]) && (uv2[0]<uv[0]) && (uv1[0]<uv2[0])) uv1[0]=uv1[0]+face->get_surface()->get_periode_u();
479     if ((uv1[1]<uv[1]) && (uv2[1]<uv[1]) && (uv1[1]<uv2[1])) uv1[1]=uv1[1]+face->get_surface()->get_periode_v();
480    
481     uv3[0]=(-uv[0]+uv1[0])*j+uv[0];
482     uv3[1]=(-uv[1]+uv1[1])*j+uv[1];
483    
484     face->evaluer(uv3,xyz3);
485     L3= ou.get_distance_curviligne(0,0.5,xyz1,xyz3,xyz2);
486     L4= ou.get_distance_curviligne(0.5,1,xyz1,xyz3,xyz2);
487    
488     }while ( 0.5-(L4/(L3+L4)) >=0.0001);
489    
490     femnoeud= new FEM_NOEUD(mgseg->get_lien_topologie(),xyz3[0],xyz3[1],xyz3[2],xori,yori,zori);
491     }
492    
493     if (L2/L1 >1.01)
494     {
495     do
496     {
497     j+=0.0001;
498     if ((uv1[0]<uv[0]) && (uv2[0]<uv[0]) && (uv2[0]<uv1[0])) uv2[0]=uv2[0]+face->get_surface()->get_periode_u();
499     if ((uv1[1]<uv[1]) && (uv2[1]<uv[1]) && (uv2[1]<uv1[1])) uv2[1]=uv2[1]+face->get_surface()->get_periode_v();
500    
501     uv3[0]=(-uv[0]+uv2[0])*j+uv[0];
502     uv3[1]=(-uv[1]+uv2[1])*j+uv[1];
503    
504     face->evaluer(uv3,xyz3);
505     L3= ou.get_distance_curviligne(0,0.5,xyz1,xyz3,xyz2);
506     L4= ou.get_distance_curviligne(0.5,1,xyz1,xyz3,xyz2);
507    
508     }while (0.5-(L3/(L3+L4)) >=0.0001);
509    
510     femnoeud= new FEM_NOEUD(mgseg->get_lien_topologie(),xyz3[0],xyz3[1],xyz3[2],xori,yori,zori);
511     }
512    
513     if (j == 0.0)
514     {
515     femnoeud= new FEM_NOEUD(mgseg->get_lien_topologie(),xyz[0],xyz[1],xyz[2],xori,yori,zori);
516     }
517    
518     }
519     if (mgseg->get_lien_topologie()!=NULL)
520     {
521     if (mgseg->get_lien_topologie()->get_dimension()==3)
522     {
523     double x=0.5*(xyz1[0]+xyz2[0]);
524     double y=0.5*(xyz1[1]+xyz2[1]);
525     double z=0.5*(xyz1[2]+xyz2[2]);
526     femnoeud= new FEM_NOEUD(mgseg->get_lien_topologie(),x,y,z,x,y,z);
527     }
528     }
529     else
530     {
531     double x=0.5*(xyz1[0]+xyz2[0]);
532     double y=0.5*(xyz1[1]+xyz2[1]);
533     double z=0.5*(xyz1[2]+xyz2[2]);
534     femnoeud= new FEM_NOEUD(mgseg->get_lien_topologie(),x,y,z,x,y,z);
535     }
536     femnoeud->change_numero(i+nbmgnoeud+1);
537     fem->ajouter_fem_noeud(femnoeud);
538     mgseg=maillage->get_suivant_segment(it_seg);
539     lstnoeuddirect.insert(lstnoeuddirect.end(),femnoeud);
540     }
541    
542     if ((mg_geometrie!=NULL) || ( (mg_geometrie==NULL) && (dimsansgeo==1)))
543     {
544    
545     mgseg=maillage->get_premier_segment(it_seg);
546     while (mgseg)
547     {
548     if (mgseg->get_lien_topologie()!=NULL)
549     if (mgseg->get_lien_topologie()->get_dimension()!=1)
550     {
551     mgseg=maillage->get_suivant_segment(it_seg);
552     continue;
553     }
554     if (mg_geometrie!=NULL)
555     if (mgseg->get_lien_topologie()==NULL)
556     if (dimsansgeo!=1)
557     {
558     mgseg=maillage->get_suivant_segment(it_seg);
559     continue;
560     }
561     FEM_NOEUD *tabnoeud[3];
562     tabnoeud[0]=lstnoeuddirect[mgseg->get_noeud1()->get_nouveau_numero()];
563     tabnoeud[1]=lstnoeuddirect[mgseg->get_nouveau_numero()+nbmgnoeud];
564     tabnoeud[2]=lstnoeuddirect[mgseg->get_noeud2()->get_nouveau_numero()];
565     FEM_SEGMENT3* seg=new FEM_SEGMENT3(mgseg,tabnoeud);
566     fem->ajouter_fem_element1(seg);
567     mgseg=maillage->get_suivant_segment(it_seg);
568     }
569     }
570     if ((mg_geometrie!=NULL) || ( (mg_geometrie==NULL) && (dimsansgeo==2)))
571     {
572     LISTE_MG_TRIANGLE::iterator it_tri;
573     MG_TRIANGLE* mgtri=maillage->get_premier_triangle(it_tri);
574     while (mgtri)
575     {
576     if (mgtri->get_lien_topologie()!=NULL) if (mgtri->get_lien_topologie()->get_dimension()!=2)
577     {
578     mgtri=maillage->get_suivant_triangle(it_tri);
579     continue;
580     }
581     if (mg_geometrie!=NULL)
582     if (mgtri->get_lien_topologie()==NULL)
583     if (dimsansgeo!=2)
584     {
585     mgtri=maillage->get_suivant_triangle(it_tri);
586     continue;
587     }
588     FEM_NOEUD *tabnoeud[6];
589     tabnoeud[0]=lstnoeuddirect[mgtri->get_noeud1()->get_nouveau_numero()];
590     tabnoeud[1]=lstnoeuddirect[mgtri->get_segment1()->get_nouveau_numero()+nbmgnoeud];
591     tabnoeud[2]=lstnoeuddirect[mgtri->get_noeud2()->get_nouveau_numero()];
592     tabnoeud[3]=lstnoeuddirect[mgtri->get_segment2()->get_nouveau_numero()+nbmgnoeud];
593     tabnoeud[4]=lstnoeuddirect[mgtri->get_noeud3()->get_nouveau_numero()];
594     tabnoeud[5]=lstnoeuddirect[mgtri->get_segment3()->get_nouveau_numero()+nbmgnoeud];
595     FEM_TRIANGLE6* tri=new FEM_TRIANGLE6(mgtri,tabnoeud);
596     fem->ajouter_fem_element2(tri);
597     mgtri=maillage->get_suivant_triangle(it_tri);
598     }
599     LISTE_MG_QUADRANGLE::iterator it_quad;
600     MG_QUADRANGLE* mgquad=maillage->get_premier_quadrangle(it_quad);
601     while (mgquad)
602     {
603     if (mgquad->get_lien_topologie()!=NULL) if (mgquad->get_lien_topologie()->get_dimension()!=2) {
604     mgquad=maillage->get_suivant_quadrangle(it_quad);
605     continue;
606     }
607     if (mg_geometrie!=NULL)
608     if (mgquad->get_lien_topologie()==NULL)
609     if (dimsansgeo!=2)
610     {
611     mgquad=maillage->get_suivant_quadrangle(it_quad);
612     continue;
613     }
614     FEM_NOEUD *tabnoeud[8];
615     tabnoeud[0]=lstnoeuddirect[mgquad->get_noeud1()->get_nouveau_numero()];
616     tabnoeud[1]=lstnoeuddirect[mgquad->get_segment1()->get_nouveau_numero()+nbmgnoeud];
617     tabnoeud[2]=lstnoeuddirect[mgquad->get_noeud2()->get_nouveau_numero()];
618     tabnoeud[3]=lstnoeuddirect[mgquad->get_segment2()->get_nouveau_numero()+nbmgnoeud];
619     tabnoeud[4]=lstnoeuddirect[mgquad->get_noeud3()->get_nouveau_numero()];
620     tabnoeud[5]=lstnoeuddirect[mgquad->get_segment3()->get_nouveau_numero()+nbmgnoeud];
621     tabnoeud[6]=lstnoeuddirect[mgquad->get_noeud4()->get_nouveau_numero()];
622     tabnoeud[7]=lstnoeuddirect[mgquad->get_segment4()->get_nouveau_numero()+nbmgnoeud];
623     FEM_QUADRANGLE8* quad=new FEM_QUADRANGLE8(mgquad,tabnoeud);
624     fem->ajouter_fem_element2(quad);
625     mgquad=maillage->get_suivant_quadrangle(it_quad);
626     }
627     }
628     if ((mg_geometrie!=NULL) || ( (mg_geometrie==NULL) && (dimsansgeo==3)))
629     {
630    
631     int nbpas=(int)param.get_valeur("Quadratisation_pas");
632     int recal=(int)param.get_valeur("Quadratisation_jmin");
633     double dismin=param.get_valeur("Quadratisation_dis");
634     LISTE_MG_TETRA::iterator it_tetra;
635     MG_TETRA* mgtetra=maillage->get_premier_tetra(it_tetra);
636     while (mgtetra)
637     {
638     FEM_NOEUD *tabnoeud[10];
639     tabnoeud[0]=lstnoeuddirect[mgtetra->get_noeud1()->get_nouveau_numero()];
640     tabnoeud[1]=lstnoeuddirect[maillage->get_mg_segment(mgtetra->get_noeud1()->get_id(),mgtetra->get_noeud2()->get_id())->get_nouveau_numero()+nbmgnoeud];
641     tabnoeud[2]=lstnoeuddirect[mgtetra->get_noeud2()->get_nouveau_numero()];
642     tabnoeud[3]=lstnoeuddirect[maillage->get_mg_segment(mgtetra->get_noeud2()->get_id(),mgtetra->get_noeud3()->get_id())->get_nouveau_numero()+nbmgnoeud];
643     tabnoeud[4]=lstnoeuddirect[mgtetra->get_noeud3()->get_nouveau_numero()];
644     tabnoeud[5]=lstnoeuddirect[maillage->get_mg_segment(mgtetra->get_noeud3()->get_id(),mgtetra->get_noeud1()->get_id())->get_nouveau_numero()+nbmgnoeud];
645     tabnoeud[6]=lstnoeuddirect[maillage->get_mg_segment(mgtetra->get_noeud1()->get_id(),mgtetra->get_noeud4()->get_id())->get_nouveau_numero()+nbmgnoeud];
646     tabnoeud[7]=lstnoeuddirect[maillage->get_mg_segment(mgtetra->get_noeud2()->get_id(),mgtetra->get_noeud4()->get_id())->get_nouveau_numero()+nbmgnoeud];
647     tabnoeud[8]=lstnoeuddirect[maillage->get_mg_segment(mgtetra->get_noeud3()->get_id(),mgtetra->get_noeud4()->get_id())->get_nouveau_numero()+nbmgnoeud];
648     tabnoeud[9]=lstnoeuddirect[mgtetra->get_noeud4()->get_nouveau_numero()];
649     FEM_TETRA10* tet=new FEM_TETRA10(mgtetra,tabnoeud);
650     fem->ajouter_fem_element3(tet);
651     if (recal==1)
652     {
653     FEM_MAILLAGE_QUADRATIQUE_OUTILS ou;
654     double dis=ou.get_distorsion2(tet);
655     if (dis<dismin-1e-12)
656     {
657     OT_VECTEUR_3D vec[10];
658     OT_VECTEUR_3D depart[10];
659     for (int i=0;i<9;i++)
660     {
661     if ((i==0) || (i==2) || (i==4)) continue;
662     double* xyzori=tet->get_fem_noeud(i)->get_coord_ori();
663     double* xyz=tet->get_fem_noeud(i)->get_coord();
664     vec[i].change_x(xyz[0]-xyzori[0]);
665     vec[i].change_y(xyz[1]-xyzori[1]);
666     vec[i].change_z(xyz[2]-xyzori[2]);
667     depart[i].change_x(xyz[0]);
668     depart[i].change_y(xyz[1]);
669     depart[i].change_z(xyz[2]);
670     }
671     int pas=0;
672     while (dis<dismin-1e-12)
673     {
674     for (int i=0;i<9;i++)
675     {
676     if ((i==0) || (i==2) || (i==4)) continue;
677     double* xyz=tet->get_fem_noeud(i)->get_coord();
678     tet->get_fem_noeud(i)->change_x(depart[i].get_x()-pas*1./1./nbpas*vec[i].get_x());
679     tet->get_fem_noeud(i)->change_y(depart[i].get_y()-pas*1./1./nbpas*vec[i].get_y());
680     tet->get_fem_noeud(i)->change_z(depart[i].get_z()-pas*1./1./nbpas*vec[i].get_z());
681     }
682     dis=ou.get_distorsion2(tet);
683     pas++;
684     }
685     char mess[255];
686     sprintf(mess," Tetra %lu replacé à %.2lf%% de sa position optimale",tet->get_id(),(pas-1.)*1./nbpas*100.);
687     affiche(mess);
688     }
689     }
690     mgtetra=maillage->get_suivant_tetra(it_tetra);
691    
692    
693     }
694     LISTE_MG_HEXA::iterator it_hexa;
695     MG_HEXA* mghex=maillage->get_premier_hexa(it_hexa);
696     while (mghex)
697     {
698     FEM_NOEUD *tabnoeud[20];
699     tabnoeud[0]=lstnoeuddirect[mghex->get_noeud1()->get_nouveau_numero()];
700     tabnoeud[1]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud1()->get_id(),mghex->get_noeud2()->get_id())->get_nouveau_numero()+nbmgnoeud];
701     tabnoeud[2]=lstnoeuddirect[mghex->get_noeud2()->get_nouveau_numero()];
702     tabnoeud[3]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud2()->get_id(),mghex->get_noeud3()->get_id())->get_nouveau_numero()+nbmgnoeud];
703     tabnoeud[4]=lstnoeuddirect[mghex->get_noeud3()->get_nouveau_numero()];
704     tabnoeud[5]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud3()->get_id(),mghex->get_noeud4()->get_id())->get_nouveau_numero()+nbmgnoeud];
705     tabnoeud[6]=lstnoeuddirect[mghex->get_noeud4()->get_nouveau_numero()];
706     tabnoeud[7]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud4()->get_id(),mghex->get_noeud1()->get_id())->get_nouveau_numero()+nbmgnoeud];
707     tabnoeud[8]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud1()->get_id(),mghex->get_noeud5()->get_id())->get_nouveau_numero()+nbmgnoeud];
708     tabnoeud[9]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud2()->get_id(),mghex->get_noeud6()->get_id())->get_nouveau_numero()+nbmgnoeud];
709     tabnoeud[10]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud3()->get_id(),mghex->get_noeud7()->get_id())->get_nouveau_numero()+nbmgnoeud];
710     tabnoeud[11]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud4()->get_id(),mghex->get_noeud8()->get_id())->get_nouveau_numero()+nbmgnoeud];
711     tabnoeud[12]=lstnoeuddirect[mghex->get_noeud5()->get_nouveau_numero()];
712     tabnoeud[13]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud5()->get_id(),mghex->get_noeud6()->get_id())->get_nouveau_numero()+nbmgnoeud];
713     tabnoeud[14]=lstnoeuddirect[mghex->get_noeud6()->get_nouveau_numero()];
714     tabnoeud[15]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud6()->get_id(),mghex->get_noeud7()->get_id())->get_nouveau_numero()+nbmgnoeud];
715     tabnoeud[16]=lstnoeuddirect[mghex->get_noeud7()->get_nouveau_numero()];
716     tabnoeud[17]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud7()->get_id(),mghex->get_noeud8()->get_id())->get_nouveau_numero()+nbmgnoeud];
717     tabnoeud[18]=lstnoeuddirect[mghex->get_noeud8()->get_nouveau_numero()];
718     tabnoeud[19]=lstnoeuddirect[maillage->get_mg_segment(mghex->get_noeud8()->get_id(),mghex->get_noeud5()->get_id())->get_nouveau_numero()+nbmgnoeud];
719     FEM_HEXA20* hex=new FEM_HEXA20(mghex,tabnoeud);
720     fem->ajouter_fem_element3(hex);
721     mghex=maillage->get_suivant_hexa(it_hexa);
722     }
723     LISTE_MG_PENTA::iterator it_pen;
724     MG_PENTA* mgpen=maillage->get_premier_penta(it_pen);
725     while (mgpen)
726     {
727     FEM_NOEUD *tabnoeud[15];
728     tabnoeud[0]=lstnoeuddirect[mgpen->get_noeud1()->get_nouveau_numero()];
729     tabnoeud[1]=lstnoeuddirect[maillage->get_mg_segment(mgpen->get_noeud1()->get_id(),mgpen->get_noeud2()->get_id())->get_nouveau_numero()+nbmgnoeud];
730     tabnoeud[2]=lstnoeuddirect[mgpen->get_noeud2()->get_nouveau_numero()];
731     tabnoeud[3]=lstnoeuddirect[maillage->get_mg_segment(mgpen->get_noeud2()->get_id(),mgpen->get_noeud3()->get_id())->get_nouveau_numero()+nbmgnoeud];
732     tabnoeud[4]=lstnoeuddirect[mgpen->get_noeud3()->get_nouveau_numero()];
733     tabnoeud[5]=lstnoeuddirect[maillage->get_mg_segment(mgpen->get_noeud3()->get_id(),mgpen->get_noeud1()->get_id())->get_nouveau_numero()+nbmgnoeud];
734     tabnoeud[6]=lstnoeuddirect[maillage->get_mg_segment(mgpen->get_noeud1()->get_id(),mgpen->get_noeud4()->get_id())->get_nouveau_numero()+nbmgnoeud];
735     tabnoeud[7]=lstnoeuddirect[maillage->get_mg_segment(mgpen->get_noeud2()->get_id(),mgpen->get_noeud5()->get_id())->get_nouveau_numero()+nbmgnoeud];
736     tabnoeud[8]=lstnoeuddirect[maillage->get_mg_segment(mgpen->get_noeud3()->get_id(),mgpen->get_noeud6()->get_id())->get_nouveau_numero()+nbmgnoeud];
737     tabnoeud[9]=lstnoeuddirect[mgpen->get_noeud4()->get_nouveau_numero()];
738     tabnoeud[10]=lstnoeuddirect[maillage->get_mg_segment(mgpen->get_noeud4()->get_id(),mgpen->get_noeud5()->get_id())->get_nouveau_numero()+nbmgnoeud];
739     tabnoeud[11]=lstnoeuddirect[mgpen->get_noeud5()->get_nouveau_numero()];
740     tabnoeud[12]=lstnoeuddirect[maillage->get_mg_segment(mgpen->get_noeud5()->get_id(),mgpen->get_noeud6()->get_id())->get_nouveau_numero()+nbmgnoeud];
741     tabnoeud[13]=lstnoeuddirect[mgpen->get_noeud6()->get_nouveau_numero()];
742     tabnoeud[14]=lstnoeuddirect[maillage->get_mg_segment(mgpen->get_noeud4()->get_id(),mgpen->get_noeud6()->get_id())->get_nouveau_numero()+nbmgnoeud];
743     FEM_PENTA15* pen=new FEM_PENTA15(mgpen,tabnoeud);
744     fem->ajouter_fem_element3(pen);
745     mgpen=maillage->get_suivant_penta(it_pen);
746     }
747     }
748     }
749    
750    
751    
752    
753     void MAILLEUR_FEM::optimise_numerotation(FEM_MAILLAGE* fem)
754     {
755     FEM_NOEUD* noeud=fem->get_fem_noeud(0);
756     FEM_GRAPH_NOEUD *graph;
757     graph=new FEM_GRAPH_NOEUD(noeud,fem);
758     int e=graph->get_excentricite();
759     TPL_MAP_ENTITE<FEM_NOEUD*> dernierniveau=graph->get_dernier_niveau();
760     int nbdernierniveau=dernierniveau.get_nb();
761     for (int i=0;i<nbdernierniveau;i++)
762     {
763     FEM_GRAPH_NOEUD* graphtmp=new FEM_GRAPH_NOEUD(dernierniveau.get(i),fem);
764     int etmp=graphtmp->get_excentricite();
765     if (etmp>e)
766     {
767     delete graph;
768     graph=graphtmp;
769     dernierniveau=graph->get_dernier_niveau();
770     nbdernierniveau=dernierniveau.get_nb();
771     i=-1;
772     e=etmp;
773     }
774     else delete graphtmp;
775     }
776     int numero=fem->get_nb_fem_noeud();
777     for (int i=0;i<e+1;i++)
778     {
779     TPL_MAP_ENTITE<FEM_NOEUD*> niveau=graph->get_niveau(i);
780     int nbnoeud=niveau.get_nb();
781     for (int j=0;j<nbnoeud;j++)
782     {
783     FEM_NOEUD* noeud=niveau.get(j);
784     noeud->change_numero_opt(numero);
785     numero--;
786     }
787     }
788    
789    
790     delete graph;
791     }
792    
793    
794     void MAILLEUR_FEM::copie_numerotation_opt(FEM_MAILLAGE* fem)
795     {
796     int nbnoeud=fem->get_nb_fem_noeud();
797     for (int i=0;i<nbnoeud;i++)
798     {
799     FEM_NOEUD* noeud=fem->get_fem_noeud(i);
800     noeud->change_numero(noeud->get_numero_opt());
801     }
802    
803     }
804