ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/geometrie/src/fem_hexa8.cpp
Revision: 757
Committed: Fri Nov 6 13:42:55 2015 UTC (9 years, 6 months ago) by francois
File size: 15429 byte(s)
Log Message:
Possibilite de calculer les comtraintes et les deformations aux points de gauss
Tester en elasticite volumique

File Contents

# User Rev Content
1 francois 310 //------------------------------------------------------------
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     // fem_triangle3.cpp
16     //
17     //------------------------------------------------------------
18     //------------------------------------------------------------
19     // COPYRIGHT 2000
20     // Version du 02/03/2006 � 11H22
21     //------------------------------------------------------------
22     //------------------------------------------------------------
23    
24    
25     #include "gestionversion.h"
26     #include "fem_hexa8.h"
27     #include "fem_maillage.h"
28     #include "fem_noeud.h"
29     #include "mg_element_maillage.h"
30     #include "ot_tenseur.h"
31 francois 754 #include "ot_quadrature_gauss.h"
32 francois 310
33    
34     FEM_HEXA8::FEM_HEXA8(unsigned long num,class MG_ELEMENT_MAILLAGE* mai,class FEM_NOEUD** tabnoeud):FEM_ELEMENT3(num,mai),FEM_TEMPLATE_ELEMENT<8>(tabnoeud)
35     {
36 francois 581 if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==3) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
37 francois 310 tab[0]->get_lien_element3()->ajouter(this);
38     tab[1]->get_lien_element3()->ajouter(this);
39     tab[2]->get_lien_element3()->ajouter(this);
40     tab[3]->get_lien_element3()->ajouter(this);
41     tab[4]->get_lien_element3()->ajouter(this);
42     tab[5]->get_lien_element3()->ajouter(this);
43     tab[6]->get_lien_element3()->ajouter(this);
44     tab[7]->get_lien_element3()->ajouter(this);
45     get_fem_noeudpetitid()->get_lien_petit_element3()->ajouter(this);
46     }
47    
48     FEM_HEXA8::FEM_HEXA8(class MG_ELEMENT_MAILLAGE* mai,FEM_NOEUD** tabnoeud):FEM_ELEMENT3(mai),FEM_TEMPLATE_ELEMENT<8>(tabnoeud)
49     {
50 francois 581 if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==3) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
51 francois 310 tab[0]->get_lien_element3()->ajouter(this);
52     tab[1]->get_lien_element3()->ajouter(this);
53     tab[2]->get_lien_element3()->ajouter(this);
54     tab[3]->get_lien_element3()->ajouter(this);
55     tab[4]->get_lien_element3()->ajouter(this);
56     tab[5]->get_lien_element3()->ajouter(this);
57     tab[6]->get_lien_element3()->ajouter(this);
58     tab[7]->get_lien_element3()->ajouter(this);
59     get_fem_noeudpetitid()->get_lien_petit_element3()->ajouter(this);
60     }
61    
62 francois 378 FEM_HEXA8::FEM_HEXA8(unsigned long num,class MG_ELEMENT_TOPOLOGIQUE* topo,class FEM_NOEUD** tabnoeud):FEM_ELEMENT3(num,topo),FEM_TEMPLATE_ELEMENT<8>(tabnoeud)
63     {
64 francois 581 if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==3) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
65 francois 378 tab[0]->get_lien_element3()->ajouter(this);
66     tab[1]->get_lien_element3()->ajouter(this);
67     tab[2]->get_lien_element3()->ajouter(this);
68     tab[3]->get_lien_element3()->ajouter(this);
69     tab[4]->get_lien_element3()->ajouter(this);
70     tab[5]->get_lien_element3()->ajouter(this);
71     tab[6]->get_lien_element3()->ajouter(this);
72     tab[7]->get_lien_element3()->ajouter(this);
73     get_fem_noeudpetitid()->get_lien_petit_element3()->ajouter(this);
74     }
75    
76     FEM_HEXA8::FEM_HEXA8(class MG_ELEMENT_TOPOLOGIQUE* topo,FEM_NOEUD** tabnoeud):FEM_ELEMENT3(topo),FEM_TEMPLATE_ELEMENT<8>(tabnoeud)
77     {
78 francois 581 if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==3) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
79 francois 378 tab[0]->get_lien_element3()->ajouter(this);
80     tab[1]->get_lien_element3()->ajouter(this);
81     tab[2]->get_lien_element3()->ajouter(this);
82     tab[3]->get_lien_element3()->ajouter(this);
83     tab[4]->get_lien_element3()->ajouter(this);
84     tab[5]->get_lien_element3()->ajouter(this);
85     tab[6]->get_lien_element3()->ajouter(this);
86     tab[7]->get_lien_element3()->ajouter(this);
87     get_fem_noeudpetitid()->get_lien_petit_element3()->ajouter(this);
88     }
89    
90     FEM_HEXA8::FEM_HEXA8(unsigned long num,class MG_ELEMENT_TOPOLOGIQUE* topo,class MG_ELEMENT_MAILLAGE* mai,class FEM_NOEUD** tabnoeud):FEM_ELEMENT3(num,topo,mai),FEM_TEMPLATE_ELEMENT<8>(tabnoeud)
91     {
92 francois 581 if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==3) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
93 francois 378 tab[0]->get_lien_element3()->ajouter(this);
94     tab[1]->get_lien_element3()->ajouter(this);
95     tab[2]->get_lien_element3()->ajouter(this);
96     tab[3]->get_lien_element3()->ajouter(this);
97     tab[4]->get_lien_element3()->ajouter(this);
98     tab[5]->get_lien_element3()->ajouter(this);
99     tab[6]->get_lien_element3()->ajouter(this);
100     tab[7]->get_lien_element3()->ajouter(this);
101     get_fem_noeudpetitid()->get_lien_petit_element3()->ajouter(this);
102     }
103    
104     FEM_HEXA8::FEM_HEXA8(class MG_ELEMENT_TOPOLOGIQUE* topo,class MG_ELEMENT_MAILLAGE* mai,FEM_NOEUD** tabnoeud):FEM_ELEMENT3(topo,mai),FEM_TEMPLATE_ELEMENT<8>(tabnoeud)
105     {
106 francois 581 if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==3) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
107 francois 378 tab[0]->get_lien_element3()->ajouter(this);
108     tab[1]->get_lien_element3()->ajouter(this);
109     tab[2]->get_lien_element3()->ajouter(this);
110     tab[3]->get_lien_element3()->ajouter(this);
111     tab[4]->get_lien_element3()->ajouter(this);
112     tab[5]->get_lien_element3()->ajouter(this);
113     tab[6]->get_lien_element3()->ajouter(this);
114     tab[7]->get_lien_element3()->ajouter(this);
115     get_fem_noeudpetitid()->get_lien_petit_element3()->ajouter(this);
116     }
117    
118    
119    
120    
121 francois 310 FEM_HEXA8::FEM_HEXA8(FEM_HEXA8& mdd):FEM_ELEMENT3(mdd),FEM_TEMPLATE_ELEMENT<8>(mdd)
122     {
123 francois 581 if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==3) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
124 francois 310 tab[0]->get_lien_element3()->ajouter(this);
125     tab[1]->get_lien_element3()->ajouter(this);
126     tab[2]->get_lien_element3()->ajouter(this);
127     tab[3]->get_lien_element3()->ajouter(this);
128     tab[4]->get_lien_element3()->ajouter(this);
129     tab[5]->get_lien_element3()->ajouter(this);
130     tab[6]->get_lien_element3()->ajouter(this);
131     tab[7]->get_lien_element3()->ajouter(this);
132     get_fem_noeudpetitid()->get_lien_petit_element3()->ajouter(this);
133     }
134     FEM_HEXA8::~FEM_HEXA8()
135     {
136 francois 663 if (liaison_topologique!=NULL)
137     if (liaison_topologique->get_dimension()==3) liaison_topologique->get_lien_fem_maillage()->supprimer(this);
138 francois 310 tab[0]->get_lien_element3()->supprimer(this);
139     tab[1]->get_lien_element3()->supprimer(this);
140     tab[2]->get_lien_element3()->supprimer(this);
141     tab[3]->get_lien_element3()->supprimer(this);
142     tab[4]->get_lien_element3()->supprimer(this);
143     tab[5]->get_lien_element3()->supprimer(this);
144     tab[6]->get_lien_element3()->supprimer(this);
145     tab[7]->get_lien_element3()->supprimer(this);
146     get_fem_noeudpetitid()->get_lien_petit_element3()->supprimer(this);
147     }
148    
149    
150     FEM_ELEMENT_MAILLAGE* FEM_HEXA8::dupliquer(FEM_MAILLAGE *femmai,long decalage)
151     {
152     FEM_NOEUD* tabnoeud[8];
153     tabnoeud[0]=femmai->get_fem_noeudid(tab[0]->get_id()+decalage);
154     tabnoeud[1]=femmai->get_fem_noeudid(tab[1]->get_id()+decalage);
155     tabnoeud[2]=femmai->get_fem_noeudid(tab[2]->get_id()+decalage);
156     tabnoeud[3]=femmai->get_fem_noeudid(tab[3]->get_id()+decalage);
157     tabnoeud[4]=femmai->get_fem_noeudid(tab[4]->get_id()+decalage);
158     tabnoeud[5]=femmai->get_fem_noeudid(tab[5]->get_id()+decalage);
159     tabnoeud[6]=femmai->get_fem_noeudid(tab[6]->get_id()+decalage);
160     tabnoeud[7]=femmai->get_fem_noeudid(tab[7]->get_id()+decalage);
161     FEM_HEXA8* quad=new FEM_HEXA8(get_id()+decalage,maillage,tabnoeud);
162     femmai->ajouter_fem_element3(quad);
163     return quad;
164     }
165    
166    
167 francois 635 bool FEM_HEXA8::valide_parametre_EF(double* uvw)
168     {
169 francois 674 if (uvw[0]>=-1.-1e-10)
170     if (uvw[1]>=-1.-1e-10)
171     if (uvw[2]>=-1.-1e-10)
172 francois 635 if (uvw[0]<=1.+1e-10)
173     if (uvw[1]<=1.+1e-10)
174     if (uvw[2]<=1.+1e-10)
175     return true;
176     return false;
177     }
178 francois 310
179 francois 684 void FEM_HEXA8::reinit_boite_3D(void)
180     {
181     FEM_TEMPLATE_ELEMENT<8>::reinit_boite_3D();
182     }
183 francois 310
184     int FEM_HEXA8::get_type_entite(void)
185     {
186     return IDFEM_HEXA8;
187     }
188    
189     int FEM_HEXA8::get_dimension(void)
190     {
191     return 3;
192     }
193    
194    
195    
196     int FEM_HEXA8::get_nb_fem_noeud(void)
197     {
198     return FEM_TEMPLATE_ELEMENT<8>::get_nb_fem_noeud();
199     }
200    
201     FEM_NOEUD* FEM_HEXA8::get_fem_noeud(int num)
202     {
203     return FEM_TEMPLATE_ELEMENT<8>::get_fem_noeud(num);
204     }
205    
206     void FEM_HEXA8::change_noeud(int num,FEM_NOEUD* noeud)
207     {
208     FEM_TEMPLATE_ELEMENT<8>::change_noeud(num,noeud);
209     }
210    
211     BOITE_3D& FEM_HEXA8::get_boite_3D(void)
212     {
213     return FEM_TEMPLATE_ELEMENT<8>::get_boite_3D();
214     }
215    
216    
217     void FEM_HEXA8::enregistrer(std::ostream& o)
218     {
219 francois 378 if (maillage!=NULL)
220     if (get_lien_topologie()!=NULL) o << "%" << get_id() << "=FEM_HEXA8($"<< get_lien_topologie()->get_id() << ",$" << maillage->get_id() << ",$" << tab[0]->get_id() << ",$" << tab[1]->get_id() << ",$" << tab[2]->get_id()<< ",$" << tab[3]->get_id() << ",$" << tab[4]->get_id() << ",$" << tab[5]->get_id() << ",$" << tab[6]->get_id()<< ",$" << tab[7]->get_id()<< ");" << std::endl;
221     else o << "%" << get_id() << "=FEM_HEXA8(NULL,$" << maillage->get_id() << ",$" << tab[0]->get_id() << ",$" << tab[1]->get_id() << ",$" << tab[2]->get_id()<< ",$" << tab[3]->get_id() << ",$" << tab[4]->get_id() << ",$" << tab[5]->get_id() << ",$" << tab[6]->get_id()<< ",$" << tab[7]->get_id()<<");" << std::endl;
222     else
223     if (get_lien_topologie()!=NULL) o << "%" << get_id() << "=FEM_HEXA8($"<< get_lien_topologie()->get_id() << ",NULL,$" << tab[0]->get_id() << ",$" << tab[1]->get_id() << ",$" << tab[2]->get_id()<< ",$" << tab[3]->get_id() << ",$" << tab[4]->get_id() << ",$" << tab[5]->get_id() << ",$" << tab[6]->get_id()<< ",$" << tab[7]->get_id()<< ");" << std::endl;
224     else o << "%" << get_id() << "=FEM_HEXA8(NULL,NULL,$" << tab[0]->get_id() << ",$" << tab[1]->get_id() << ",$" << tab[2]->get_id()<< ",$" << tab[3]->get_id() << ",$" << tab[4]->get_id() << ",$" << tab[5]->get_id() << ",$" << tab[6]->get_id()<< ",$" << tab[7]->get_id()<<");" << std::endl;
225    
226 francois 310 }
227 francois 637
228    
229     int FEM_HEXA8::get_nb_pt_gauss(int degre)
230     {
231 francois 754 return OT_POINTS_GAUSS::get_nb_point_hexa_prod(degre);
232 francois 637 }
233    
234     void FEM_HEXA8::get_pt_gauss(int degre,int num,double &w,double *uvw)
235     {
236 francois 754 return OT_POINTS_GAUSS::get_pt_gauss_hex_prod(degre,num,w,uvw);
237 francois 637 }
238 francois 754
239 francois 757 int FEM_HEXA8::get_degre_gauss(int num)
240     {
241     return OT_POINTS_GAUSS::get_degre_gauss_hexa_prod(num);
242     }
243    
244 francois 638 int FEM_HEXA8::get_nb_fonction_interpolation(void)
245 francois 310 {
246     return 8;
247     }
248    
249 francois 638 int FEM_HEXA8::get_degremax_fonction_interpolation(void)
250     {
251     return 3;
252     }
253 francois 310 double FEM_HEXA8::get_fonction_interpolation(int num,double *uv)
254     {
255     double val;
256     double a1=1.+uv[0],a2=1.-uv[0];
257     double b1=1.+uv[1],b2=1.-uv[1];
258     double c1=1.+uv[2],c2=1.-uv[2];
259     switch (num)
260     {
261     case 1:
262     val=0.125*a2*b2*c2;
263     break;
264     case 2:
265     val=0.125*a1*b2*c2;
266     break;
267     case 3:
268     val=0.125*a1*b1*c2;
269     break;
270     case 4:
271     val=0.125*a2*b1*c2;
272     break;
273     case 5:
274     val=0.125*a2*b2*c1;
275     break;
276     case 6:
277     val=0.125*a1*b2*c1;
278     break;
279     case 7:
280     val=0.125*a1*b1*c1;
281     break;
282     case 8:
283     val=0.125*a2*b1*c1;
284     break;
285     }
286     return val;
287     }
288    
289     double FEM_HEXA8::get_fonction_derive_interpolation(int num,int num_variable,double *uv)
290     {
291     double val;
292     double a1=1.+uv[0],a2=1.-uv[0];
293     double b1=1.+uv[1],b2=1.-uv[1];
294     double c1=1.+uv[2],c2=1.-uv[2];
295     switch (num)
296     {
297     case 1:
298     switch (num_variable)
299     {
300     case 1:
301     val=-0.125*b2*c2;
302     break;
303     case 2:
304     val=-0.125*a2*c2;
305     break;
306     case 3:
307     val=-0.125*a2*b2;
308     break;
309     } break;
310     case 2:
311     switch (num_variable)
312     {
313     case 1:
314     val=0.125*b2*c2;
315     break;
316     case 2:
317     val=-0.125*a1*c2;
318     break;
319     case 3:
320     val=-0.125*a1*b2;
321     break;
322     } break;
323     case 3:
324     switch (num_variable)
325     {
326     case 1:
327     val=0.125*b1*c2;
328     break;
329     case 2:
330     val=0.125*a1*c2;
331     break;
332     case 3:
333     val=-0.125*a1*b1;
334     break;
335     } break;
336     case 4:
337     switch (num_variable)
338     {
339     case 1:
340     val=-0.125*b1*c2;
341     break;
342     case 2:
343     val=0.125*a2*c2;
344     break;
345     case 3:
346     val=-0.125*a2*b1;
347     break;
348     } break;
349     case 5:
350     switch (num_variable)
351     {
352     case 1:
353     val=-0.125*b2*c1;
354     break;
355     case 2:
356     val=-0.125*a2*c1;
357     break;
358     case 3:
359     val=0.125*a2*b2;
360     break;
361     } break;
362     case 6:
363     switch (num_variable)
364     {
365     case 1:
366     val=0.125*b2*c1;
367     break;
368     case 2:
369     val=-0.125*a1*c1;
370     break;
371     case 3:
372     val=0.125*a1*b2;
373     break;
374     } break;
375     case 7:
376     switch (num_variable)
377     {
378     case 1:
379     val=0.125*b1*c1;
380     break;
381     case 2:
382     val=0.125*a1*c1;
383     break;
384     case 3:
385     val=0.125*a1*b1;
386     break;
387     } break;
388     case 8:
389     switch (num_variable)
390     {
391     case 1:
392     val=-0.125*b1*c1;
393     break;
394     case 2:
395     val=0.125*a2*c1;
396     break;
397     case 3:
398     val=0.125*a2*b1;
399     break;
400     } break;
401     }
402     return val;
403     }
404    
405     double FEM_HEXA8::get_jacobien(double* jac,double *uv,int& li,int& col,double unite)
406     {
407     li=3;col=3;
408     OT_TENSEUR Ni(3,8),Xi(8,3);
409     for (int i=0;i<3;i++)
410     for (int j=0;j<8;j++)
411     Ni(i,j)=get_fonction_derive_interpolation(j+1,i+1,uv);
412     for (int i=0;i<8;i++)
413     {
414     double *xyz=tab[i]->get_coord();
415     Xi(i,0)=xyz[0]*unite;
416     Xi(i,1)=xyz[1]*unite;
417     Xi(i,2)=xyz[2]*unite;
418     }
419     OT_TENSEUR jacobi=Ni*Xi;
420     jac[0]=jacobi(0,0).get_x();
421     jac[1]=jacobi(0,1).get_x();
422     jac[2]=jacobi(0,2).get_x();
423     jac[3]=jacobi(1,0).get_x();
424     jac[4]=jacobi(1,1).get_x();
425     jac[5]=jacobi(1,2).get_x();
426     jac[6]=jacobi(2,0).get_x();
427     jac[7]=jacobi(2,1).get_x();
428     jac[8]=jacobi(2,2).get_x();
429     double det=jac[0]*jac[4]*jac[8]+jac[3]*jac[7]*jac[2]+jac[6]*jac[1]*jac[5]-jac[2]*jac[4]*jac[6]-jac[5]*jac[7]*jac[0]-jac[8]*jac[1]*jac[3];
430     return det;
431    
432     }
433    
434     void FEM_HEXA8::get_inverse_jacob(double* j,double *uv,double unite)
435     {
436     double jac[9];
437     int li,col;
438     double detj=get_jacobien(jac,uv,li,col,unite);
439     j[0*3+0] = (jac[1*3+1]*jac[2*3+2]-jac[1*3+2]*jac[2*3+1])/(detj);
440     j[0*3+1] = -(jac[0*3+1]*jac[2*3+2]-jac[0*3+2]*jac[2*3+1])/(detj);
441     j[0*3+2] =-(-jac[0*3+1]*jac[1*3+2]+jac[0*3+2]*jac[1*3+1])/(detj);
442    
443     j[1*3+0] = -(jac[1*3+0]*jac[2*3+2]-jac[1*3+2]*jac[2*3+0])/(detj);
444     j[1*3+1] = (jac[0*3+0]*jac[2*3+2]-jac[0*3+2]*jac[2*3+0])/(detj);
445     j[1*3+2] = -(jac[0*3+0]*jac[1*3+2]-jac[0*3+2]*jac[1*3+0])/(detj);
446    
447     j[2*3+0] =-(-jac[1*3+0]*jac[2*3+1]+jac[1*3+1]*jac[2*3+0])/(detj);
448     j[2*3+1] = -(jac[0*3+0]*jac[2*3+1]-jac[0*3+1]*jac[2*3+0])/(detj);
449     j[2*3+2] = (jac[0*3+0]*jac[1*3+1]-jac[0*3+1]*jac[1*3+0])/(detj);
450    
451 francois 405 }
452    
453    
454    
455 francois 406 int FEM_HEXA8::verifie_validite_decoupage_xfem(double *vol)
456 francois 405 {
457    
458 francois 310 }