ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mtu/src/mg_face.cpp
Revision: 1180
Committed: Fri Sep 20 20:10:22 2024 UTC (8 months ago) by francois
File size: 14870 byte(s)
Log Message:
Correction dans la quadratisation sans geometrie

File Contents

# User Rev Content
1 francois 1158 //####//------------------------------------------------------------
2     //####//------------------------------------------------------------
3     //####// MAGiC
4     //####// Jean Christophe Cuilliere et Vincent FRANCOIS
5     //####// Departement de Genie Mecanique - UQTR
6     //####//------------------------------------------------------------
7     //####// MAGIC est un projet de recherche de l equipe ERICCA
8     //####// du departement de genie mecanique de l Universite du Quebec a Trois Rivieres
9     //####// http://www.uqtr.ca/ericca
10     //####// http://www.uqtr.ca/
11     //####//------------------------------------------------------------
12     //####//------------------------------------------------------------
13     //####//
14     //####// mg_face.cpp
15     //####//
16     //####//------------------------------------------------------------
17     //####//------------------------------------------------------------
18     //####// COPYRIGHT 2000-2024
19     //####// jeu 13 jun 2024 11:58:54 EDT
20     //####//------------------------------------------------------------
21     //####//------------------------------------------------------------
22 francois 283
23    
24     #include "gestionversion.h"
25     #include <math.h>
26     #include "mg_face.h"
27     #include "vct_face.h"
28 francois 375 #include "mg_definition.h"
29 francois 283 #include "ot_mathematique.h"
30    
31 couturad 814 MG_FACE::MG_FACE(std::string idori,unsigned long num,MG_SURFACE* srf,int sens):MG_ELEMENT_TOPOLOGIQUE(num,idori),surface(srf),orientation(sens),vect(NULL),nb_pole(-1)
32 francois 283 {
33     }
34    
35 couturad 814 MG_FACE::MG_FACE(std::string idori,MG_SURFACE* srf,int sens):MG_ELEMENT_TOPOLOGIQUE(idori),surface(srf),orientation(sens),vect(NULL),nb_pole(-1)
36 francois 283 {
37     }
38    
39 couturad 814 MG_FACE::MG_FACE(MG_FACE& mdd):MG_ELEMENT_TOPOLOGIQUE(mdd),lst_boucle(mdd.lst_boucle),surface(mdd.surface),orientation(mdd.orientation),vect(NULL),nb_pole(mdd.nb_pole)
40 francois 283 {
41     }
42    
43     MG_FACE::~MG_FACE()
44     {
45     if (vect!=NULL) delete vect;
46     }
47    
48     void MG_FACE::ajouter_mg_boucle(class MG_BOUCLE* mgbou)
49     {
50     lst_boucle.insert(lst_boucle.end(),mgbou);
51     }
52    
53     void MG_FACE::supprimer_mg_boucle(class MG_BOUCLE* mgbou)
54     {
55     std::vector<MG_BOUCLE*>::iterator i;
56     for (i=lst_boucle.begin();i!=lst_boucle.end();i++)
57     {
58     if ((*i)==mgbou)
59     {
60     lst_boucle.erase(i);
61     return;
62     }
63     }
64     }
65    
66    
67     int MG_FACE::get_nb_mg_boucle(void)
68     {
69     return lst_boucle.size();
70     }
71    
72     MG_BOUCLE* MG_FACE::get_mg_boucle(int num)
73     {
74     return lst_boucle[num];
75     }
76    
77    
78     void MG_FACE::ajouter_mg_coface(class MG_COFACE* coface)
79     {
80     lst_coface.insert(lst_coface.end(),coface);
81     }
82    
83     int MG_FACE::get_nb_mg_coface(void)
84     {
85     return lst_coface.size();
86     }
87    
88    
89     void MG_FACE::supprimer_mg_coface(class MG_COFACE* coface)
90     {
91     std::vector<MG_COFACE*>::iterator i;
92     for (i=lst_coface.begin();i!=lst_coface.end();i++)
93     {
94     if ((*i)==coface)
95     {
96     lst_coface.erase(i);
97     return;
98     }
99     }
100     }
101    
102    
103    
104     MG_COFACE* MG_FACE::get_mg_coface(int num)
105     {
106     return lst_coface[num];
107     }
108    
109     MG_SURFACE* MG_FACE::get_surface(void)
110     {
111     return surface;
112     }
113    
114     void MG_FACE::get_topologie_sousjacente(TPL_MAP_ENTITE<MG_ELEMENT_TOPOLOGIQUE*> *lst)
115     {
116     int nb=lst_boucle.size();
117     for (int i=0;i<nb;i++)
118     {
119     MG_BOUCLE* bou=lst_boucle[i];
120     int nb2=bou->get_nb_mg_coarete();
121     for (int j=0;j<nb2;j++)
122     {
123     MG_ARETE* are=bou->get_mg_coarete(j)->get_arete();
124     lst->ajouter(are);
125     are->get_topologie_sousjacente(lst);
126     }
127     }
128     }
129    
130     int MG_FACE::get_dimension(void)
131     {
132     return 2;
133     }
134    
135 couturad 906 int MG_FACE::get_type(void)
136     {
137     return TYPE_ELEMENT_TOPOLOGIQUE::FACE;
138     }
139    
140 francois 283 int MG_FACE::get_orientation(void)
141     {
142     return orientation;
143     }
144    
145 francois 1180 bool MG_FACE::est_une_topo_element(void)
146 francois 576 {
147     return false;
148     }
149    
150 francois 283 int MG_FACE::valide_parametre_u(double& u)
151     {
152     if (surface->est_periodique_u()) return 1;
153     double param=u;
154     if (orientation!=MEME_SENS) param=surface->get_umin()+surface->get_umax()-param;
155     double param_min,param_max;
156     param_min=surface->get_umin();
157     param_max=surface->get_umax();
158     if (param<param_min)
159     {
160     u=param_min;
161     if (orientation!=MEME_SENS) u=surface->get_umin()+surface->get_umax()-u;
162     return 0;
163     }
164     if (param>param_max)
165     {
166     u=param_max;
167     if (orientation!=MEME_SENS) u=surface->get_umin()+surface->get_umax()-u;
168     return 0;
169     }
170     return 1;
171     }
172     int MG_FACE::valide_parametre_v(double& v)
173     {
174     if (surface->est_periodique_v()) return 1;
175     double param=v;
176     double param_min,param_max;
177     param_min=surface->get_vmin();
178     param_max=surface->get_vmax();
179     if (param<param_min)
180     {
181     v=param_min;
182     return 0;
183     }
184     if (param>param_max)
185     {
186     v=param_max;
187     return 0;
188     }
189     return 1;
190     }
191    
192     void MG_FACE::evaluer(double *uv,double *xyz)
193     {
194     double param[2]={uv[0],uv[1]};
195     if (orientation!=MEME_SENS) param[0]=surface->get_umin()+surface->get_umax()-param[0];
196     surface->evaluer(param,xyz);
197     }
198    
199     void MG_FACE::deriver(double *uv,double *xyzdu, double *xyzdv)
200     {
201     double param[2]={uv[0],uv[1]};
202     if (orientation!=MEME_SENS) param[0]=surface->get_umin()+surface->get_umax()-param[0];
203     surface->deriver(param,xyzdu,xyzdv);
204     if (orientation!=MEME_SENS)
205     {
206     xyzdu[0]=-xyzdu[0];
207     xyzdu[1]=-xyzdu[1];
208     xyzdu[2]=-xyzdu[2];
209     }
210     }
211    
212     void MG_FACE::deriver_seconde(double *uv,double* xyzduu,double* xyzduv,double* xyzdvv,double *xyz , double *xyzdu, double *xyzdv)
213     {
214     double param[2]={uv[0],uv[1]};
215     if (orientation!=MEME_SENS) param[0]=surface->get_umin()+surface->get_umax()-param[0];
216     surface->deriver_seconde(param,xyzduu,xyzduv,xyzdvv,xyz,xyzdu,xyzdv);
217     if (orientation!=MEME_SENS)
218     {
219     xyzdu[0]=-xyzdu[0];
220     xyzdu[1]=-xyzdu[1];
221     xyzdu[2]=-xyzdu[2];
222     xyzduv[0]=-xyzduv[0];
223     xyzduv[1]=-xyzduv[1];
224     xyzduv[2]=-xyzduv[2];
225     }
226     }
227    
228     void MG_FACE::inverser(double *uv,double *xyz,double precision)
229     {
230     surface->inverser(uv,xyz,precision);
231     if (orientation!=MEME_SENS) uv[0]=surface->get_umin()+surface->get_umax()-uv[0];
232     }
233    
234     void MG_FACE::calcul_normale(double *uv,double *normale)
235     {
236     double xyzdu[3];
237     double xyzdv[3];
238    
239     deriver(uv,xyzdu,xyzdv);
240     OT_VECTEUR_3D xu(xyzdu);
241     OT_VECTEUR_3D xv(xyzdv);
242     OT_VECTEUR_3D n=xu&xv;
243     normale[0]=n.get_x();
244     normale[1]=n.get_y();
245     normale[2]=n.get_z();
246     }
247    
248     void MG_FACE::calcul_normale_unitaire(double *uv,double *normale)
249     {
250     double xyzdu[3];
251     double xyzdv[3];
252    
253     deriver(uv,xyzdu,xyzdv);
254     OT_VECTEUR_3D xu(xyzdu);
255     OT_VECTEUR_3D xv(xyzdv);
256     OT_VECTEUR_3D n=xu&xv;
257     n.norme();
258     normale[0]=n.get_x();
259     normale[1]=n.get_y();
260     normale[2]=n.get_z();
261     }
262    
263    
264     void MG_FACE::get_EFG(double *uv,double& E,double& F,double& G)
265     {
266     double xyzdu[3];
267     double xyzdv[3];
268     deriver(uv,xyzdu,xyzdv);
269     E=xyzdu[0]*xyzdu[0]+xyzdu[1]*xyzdu[1]+xyzdu[2]*xyzdu[2];
270     F=xyzdu[0]*xyzdv[0]+xyzdu[1]*xyzdv[1]+xyzdu[2]*xyzdv[2];
271     G=xyzdv[0]*xyzdv[0]+xyzdv[1]*xyzdv[1]+xyzdv[2]*xyzdv[2];
272     }
273    
274     void MG_FACE::get_M(double *uv,double& M1,double& M2,double& M3)
275     {
276     double E,F,G;
277     double xyz[3],xyzdu[3],xyzdv[3],xyzduu[3],xyzduv[3],xyzdvv[3];
278    
279     deriver_seconde(uv,xyzduu,xyzduv,xyzdvv,xyz,xyzdu,xyzdv);
280     E=xyzdu[0]*xyzdu[0]+xyzdu[1]*xyzdu[1]+xyzdu[2]*xyzdu[2];
281     G=xyzdv[0]*xyzdv[0]+xyzdv[1]*xyzdv[1]+xyzdv[2]*xyzdv[2];
282     double Edu=2.*(xyzdu[0]*xyzduu[0]+xyzdu[1]*xyzduu[1]+xyzdu[2]*xyzduu[2]);
283     double Gdv=2.*(xyzdv[0]*xyzdvv[0]+xyzdv[1]*xyzdvv[1]+xyzdv[2]*xyzdvv[2]);
284     double Edv=2.*(xyzdu[0]*xyzduv[0]+xyzdu[1]*xyzduv[1]+xyzdu[2]*xyzduv[2]);
285     double Gdu=2.*(xyzdv[0]*xyzduv[0]+xyzdv[1]*xyzduv[1]+xyzdv[2]*xyzduv[2]);
286     double m1[3],m2[3],m3[3];
287     m1[0]=xyzduu[0]/E-0.5*Edu*xyzdu[0]/E/E;
288     m1[1]=xyzduu[1]/E-0.5*Edu*xyzdu[1]/E/E;
289     m1[2]=xyzduu[2]/E-0.5*Edu*xyzdu[2]/E/E;
290     m2[0]=xyzduv[0]/sqrt(E*G)-0.5*xyzdu[0]*Edv/E/sqrt(E*G)-0.5*xyzdv[0]*Gdu/G/sqrt(E*G);
291     m2[1]=xyzduv[1]/sqrt(E*G)-0.5*xyzdu[1]*Edv/E/sqrt(E*G)-0.5*xyzdv[1]*Gdu/G/sqrt(E*G);
292     m2[2]=xyzduv[2]/sqrt(E*G)-0.5*xyzdu[2]*Edv/E/sqrt(E*G)-0.5*xyzdv[2]*Gdu/G/sqrt(E*G);
293     m3[0]=xyzdvv[0]/G-0.5*Gdv*xyzdv[0]/G/G;
294     m3[1]=xyzdvv[1]/G-0.5*Gdv*xyzdv[1]/G/G;
295     m3[2]=xyzdvv[2]/G-0.5*Gdv*xyzdv[2]/G/G;
296     M1=sqrt(m1[0]*m1[0]+m1[1]*m1[1]+m1[2]*m1[2]);
297     M2=sqrt(m2[0]*m2[0]+m2[1]*m2[1]+m2[2]*m2[2]);
298     M3=sqrt(m3[0]*m3[0]+m3[1]*m3[1]+m3[2]*m3[2]);
299     }
300    
301     void MG_FACE::get_LMN(double *uv,double& L,double& M,double &N)
302     {
303    
304     double xyz[3],xyzdu[3],xyzdv[3],xyzduu[3],xyzduv[3],xyzdvv[3],normal[3];
305     deriver_seconde(uv,xyzduu,xyzduv,xyzdvv,xyz,xyzdu,xyzdv);
306     calcul_normale_unitaire(uv,normal);
307     L=xyzduu[0]*normal[0]+xyzduu[1]*normal[1]+xyzduu[2]*normal[2];
308     M=xyzduv[0]*normal[0]+xyzduv[1]*normal[1]+xyzduv[2]*normal[2];
309     N=xyzdvv[0]*normal[0]+xyzdvv[1]*normal[1]+xyzdvv[2]*normal[2];
310     }
311    
312    
313     void MG_FACE::get_courbure(double *uv,double& cmax,double& cmin)
314     {
315     double E,F,G,L,M,N;
316     get_EFG(uv,E,F,G);
317     get_LMN(uv,L,M,N);
318     double a=E*G-F*F;
319     double b=-E*N-G*L+2*F*M;
320     double c=L*N-M*M;
321     double delta=b*b-4*a*c;
322     if (delta<0.00000001) delta=0.;
323     double x1=(-b+sqrt(delta))/2./a;
324     double x2=(-b-sqrt(delta))/2./a;
325     if (fabs(x1)>fabs(x2)) {
326     cmax=x1;
327     cmin=x2;
328     } else {
329     cmax=x2;
330     cmin=x1;
331     }
332     }
333    
334 francois 632
335    
336 francois 1095 BOITE_3D MG_FACE::get_boite_3D(void)
337 francois 632 {
338     double xmin=1e308,ymin=1e308,zmin=1e308;
339     double xmax=-1e308,ymax=-1e308,zmax=-1e308;
340     double umin=1e308,vmin=1e308;
341     double umax=-1e308,vmax=-1e308;
342     int nb_boucle=get_nb_mg_boucle();
343 francois 824 if (surface->est_periodique_u())
344     {
345     umin=0.;
346     umax=surface->get_periode_u();
347     }
348     if (surface->est_periodique_v())
349     {
350     vmin=0.;
351     vmax=surface->get_periode_v();
352     }
353 francois 632 for (int i=0;i<nb_boucle;i++)
354     {
355     MG_BOUCLE* bou=get_mg_boucle(i);
356     int nb_arete=bou->get_nb_mg_coarete();
357     for (int j=0;j<nb_arete;j++)
358     {
359     MG_ARETE* arete=bou->get_mg_coarete(j)->get_arete();
360     double tmin=arete->get_tmin();
361     double tmax=arete->get_tmax();
362     double tdemi=0.5*(tmax+tmin);
363     double xyz1[3],xyz2[3];
364     for (int k=0;k<pas_echantillon+1;k++)
365     {
366     double t=tmin+1.0*k/pas_echantillon*(tmax-tmin);
367     double xyz[3];
368     arete->evaluer(t,xyz);
369     if (xyz[0]<xmin) xmin=xyz[0];
370     if (xyz[0]>xmax) xmax=xyz[0];
371     if (xyz[1]<ymin) ymin=xyz[1];
372     if (xyz[1]>ymax) ymax=xyz[1];
373     if (xyz[2]<zmin) zmin=xyz[2];
374     if (xyz[2]>zmax) zmax=xyz[2];
375     double uv[2];
376 francois 1095 inverser(uv,xyz);
377 francois 632 if (uv[0]>umax) umax=uv[0];
378     if (uv[0]<umin) umin=uv[0];
379     if (uv[1]>vmax) vmax=uv[1];
380     if (uv[1]<vmin) vmin=uv[1];
381     }
382     }
383 francois 1095 std::vector<double> lstpoint;
384     get_echantillonnage(pas_echantillon*pas_echantillon,lstpoint);
385     if (lstpoint.size()<1)
386     {
387     for (int k=0;k<pas_echantillon+1;k++)
388     for (int l=0;l<pas_echantillon+1;l++)
389 francois 632 {
390     double uv[2];
391     uv[0]=umin+1.0*k/pas_echantillon*(umax-umin);
392     uv[1]=vmin+1.0*l/pas_echantillon*(vmax-vmin);
393     double xyz[3];
394     evaluer(uv,xyz);
395     if (xyz[0]<xmin) xmin=xyz[0];
396     if (xyz[1]<ymin) ymin=xyz[1];
397     if (xyz[2]<zmin) zmin=xyz[2];
398     if (xyz[0]>xmax) xmax=xyz[0];
399     if (xyz[1]>ymax) ymax=xyz[1];
400     if (xyz[2]>zmax) zmax=xyz[2];
401     }
402 francois 1095 }
403     else
404     for (int k=0;k<lstpoint.size()/5;k++)
405     {
406     double xyz[3]={lstpoint[5*k],lstpoint[5*k+1],lstpoint[5*k+2]};
407     if (xyz[0]<xmin) xmin=xyz[0];
408     if (xyz[1]<ymin) ymin=xyz[1];
409     if (xyz[2]<zmin) zmin=xyz[2];
410     if (xyz[0]>xmax) xmax=xyz[0];
411     if (xyz[1]>ymax) ymax=xyz[1];
412     if (xyz[2]>zmax) zmax=xyz[2];
413     }
414 francois 632 }
415     BOITE_3D boite(xmin,ymin,zmin,xmax,ymax,zmax);
416     return boite;
417    
418     }
419    
420    
421 couturad 814 int MG_FACE::get_nb_pole(void)
422     {
423     return nb_pole;
424     }
425 francois 632
426    
427 francois 820 void MG_FACE::get_liste_pole_uv(std::vector<double> *liste_pole_uv,double eps)
428 couturad 814 {
429 francois 820 surface->get_liste_pole(liste_pole_uv,eps);
430     nb_pole = liste_pole_uv->size()/2;
431 couturad 814 }
432 francois 632
433 francois 1095 void MG_FACE::get_xyz_min_max(double* xyzmin, double* xyzmax)
434 couturad 919 {
435 francois 984 /*xyzmin[0] = std::numeric_limits<double>::max();
436 couturad 951 xyzmin[1] = std::numeric_limits<double>::max();
437     xyzmin[2] = std::numeric_limits<double>::max();
438     xyzmax[0] = std::numeric_limits<double>::min();
439     xyzmax[1] = std::numeric_limits<double>::min();
440     xyzmax[2] = std::numeric_limits<double>::min();
441 couturad 919 double umin = get_surface()->get_umin();
442     double umax = get_surface()->get_umax();
443     if(get_surface()->est_periodique_u())
444     {
445     umin=0.0;
446     umax=get_surface()->get_periode_u();
447     }
448     double pas_u = (umax-umin)/(double)nb_pas;
449     double vmin = get_surface()->get_vmin();
450     double vmax = get_surface()->get_vmax();
451     if(get_surface()->est_periodique_v())
452     {
453     vmin=0.0;
454     vmax=get_surface()->get_periode_v();
455     }
456 francois 984 double pas_v = (vmax-vmin)/(double)nb_pas;
457 couturad 919 double uv[2];
458     for(int i=0;i<nb_pas;i++)
459     {
460     uv[0] = umin+pas_u*i;
461     for(int j=0;j<nb_pas;j++)
462     {
463     uv[1] = vmin+j*pas_v;
464     double xyz[3];
465     evaluer(uv,xyz);
466     if(xyz[0]>xyzmax[0]) xyzmax[0]=xyz[0];
467     else if(xyz[0]<xyzmin[0]) xyzmin[0]=xyz[0];
468     if(xyz[1]>xyzmax[1]) xyzmax[1]=xyz[1];
469     else if(xyz[1]<xyzmin[1]) xyzmin[1]=xyz[1];
470     if(xyz[2]>xyzmax[2]) xyzmax[2]=xyz[2];
471     else if(xyz[2]<xyzmin[2]) xyzmin[2]=xyz[2];
472     }
473 francois 984 }*/
474 francois 1095 BOITE_3D boite=get_boite_3D();
475 francois 984 xyzmin[0] = boite.get_xmin();
476     xyzmin[1] = boite.get_ymin();
477     xyzmin[2] = boite.get_zmin();
478     xyzmax[0] = boite.get_xmax();
479     xyzmax[1] = boite.get_ymax();
480     xyzmax[2] = boite.get_zmax();
481    
482 couturad 919 }
483    
484 francois 1095 void MG_FACE::get_echantillonnage(int numechantillon, std::vector<double> &tab)
485     {
486     surface->get_echantillonnage(numechantillon,tab,epsilon_echantillon,angle_dev_echantillon);
487    
488     }
489 couturad 919
490 francois 1095
491 couturad 814 void MG_FACE::change_nb_pole(int val)
492     {
493     nb_pole=val;
494     }
495 francois 632
496 francois 283 VCT& MG_FACE::get_vectorisation(void)
497     {
498     if (vect==NULL) vect=new VCT_FACE(this);
499     return *vect;
500     }
501    
502 francois 763 void MG_FACE::enregistrer(std::ostream& o,double version)
503 francois 283 {
504     o << "%" << get_id() << "=FACE("<< get_idoriginal() << ",$" << surface->get_id() << ",(";
505     for (unsigned int i=0;i<lst_boucle.size();i++)
506     {
507     o << "$" << lst_boucle[i]->get_id();
508     if (i!=lst_boucle.size()-1) o << ",";
509     else o << ")";
510     }
511     int nb=get_nb_ccf();
512 couturad 814 o << "," << orientation << "," << nb_pole << ",";
513 francois 763 if (version<2)
514 francois 283 {
515 francois 763 o << nb;
516     if (nb!=0)
517     {
518 francois 283 o << ",(";
519     for (int i=0;i<nb;i++)
520     {
521     char nom[3];
522     get_type_ccf(i,nom);
523     o << "(" << nom << "," << get_valeur_ccf(i) << ")";
524     if (i!=nb-1) o << "," ;
525     }
526     o << ")";
527 francois 763
528     }
529 francois 283 }
530 francois 763 else
531     enregistrer_ccf(o,version);
532 francois 283 o << ");" << std::endl;
533     }
534    
535