ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/geometrie/src/fem_hexa8.cpp
Revision: 637
Committed: Mon Jan 26 21:26:37 2015 UTC (10 years, 7 months ago) by francois
File size: 16858 byte(s)
Log Message:
entree des points de gauss pour les différents éléments finis depuis reference JCC

File Contents

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