ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/geometrie/src/fem_tetra10.cpp
Revision: 309
Committed: Tue Jan 24 15:25:16 2012 UTC (13 years, 3 months ago) by francois
File size: 11972 byte(s)
Log Message: 
Preparation du maillage FEM pour inclure de nouveaux élements non necessairement triangle/tetra

File Contents

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