ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/mailleur2d_structure.cpp
Revision: 1053
Committed: Fri Oct 16 16:17:23 2020 UTC (4 years, 6 months ago) by francois
File size: 9944 byte(s)
Log Message:
mailleur2d_structure de quadrangle

File Contents

# User Rev Content
1 francois 1053 #include <stdio.h>
2     #include <math.h>
3     #include "mg_gestionnaire.h"
4     #include "mailleur2d_structure.h"
5     #include "mg_maillage.h"
6     #include "fem_maillage.h"
7     #include "fem_hexa20.h"
8     #include "mg_geometrie.h"
9     #include "mg_sommet_noeud.h"
10     #include "mg_arete_element.h"
11     #include "mg_face_element.h"
12     #include "mg_volume_element.h"
13     #include "tpl_octree.h"
14     #include "mg_maillage_outils.h"
15     #include <vector>
16    
17    
18     struct integrale
19     {
20     double ti;
21     double li;
22     };
23    
24     MAILLEUR2D_STRUCTURE::MAILLEUR2D_STRUCTURE(double xmintmp, double xmaxtmp, double ymintmp, double ymaxtmp, int ntmp, int mtmp, MG_GESTIONNAIRE* gesttmp,int avecgeo,double unit):MAILLEUR(false),xmin(xmintmp),xmax(xmaxtmp),ymin(ymintmp),ymax(ymaxtmp),n(ntmp),m(mtmp),gest(gesttmp),geovirtuel(avecgeo),densitevariable(false),unite(unit),nbfonction(0)
25     {
26     }
27     MAILLEUR2D_STRUCTURE::MAILLEUR2D_STRUCTURE( double xmintmp, double xmaxtmp, double ymintmp, double ymaxtmp, double (*fx)(double),double (*fy)(double),MG_GESTIONNAIRE* gesttmp,int avecgeo,double unit):MAILLEUR(false),xmin(xmintmp),xmax(xmaxtmp),ymin(ymintmp),ymax(ymaxtmp),gest(gesttmp),geovirtuel(avecgeo),enx(fx),eny(fy),densitevariable(true),unite(unit),nbfonction(0)
28     {
29     }
30    
31     MAILLEUR2D_STRUCTURE::MAILLEUR2D_STRUCTURE( MAILLEUR2D_STRUCTURE& mdd):MAILLEUR(mdd),
32     xmin(mdd.xmin),xmax(mdd.xmax),ymin(mdd.ymin),ymax(mdd.ymax),n(mdd.n),m(mdd.m),gest(mdd.gest),nbfonction(mdd.nbfonction)
33     {
34     }
35    
36     MAILLEUR2D_STRUCTURE::~MAILLEUR2D_STRUCTURE()
37     {
38     }
39    
40    
41     void MAILLEUR2D_STRUCTURE::discretise(std::vector<double> &tab,double min,double max,int &nb,double (*en)(double))
42     {
43     double t1=0;
44     double t2=1;
45     std::vector<integrale> tab_integrale;
46     double ti,tii,t;
47     double longueur_metrique=0.;
48     tii=t1;
49     double dt=1./pas;
50     do
51     {
52     ti=tii;
53     tii=ti+(t2-t1)*dt;
54     if (tii>t2) tii=t2;
55     t=0.7886751345*ti+0.2113248654*tii;
56     double facteur1=(max-min)/(*en)(t);
57     t=0.7886751345*tii+0.2113248654*ti;
58     double facteur2=(max-min)/(*en)(t);
59     longueur_metrique=longueur_metrique+0.5*(tii-ti)*(facteur1+facteur2);
60     integrale pas;
61     pas.ti=ti;
62     pas.li=longueur_metrique;
63     tab_integrale.insert(tab_integrale.end(),pas);
64     }
65     while (tii<t2);
66     nb=(int)floor(longueur_metrique);
67     if (longueur_metrique-floor(longueur_metrique)>0.5) nb++;
68     if (nb<1) nb=1;
69     double valeur_cible=longueur_metrique/nb;
70     tab.push_back(min);
71     int nbcree=0;
72     int nb_pas_integrale=tab_integrale.size();
73     for (int i=0;i<nb_pas_integrale;i++)
74     {
75     while ((tab_integrale[i].li>(nbcree+1)*valeur_cible) && (nbcree<nb-1))
76     {
77     double ti=tab_integrale[i].ti;
78     double tii=t2;
79     double li=0.;
80     if (i!=0) li=tab_integrale[i-1].li;
81     double lii=tab_integrale[i].li;
82     if (i!=nb_pas_integrale-1) tii=tab_integrale[i+1].ti;
83     t=ti+(tii-ti)/(lii-li)*((nbcree+1)*valeur_cible-li);
84     double x=min+t*(max-min);
85     tab.push_back(x);
86     nbcree++;
87     }
88     }
89     tab.push_back(max);
90    
91     }
92     void MAILLEUR2D_STRUCTURE::ajouter_fonction_geometrie(int nb)
93     {
94     nbfonction=nb;
95     }
96    
97     int MAILLEUR2D_STRUCTURE::maille(MG_GROUPE_TOPOLOGIQUE* mggt)
98     {
99     affiche((char*)"Création d'un maillage MG structuré");
100     double intx, inty;
101     std::vector<double> tabx,taby;
102     if (densitevariable==false)
103     {
104     intx=fabs(xmax-xmin)/n;
105     inty=fabs(ymax-ymin)/m;
106     }
107     else
108     {
109     discretise(tabx,xmin,xmax,n,enx);
110     discretise(taby,ymin,ymax,m,eny);
111     }
112    
113     MG_GEOMETRIE* geo=NULL;
114     MG_SOMMET_NOEUD *som1=NULL;
115     MG_SOMMET_NOEUD *som2=NULL;
116     MG_SOMMET_NOEUD *som3=NULL;
117     MG_SOMMET_NOEUD *som4=NULL;
118     MG_ARETE_ELEMENT* are1=NULL;
119     MG_ARETE_ELEMENT* are2=NULL;
120     MG_ARETE_ELEMENT* are3=NULL;
121     MG_ARETE_ELEMENT* are4=NULL;
122     MG_FACE_ELEMENT* face1=NULL;
123     if (geovirtuel==true)
124     {
125     geo=new MG_GEOMETRIE((char*)"VIRTUEL",(char*)"VIRTUEL");
126     geo->change_valeur_unite(unite);
127     gest->ajouter_mg_geometrie(geo);
128     som1=new MG_SOMMET_NOEUD(NULL);
129     geo->ajouter_mg_sommet(som1);
130     som2=new MG_SOMMET_NOEUD(NULL);
131     geo->ajouter_mg_sommet(som2);
132     som3=new MG_SOMMET_NOEUD(NULL);
133     geo->ajouter_mg_sommet(som3);
134     som4=new MG_SOMMET_NOEUD(NULL);
135     geo->ajouter_mg_sommet(som4);
136     are1=new MG_ARETE_ELEMENT();
137     geo->ajouter_mg_arete(are1);
138     are2=new MG_ARETE_ELEMENT();
139     geo->ajouter_mg_arete(are2);
140     are3=new MG_ARETE_ELEMENT();
141     geo->ajouter_mg_arete(are3);
142     are4=new MG_ARETE_ELEMENT();
143     geo->ajouter_mg_arete(are4);
144     face1=new MG_FACE_ELEMENT();
145     geo->ajouter_mg_face(face1);
146     for (int i=0;i<nbfonction;i++)
147     {
148     MG_GEOM_FONCTION *fonc=new MG_GEOM_FONCTION(0);
149     geo->ajouter_mg_geom_fonction(fonc);
150     }
151     }
152     MG_MAILLAGE* mai=new MG_MAILLAGE(geo);
153     gest->ajouter_mg_maillage(mai);
154     std::vector<MG_NOEUD*> vecnod;
155     for (int j=0; j<=m;j++)
156     {
157     for (int i=0; i<=n;i++)
158     {
159     MG_ELEMENT_TOPOLOGIQUE *topo=NULL;
160     if (geovirtuel)
161     {
162     if ((i==0) && (j==0)) topo=som1;
163     else if ((i==0) && (j==m)) topo=som2;
164     else if ((i==0) && (j==m)) topo=som3;
165     else if ((i==0) && (j==0)) topo=som4;
166     else if (j==0) topo=are1;
167     else if (i==n) topo=are2;
168     else if (j==m) topo=are3;
169     else if (j==0) topo=are4;
170     else topo=face1;
171     }
172     double x,y;
173     if (densitevariable==false)
174     {
175     y=ymin+(j*inty);
176     x=xmin+(i*intx);
177     }
178     else
179     {
180     x=tabx[i];
181     y=taby[j];
182     }
183     MG_NOEUD* nou=mai->ajouter_mg_noeud(topo,x,y,0.,MAGIC::ORIGINE::MAILLEUR_AUTO);
184     vecnod.push_back(nou);
185     if (geovirtuel)
186     {
187     if ((i==0) && (j==0)) som1->change_mg_noeud(nou);
188     else if ((i==n) && (j==0) ) som2->change_mg_noeud(nou);
189     else if ((i==n) && (j==m) ) som3->change_mg_noeud(nou);
190     else if ((i==0) && (j==m) ) som4->change_mg_noeud(nou);
191     }
192     }
193     }
194     if (geovirtuel)
195     {
196     for(int i=0; i<n;i++)
197     {
198     int j=0;
199     MG_NOEUD* nod1=vecnod[i+((n+1)*j)];
200     MG_NOEUD* nod2=vecnod[(i+1)+((n+1)*j)];
201     mai->ajouter_mg_segment(are1,nod1,nod2,MAGIC::ORIGINE::MAILLEUR_AUTO);
202     j=m;
203     MG_NOEUD* nod3=vecnod[i+((n+1)*j)];
204     MG_NOEUD* nod4=vecnod[(i+1)+((n+1)*j)];
205     mai->ajouter_mg_segment(are3,nod3,nod4,MAGIC::ORIGINE::MAILLEUR_AUTO);
206     }
207     for(int j=0; j<m;j++)
208     {
209     int i=0;
210     MG_NOEUD* nod1=vecnod[i+((n+1)*j)];
211     MG_NOEUD* nod2=vecnod[i+((n+1)*(j+1))];
212     mai->ajouter_mg_segment(are4,nod1,nod2,MAGIC::ORIGINE::MAILLEUR_AUTO);
213     i=n;
214     MG_NOEUD* nod3=vecnod[i+((n+1)*j)];
215     MG_NOEUD* nod4=vecnod[i+((n+1)*(j+1))];
216     mai->ajouter_mg_segment(are2,nod3,nod4,MAGIC::ORIGINE::MAILLEUR_AUTO);
217     }
218     }
219     for (int i=0;i<n;i++)
220     for (int j=0;j<m;j++)
221     {
222     MG_NOEUD* nod1=vecnod[i+((n+1)*j)];
223     MG_NOEUD* nod2=vecnod[(i+1)+((n+1)*j)];
224     MG_NOEUD* nod3=vecnod[(i+1)+((n+1)*(j+1))];
225     MG_NOEUD* nod4=vecnod[i+((n+1)*(j+1))];
226     mai->ajouter_mg_quadrangle(face1,nod1,nod4,nod3,nod2,MAGIC::ORIGINE::MAILLEUR_AUTO);
227     }
228    
229     BOITE_3D b(xmin,ymin,0.,xmax,ymax,0.);
230     mai->change_param_structure(b,n,m,0);
231     /*affiche((char*)"Création d'un maillage FEM structuré");
232     FEM_MAILLAGE* fem = new FEM_MAILLAGE(geo, mai,1);
233     gest->ajouter_fem_maillage(fem);
234     fem->construire();*/
235     return OK;
236     }
237    
238    
239     void MAILLEUR2D_STRUCTURE::ajuste(class MG_MAILLAGE* mai, MG_MAILLAGE* maibase, int nb)
240     {
241     LISTE_MG_NOEUD::iterator it;
242     double xmin=1e300,ymin=1e300,zmin=1e300;
243     double xmax=-1e300,ymax=-1e300,zmax=-1e300;
244     TPL_LISTE_ENTITE<MG_NOEUD*> lst;
245     for (MG_NOEUD* no=maibase->get_premier_noeud(it);no!=NULL;no=maibase->get_suivant_noeud(it))
246     {
247     if (no->get_x()<xmin) xmin=no->get_x();
248     if (no->get_y()<ymin) ymin=no->get_y();
249     if (no->get_z()<zmin) zmin=no->get_z();
250     if (no->get_x()>xmax) xmax=no->get_x();
251     if (no->get_y()>ymax) ymax=no->get_y();
252     if (no->get_z()>zmax) zmax=no->get_z();
253     lst.ajouter(no);
254     }
255     BOITE_3D b(xmin,ymin,zmin,xmax,ymax,zmax);
256     b.change_grosseur(1.05);
257     TPL_OCTREE<MG_TRIANGLE*,MG_NOEUD*> oc;
258     oc.initialiser(&lst,1,b.get_xmin(),b.get_xmin(),b.get_zmin(),b.get_xmax(),b.get_xmax(),b.get_zmax());
259     LISTE_MG_TRIANGLE::iterator itt;
260     for (MG_TRIANGLE* tri=maibase->get_premier_triangle(itt);tri!=NULL;tri=maibase->get_suivant_triangle(itt))
261     oc.inserer(tri);
262    
263     for (MG_NOEUD* no=mai->get_premier_noeud(it);no!=NULL;no=mai->get_suivant_noeud(it))
264     {
265     double x=no->get_x();
266     double y=no->get_y();
267     double z=no->get_z();
268     TPL_MAP_ENTITE<MG_TRIANGLE*> lsttrouve;
269     oc.rechercher(x,y,z,1e-7,lsttrouve);
270     TPL_MAP_ENTITE<MG_TRIANGLE*>::ITERATEUR itf;
271     bool trouve=false;
272     for (MG_TRIANGLE* tri=lsttrouve.get_premier(itf);tri!=NULL;tri=lsttrouve.get_suivant(itf))
273     {
274     MG_MAILLAGE_OUTILS ot;
275     int res=ot.estdansletriangle(tri,x,y,z);
276     int val=ot.compare_etat_triangle(res,MG_MAILLAGE_OUTILS::INTERIEUR);
277     if (val>0) trouve=true;
278     }
279     if (trouve) no->change_nouveau_numero(1); else no->change_nouveau_numero(0);
280     }
281     LISTE_MG_QUADRANGLE::iterator itq;
282     std::vector<unsigned long> listeasupprimer;
283     for (MG_QUADRANGLE* quad=mai->get_premier_quadrangle(itq);quad!=NULL;quad=mai->get_suivant_quadrangle(itq))
284     {
285     MG_NOEUD* no1=quad->get_noeud1();
286     MG_NOEUD* no2=quad->get_noeud2();
287     MG_NOEUD* no3=quad->get_noeud3();
288     MG_NOEUD* no4=quad->get_noeud4();
289     int val=no1->get_nouveau_numero()+no2->get_nouveau_numero()+no3->get_nouveau_numero()+no4->get_nouveau_numero();
290     if (val<nb)
291     listeasupprimer.push_back(quad->get_id());
292     }
293     for (int i=0;i<listeasupprimer.size();i++)
294     mai->supprimer_mg_quadrangleid(listeasupprimer[i]);
295    
296    
297    
298     }