ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/solveur/src/sl_solvelastoplastique.cpp
Revision: 166
Committed: Thu Oct 23 21:23:28 2008 UTC (16 years, 6 months ago) by souaissa
Original Path: magic/lib/solveur/solveur/src/sl_solvelastoplastique.cpp
File size: 4170 byte(s)
Log Message:
mise a jour du solveur : ajout de la resolution partie plastique

File Contents

# User Rev Content
1 souaissa 166 #include "gestionversion.h"
2     //---------------------------------------------------------------------------
3    
4    
5     #pragma hdrstop
6    
7     #include "sl_solvelastoplastique.h"
8     #include "sl_critere_plastique.h"
9     #include "sl_operateurs_plastique.h"
10     //---------------------------------------------------------------------------
11    
12     #pragma package(smart_init)
13    
14    
15     SL_SOLVE_ELASTOPLASTIQUE::SL_SOLVE_ELASTOPLASTIQUE(int ele_type,char* file_name,unsigned long mai_id,long nbincr,char*out_file)
16     {
17     Nb_Charge=nbincr;
18     EL_TYPE=ele_type;
19     strcpy(nom_fichier,file_name);
20     strcpy(out_fichier,out_file);
21     mggest=new MG_FILE(nom_fichier);
22     mai=mggest->get_fem_maillageid(mai_id);
23     id_mai=mai_id;
24    
25     }
26    
27     SL_SOLVE_ELASTOPLASTIQUE::~SL_SOLVE_ELASTOPLASTIQUE()
28     {
29     delete mggest;
30     delete mai;
31     }
32    
33     void SL_SOLVE_ELASTOPLASTIQUE::resoudre(ostream& os,double*kg,Matrice& contrainte_elast,Matrice& def_elast)
34     {
35     int nb_tetra=mai->get_nb_fem_tetra();
36     Matrice def_plastique (nb_tetra,6,0.);
37     Matrice cntr_plastque(nb_tetra,6,0.);
38     Matrice cntr_elastique(nb_tetra,6,0.);
39     Matrice def_elastique=def_elast*(1./Nb_Charge);
40     Matrice d_sig_elast(nb_tetra,6,0.);
41     d_sig_elast=contrainte_elast*(1./Nb_Charge);
42     double max= d_sig_elast.Get_AbMax();
43    
44     vecteur lamda(nb_tetra);
45     SL_OPERATEURS_PLASTIQUE op_plastique(mggest);
46     double limte_elast=op_plastique.get_limite_elastque();
47     int N_plas= limte_elast/max;
48     Matrice cont_chargement=d_sig_elast*N_plas;
49     Matrice cntr_tot=d_sig_elast*N_plas;
50     def_elastique=def_elastique*N_plas;
51     Matrice S=op_plastique.get_operateur_s(os,kg);
52     Matrice K_1=op_plastique.modifier_matrice_rigidite(os,kg);
53     int param_plastique=0;
54     for(int i=N_plas;i<N_plas+100;i++)
55     {
56     os<<"INCR_CHARGE : "<<i<<endl;
57     Matrice d_sig(nb_tetra,6,0.);
58     Matrice d_sig_plast(nb_tetra,6,0.);
59     Matrice d_eps_plast(nb_tetra,6,0.);
60     vecteur d_lamda(nb_tetra);
61     op_plastique.get_operateur_m(os,S,cont_chargement,d_sig_elast, lamda,d_lamda,param_plastique);
62    
63     if (param_plastique) d_eps_plast=get_def_plastique(cntr_tot,d_lamda) ;
64     if (param_plastique) d_sig_plast=get_cntr_plastique(S,cntr_tot,d_eps_plast);
65    
66     cont_chargement+= d_sig_elast;
67    
68     if (param_plastique) cntr_plastque+= d_sig_plast;
69     def_elastique+= d_eps_plast;
70     d_sig= d_sig_elast+d_sig_plast ;
71     lamda+=d_lamda;
72     cntr_tot+= d_sig;
73     if (param_plastique)
74     {
75     os<<"DEF_PLASTIQUE"<<endl;
76     os<<d_eps_plast<<endl;
77     os<<"SIG_PLASTIQUE"<<endl;
78     os<<d_sig_plast<<endl;
79     os<<"CONT_CHARGEMENT"<<endl;
80     os<<cont_chargement<<endl;
81     os<<"CONTR_TOT"<<endl;
82     os<<cntr_tot<<endl;
83     }
84     }
85    
86     }
87    
88     Matrice SL_SOLVE_ELASTOPLASTIQUE::get_def_plastique(Matrice& cntr_tot,vecteur&lamda)
89     {
90     int nb_tetra=mai->get_nb_fem_tetra();
91     Matrice def_plastique(0,6,0.);
92     SL_CRITERE_PLASTICITE critere;
93    
94     for(int i=0;i<nb_tetra;i++)
95     {
96     vecteur df_sig=critere.get_derivee_du_critere(cntr_tot.Get_Row(i));
97     vecteur df_eps=lamda[i]*df_sig;
98     def_plastique.Add_Row(df_eps,i);
99     }
100    
101     return def_plastique;
102    
103     }
104    
105    
106    
107     Matrice SL_SOLVE_ELASTOPLASTIQUE :: get_cntr_plastique(Matrice S,Matrice& cntr_tot,Matrice& def_plastque)
108     {
109     int nb_tetra=mai->get_nb_fem_tetra();
110     Matrice cntr_plastique(0,6,0.);
111    
112     for (int a=0;a<nb_tetra;a++)
113     {
114     vecteur res(6);
115    
116     for(int b=0;b<nb_tetra;b++)
117     {
118     Matrice sab=S.Get_SubMatrice(6*a,6*b,6,6);
119     vecteur d_epsp_b=def_plastque.Get_Row(b);
120     res+=sab* d_epsp_b;
121     }
122     cntr_plastique.Add_Row(res,a);
123     }
124    
125     return cntr_plastique;
126     }
127    
128    
129    
130     void SL_SOLVE_ELASTOPLASTIQUE::resoud()
131     {
132     SL_SOLVE_ELASTIQUE solve_elst(EL_TYPE,nom_fichier,id_mai,Nb_Charge);
133     int taille=3*mai->get_nb_fem_noeud();
134     double *fg=new double[taille];
135     double *kg=new double[taille];
136     double *du=new double[taille];
137     Matrice cont(0,6,0.0);
138     Matrice def(0,6,0.0);
139     ofstream os(out_fichier,ios_base::app);
140     solve_elst.resoud(du,kg,fg,cont,def);
141     resoudre(os,kg,cont,def) ;
142    
143     delete [] kg;
144     delete [] fg;
145     delete [] du;
146     }