ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/geometrie/src/fem_triangle3.cpp
Revision: 786
Committed: Wed Mar 9 21:10:49 2016 UTC (9 years, 2 months ago) by francois
File size: 11028 byte(s)
Log Message: 
Possibilité d'obtenir la matrice de rigidite de code aster sous form de matrice bande.
Pour cela ajouter -opt dans l'operation -fem à la création du maillage FEM et dans -calculaster utiliser les nouveaux types d'étude rigidite.
Attention pour cela il faut une version de code aster modifiée
De manière automatique elle peut etre optenu en faisant ./instal_aster 12.5mod
Uniquement implantée avec la version 12.5 de code aster.

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     // 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_triangle3.h"
27     #include "fem_maillage.h"
28     #include "fem_noeud.h"
29     #include "mg_element_maillage.h"
30 francois 754 #include "ot_quadrature_gauss.h"
31 francois 283
32    
33 francois 309 FEM_TRIANGLE3::FEM_TRIANGLE3(unsigned long num,class MG_ELEMENT_MAILLAGE* mai,class FEM_NOEUD** tabnoeud):FEM_ELEMENT2(num,mai),FEM_TEMPLATE_ELEMENT<3>(tabnoeud)
34 francois 283 {
35     if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==2) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
36 francois 309 tab[0]->get_lien_element2()->ajouter(this);
37     tab[1]->get_lien_element2()->ajouter(this);
38     tab[2]->get_lien_element2()->ajouter(this);
39     get_fem_noeudpetitid()->get_lien_petit_element2()->ajouter(this);
40 francois 283 }
41    
42 francois 309 FEM_TRIANGLE3::FEM_TRIANGLE3(class MG_ELEMENT_MAILLAGE* mai,FEM_NOEUD** tabnoeud):FEM_ELEMENT2(mai),FEM_TEMPLATE_ELEMENT<3>(tabnoeud)
43 francois 283 {
44     if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==2) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
45 francois 309 tab[0]->get_lien_element2()->ajouter(this);
46     tab[1]->get_lien_element2()->ajouter(this);
47     tab[2]->get_lien_element2()->ajouter(this);
48     get_fem_noeudpetitid()->get_lien_petit_element2()->ajouter(this);
49 francois 283 }
50 francois 378 FEM_TRIANGLE3::FEM_TRIANGLE3(unsigned long num,class MG_ELEMENT_TOPOLOGIQUE* topo,class FEM_NOEUD** tabnoeud):FEM_ELEMENT2(num,topo),FEM_TEMPLATE_ELEMENT<3>(tabnoeud)
51     {
52     if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==2) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
53     tab[0]->get_lien_element2()->ajouter(this);
54     tab[1]->get_lien_element2()->ajouter(this);
55     tab[2]->get_lien_element2()->ajouter(this);
56     get_fem_noeudpetitid()->get_lien_petit_element2()->ajouter(this);
57     }
58 francois 283
59 francois 378 FEM_TRIANGLE3::FEM_TRIANGLE3(class MG_ELEMENT_TOPOLOGIQUE* topo,FEM_NOEUD** tabnoeud):FEM_ELEMENT2(topo),FEM_TEMPLATE_ELEMENT<3>(tabnoeud)
60     {
61     if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==2) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
62     tab[0]->get_lien_element2()->ajouter(this);
63     tab[1]->get_lien_element2()->ajouter(this);
64     tab[2]->get_lien_element2()->ajouter(this);
65     get_fem_noeudpetitid()->get_lien_petit_element2()->ajouter(this);
66     }
67     FEM_TRIANGLE3::FEM_TRIANGLE3(unsigned long num,class MG_ELEMENT_TOPOLOGIQUE* topo,class MG_ELEMENT_MAILLAGE* mai,class FEM_NOEUD** tabnoeud):FEM_ELEMENT2(num,topo,mai),FEM_TEMPLATE_ELEMENT<3>(tabnoeud)
68     {
69     if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==2) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
70     tab[0]->get_lien_element2()->ajouter(this);
71     tab[1]->get_lien_element2()->ajouter(this);
72     tab[2]->get_lien_element2()->ajouter(this);
73     get_fem_noeudpetitid()->get_lien_petit_element2()->ajouter(this);
74     }
75    
76     FEM_TRIANGLE3::FEM_TRIANGLE3(class MG_ELEMENT_TOPOLOGIQUE* topo,class MG_ELEMENT_MAILLAGE* mai,FEM_NOEUD** tabnoeud):FEM_ELEMENT2(topo,mai),FEM_TEMPLATE_ELEMENT<3>(tabnoeud)
77     {
78     if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==2) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
79     tab[0]->get_lien_element2()->ajouter(this);
80     tab[1]->get_lien_element2()->ajouter(this);
81     tab[2]->get_lien_element2()->ajouter(this);
82     get_fem_noeudpetitid()->get_lien_petit_element2()->ajouter(this);
83     }
84 francois 309 FEM_TRIANGLE3::FEM_TRIANGLE3(FEM_TRIANGLE3& mdd):FEM_ELEMENT2(mdd),FEM_TEMPLATE_ELEMENT<3>(mdd)
85 francois 283 {
86     if (liaison_topologique!=NULL) if (liaison_topologique->get_dimension()==2) liaison_topologique->get_lien_fem_maillage()->ajouter(this);
87 francois 309 tab[0]->get_lien_element2()->ajouter(this);
88     tab[1]->get_lien_element2()->ajouter(this);
89     tab[2]->get_lien_element2()->ajouter(this);
90     get_fem_noeudpetitid()->get_lien_petit_element2()->ajouter(this);
91 francois 283 }
92     FEM_TRIANGLE3::~FEM_TRIANGLE3()
93     {
94 francois 663 if (liaison_topologique!=NULL)
95     if (liaison_topologique->get_dimension()==0) liaison_topologique->get_lien_fem_maillage()->supprimer(this);
96 francois 309 tab[0]->get_lien_element2()->supprimer(this);
97     tab[1]->get_lien_element2()->supprimer(this);
98     tab[2]->get_lien_element2()->supprimer(this);
99     get_fem_noeudpetitid()->get_lien_petit_element2()->supprimer(this);
100 francois 283 }
101    
102    
103     FEM_ELEMENT_MAILLAGE* FEM_TRIANGLE3::dupliquer(FEM_MAILLAGE *femmai,long decalage)
104     {
105     FEM_NOEUD* tabnoeud[3];
106     tabnoeud[0]=femmai->get_fem_noeudid(tab[0]->get_id()+decalage);
107     tabnoeud[1]=femmai->get_fem_noeudid(tab[1]->get_id()+decalage);
108     tabnoeud[2]=femmai->get_fem_noeudid(tab[2]->get_id()+decalage);
109     FEM_TRIANGLE3* tri=new FEM_TRIANGLE3(get_id()+decalage,maillage,tabnoeud);
110 francois 309 femmai->ajouter_fem_element2(tri);
111 francois 283 return tri;
112     }
113    
114 francois 786 void FEM_TRIANGLE3::get_voisin_noeud(class FEM_NOEUD* no,TPL_LISTE_ENTITE<FEM_NOEUD*> &voisin)
115     {
116     voisin.vide();
117     if (no==tab[0])
118     {
119     voisin.ajouter(tab[1]);
120     voisin.ajouter(tab[2]);
121     }
122     if (no==tab[1])
123     {
124     voisin.ajouter(tab[0]);
125     voisin.ajouter(tab[2]);
126     }
127     if (no==tab[2])
128     {
129     voisin.ajouter(tab[0]);
130     voisin.ajouter(tab[1]);
131     }
132     }
133 francois 283
134 francois 684 void FEM_TRIANGLE3::reinit_boite_3D(void)
135     {
136     FEM_TEMPLATE_ELEMENT<3>::reinit_boite_3D();
137     }
138 francois 283
139    
140     int FEM_TRIANGLE3::get_type_entite(void)
141     {
142     return IDFEM_TRIANGLE3;
143     }
144    
145     int FEM_TRIANGLE3::get_dimension(void)
146     {
147     return 2;
148     }
149    
150    
151     int FEM_TRIANGLE3::get_nb_fem_noeud(void)
152     {
153     return FEM_TEMPLATE_ELEMENT<3>::get_nb_fem_noeud();
154     }
155    
156     FEM_NOEUD* FEM_TRIANGLE3::get_fem_noeud(int num)
157     {
158     return FEM_TEMPLATE_ELEMENT<3>::get_fem_noeud(num);
159     }
160    
161     void FEM_TRIANGLE3::change_noeud(int num,FEM_NOEUD* noeud)
162     {
163     FEM_TEMPLATE_ELEMENT<3>::change_noeud(num,noeud);
164     }
165    
166     BOITE_3D& FEM_TRIANGLE3::get_boite_3D(void)
167     {
168     return FEM_TEMPLATE_ELEMENT<3>::get_boite_3D();
169     }
170    
171    
172 francois 763 void FEM_TRIANGLE3::enregistrer(std::ostream& o,double version)
173 francois 283 {
174 francois 378 if (maillage!=NULL)
175     if (get_lien_topologie()!=NULL) o << "%" << get_id() << "=FEM_TRIANGLE3($"<< get_lien_topologie()->get_id() << ",$" << maillage->get_id() << ",$" << tab[0]->get_id() << ",$" << tab[1]->get_id() << ",$" << tab[2]->get_id()<< ");" << std::endl;
176     else o << "%" << get_id() << "=FEM_TRIANGLE3(NULL,$" << maillage->get_id() << ",$" << tab[0]->get_id() << ",$" << tab[1]->get_id() << ",$" << tab[2]->get_id()<< ");" << std::endl;
177     else
178     if (get_lien_topologie()!=NULL) o << "%" << get_id() << "=FEM_TRIANGLE3($"<< get_lien_topologie()->get_id() << ",NULL,$" << tab[0]->get_id() << ",$" << tab[1]->get_id() << ",$" << tab[2]->get_id()<< ");" << std::endl;
179     else o << "%" << get_id() << "=FEM_TRIANGLE3(NULL,NULL,$" << tab[0]->get_id() << ",$" << tab[1]->get_id() << ",$" << tab[2]->get_id()<< ");" << std::endl;
180    
181 francois 637 }
182 francois 310
183 francois 637 int FEM_TRIANGLE3::get_nb_pt_gauss(int degre)
184     {
185 francois 754 return OT_POINTS_GAUSS::get_nb_point_tri(degre);
186 francois 637 }
187    
188     void FEM_TRIANGLE3::get_pt_gauss(int degre,int num,double &w,double *uv)
189     {
190 francois 754 return OT_POINTS_GAUSS::get_pt_gauss_tri(degre,num,w,uv);
191 francois 637 }
192 francois 757
193     int FEM_TRIANGLE3::get_degre_gauss(int num)
194     {
195     return OT_POINTS_GAUSS::get_degre_gauss_tri(num);
196     }
197 francois 637
198 francois 638 int FEM_TRIANGLE3::get_nb_fonction_interpolation(void)
199 francois 310 {
200     return 3;
201     }
202 francois 638 int FEM_TRIANGLE3::get_degremax_fonction_interpolation(void)
203     {
204     return 1;
205     }
206 francois 310 double FEM_TRIANGLE3::get_fonction_interpolation(int num,double *uv)
207     {
208     double val;
209     switch (num)
210     {
211     case 1:
212     val=1-uv[0]-uv[1];
213     break;
214     case 2:
215     val=uv[0];
216     break;
217     case 3:
218     val=uv[1];
219     break;
220    
221     }
222     return val;
223     }
224    
225     double FEM_TRIANGLE3::get_fonction_derive_interpolation(int num,int num_variable,double *uv)
226     {
227     double val;
228    
229     switch (num)
230     {
231     case 1:
232     switch (num_variable)
233     {
234     case 1:
235     val=-1;
236     break;
237     case 2:
238     val=-1;
239     break;
240     } break;
241     case 2:
242     switch (num_variable)
243     {
244     case 1:
245     val=1;
246     break;
247     case 2:
248     val=0.;
249     break;
250     }break;
251     case 3:
252     switch (num_variable)
253     {
254     case 1:
255     val=0.;
256     break;
257     case 2:
258     val=1.;
259     break;
260     }break;
261    
262     }
263     return val;
264     }
265    
266     double FEM_TRIANGLE3::get_jacobien(double* jac,double *uv,int& li,int& col,double unite)
267 gervaislavoie 382 {
268 francois 647 int nb=get_nb_fem_noeud();
269     OT_VECTEUR_3D *point= new OT_VECTEUR_3D[nb];
270 gervaislavoie 382
271 francois 647 for (int i=0;i<nb;i++)
272     {
273     point[i].change_x(get_fem_noeud(i)->get_x());
274     point[i].change_y(get_fem_noeud(i)->get_y());
275     point[i].change_z(get_fem_noeud(i)->get_z());
276     }
277    
278     OT_MATRICE_3D J;
279     for (int i=0;i<2;i++)
280     {
281     J(i,0)=0.;J(i,1)=0.;J(i,2)=0.;
282     for (int k=0;k<nb;k++)
283     {
284     J(i,0)=J(i,0)+get_fonction_derive_interpolation(k+1,i+1,uv)*point[k].get_x();
285     J(i,1)=J(i,1)+get_fonction_derive_interpolation(k+1,i+1,uv)*point[k].get_y();
286     J(i,2)=J(i,2)+get_fonction_derive_interpolation(k+1,i+1,uv)*point[k].get_z();
287     }
288     }
289     OT_VECTEUR_3D J1(J(0,0),J(0,1),J(0,2));
290     OT_VECTEUR_3D J2(J(1,0),J(1,1),J(1,2));
291     OT_VECTEUR_3D J3=J1&J2;
292     J3.norme();
293     J(2,0)=J3(0);
294     J(2,1)=J3(1);
295     J(2,2)=J3(2);
296 gervaislavoie 389
297 francois 647 jac[0]=J(0,0);
298     jac[1]=J(0,1);
299     jac[2]=J(0,2);
300 gervaislavoie 389
301 francois 647 jac[3]=J(1,0);
302     jac[4]=J(1,1);
303     jac[5]=J(1,2);
304 gervaislavoie 389
305 francois 647 jac[6]=J(2,0);
306     jac[7]=J(2,1);
307     jac[8]=J(2,2);
308     li=3;col=3;
309 gervaislavoie 389
310 francois 647 return J.get_determinant();
311 francois 310 }
312    
313 francois 635
314     bool FEM_TRIANGLE3::valide_parametre_EF(double* uvw)
315     {
316     if (uvw[0]>=-1e-10)
317     if (uvw[1]>=-1e-10)
318     if (uvw[0]+uvw[1]<=1.+1e-10)
319     return true;
320     return false;
321     }
322    
323 francois 310 void FEM_TRIANGLE3::get_inverse_jacob(double* j,double *uv,double unite)
324     {
325 francois 647 double jac[9];
326     int li,col;
327     get_jacobien(jac,uv,li,col);
328     OT_MATRICE_3D J;
329     J(0,0)=jac[0];
330     J(0,1)=jac[1];
331     J(0,2)=jac[2];
332     J(1,0)=jac[3];
333     J(1,1)=jac[4];
334     J(1,2)=jac[5];
335     J(2,0)=jac[6];
336     J(2,1)=jac[7];
337     J(2,2)=jac[8];
338     OT_MATRICE_3D j_i=J.inverse();
339    
340 gervaislavoie 382
341     j[0]=j_i(0,0);
342 francois 647 j[1]=j_i(0,1);
343     j[2]=j_i(0,2);
344 gervaislavoie 382
345 francois 647 j[3]=j_i(1,0);
346 gervaislavoie 382 j[4]=j_i(1,1);
347 francois 647 j[5]=j_i(1,2);
348 gervaislavoie 382
349 francois 647 j[6]=j_i(2,0);
350     j[7]=j_i(2,1);
351 gervaislavoie 382 j[8]=j_i(2,2);
352 francois 310 }