rio/Herezh_dev
1
0
Fork 0
Code Issues Pull requests Projects Releases Wiki Activity
Herezh_dev/herezh_pp/comportement/Hyper_elastique/HyperDN.cc

228 lines
11 KiB
C++
Executable file
Raw Blame History

// This file is part of the Herezh++ application.
//
// The finite element software Herezh++ is dedicated to the field
// of mechanics for large transformations of solid structures.
// It is developed by Gérard Rio (APP: IDDN.FR.010.0106078.000.R.P.2006.035.20600)
// INSTITUT DE RECHERCHE DUPUY DE LÔME (IRDL) <https://www.irdl.fr/>.
//
// Herezh++ is distributed under GPL 3 license ou ultérieure.
//
// Copyright (C) 1997-2021 Université Bretagne Sud (France)
// AUTHOR : Gérard Rio
// E-MAIL : gerardrio56@free.fr
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License,
// or (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty
// of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
// See the GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program. If not, see <https://www.gnu.org/licenses/>.
//
// For more information, please consult: <https://herezh.irdl.fr/>.
//#include "HyperDN.h"
// Constructeur par defaut
template <class TensHH,class TensBB,class TensBH,class TensHB,
class Tens_nHH,class Tens_nBB>
HyperDN<TensHH,TensBB,TensBH,TensHB,Tens_nHH,Tens_nBB>::HyperDN () :
Loi_comp_abstraite()
{}
// Constructeur utile si l'identificateur du nom de la loi
// de comportement et la dimension sont connus
template <class TensHH,class TensBB,class TensBH,class TensHB,
class Tens_nHH,class Tens_nBB>
HyperDN<TensHH,TensBB,TensBH,TensHB,Tens_nHH,Tens_nBB>
::HyperDN (Enum_comp id_compor,Enum_categorie_loi_comp categorie_comp,int dimension) :
Loi_comp_abstraite(id_compor,categorie_comp,dimension)
{}
// Constructeur utile si l'identificateur du nom de la loi
// de comportement et la dimension sont connus
template <class TensHH,class TensBB,class TensBH,class TensHB,
class Tens_nHH,class Tens_nBB>
HyperDN<TensHH,TensBB,TensBH,TensHB,Tens_nHH,Tens_nBB>
::HyperDN (char* nom,Enum_categorie_loi_comp categorie_comp,int dimension) :
Loi_comp_abstraite(nom,categorie_comp,dimension)
{} // DESTRUCTEUR :
template <class TensHH,class TensBB,class TensBH,class TensHB,
class Tens_nHH,class Tens_nBB>
HyperDN<TensHH,TensBB,TensBH,TensHB,Tens_nHH,Tens_nBB>
::~HyperDN ()
{}
// constructeur de copie
template <class TensHH,class TensBB,class TensBH,class TensHB,
class Tens_nHH,class Tens_nBB>
HyperDN<TensHH,TensBB,TensBH,TensHB,Tens_nHH,Tens_nBB>
::HyperDN (const HyperDN & a) :
Loi_comp_abstraite (a)
{}
// ========== codage des METHODES VIRTUELLES protegees:================
// calcul des contraintes a tdt
template <class TensHH,class TensBB,class TensBH,class TensHB,
class Tens_nHH,class Tens_nBB>
void HyperDN<TensHH,TensBB,TensBH,TensHB,Tens_nHH,Tens_nBB>
::Calcul_SigmaHH
(TenseurHH& ,TenseurBB& ,DdlElement & tab_ddl,
TenseurBB & ,TenseurHH & ,BaseB& ,BaseH& ,TenseurBB & epsBB_,TenseurBB & ,
TenseurBB & gijBB_,TenseurHH & gijHH_,Tableau <TenseurBB *>& d_gijBB_,
double& jacobien_0,double& jacobien_,TenseurHH & sigHH
,EnergieMeca & energ,const EnergieMeca & energ_t,double& module_compressibilite,double& module_cisaillement
,const Met_abstraite::Expli_t_tdt& )
{
#ifdef MISE_AU_POINT
if (tab_ddl.NbDdl() != d_gijBB_.Taille())
{ cout << "\nErreur : le nb de ddl est != de la taille de d_gijBB_ !\n";
cout << " HyperDN::Calcul_SigmaHH\n";
Sortie(1);
};
#endif
TensBH epsBH; // la déformation en mixte
TensBH * IdGBH; // l'identitée
double alpha_0,alpha_1,alpha_2; // coefficients pour le calcul de la loi de comportement
// Calcul des 3 invariants , de epsBH,
// retour de IdGBH qui pointe sur le bon tenseur
double Ieps,V,bIIb;
IdGBH = Invariants (epsBB_,gijBB_,gijHH_,jacobien_0,jacobien_,Ieps,V,bIIb,epsBH);
// calcul du potentiel et de ses dérivées premières
double E,EV,EbIIb,EIeps; // potentiel et dérivées
Potentiel(jacobien_0,Ieps,V,bIIb,E,EV,EbIIb,EIeps);
// calcul des alphas
Alpha(E,EV,EbIIb,EIeps,jacobien_0,Ieps,V,bIIb,alpha_0,alpha_1,alpha_2);
// calcul des contraintes finales
TensBH sigBH = alpha_0 * (*IdGBH) + alpha_1 * epsBH
+ alpha_2 * (epsBH * epsBH);
sigHH = gijHH_ * sigBH; // en deux fois contravariants
}
// calcul des contraintes et de ses variations a t+dt
template <class TensHH,class TensBB,class TensBH,class TensHB,
class Tens_nHH,class Tens_nBB>
void HyperDN<TensHH,TensBB,TensBH,TensHB,Tens_nHH,Tens_nBB>
::Calcul_DsigmaHH_tdt
(TenseurHH& ,TenseurBB& ,DdlElement & tab_ddl
,BaseB& ,TenseurBB & ,TenseurHH & ,
BaseB& ,Tableau <BaseB> & ,BaseH& ,Tableau <BaseH> & ,
TenseurBB & epsBB_tdt,Tableau <TenseurBB *>& d_epsBB,
TenseurBB & ,TenseurBB & gijBB_tdt,TenseurHH & gijHH_tdt,
Tableau <TenseurBB *>& d_gijBB_tdt,
Tableau <TenseurHH *>& d_gijHH_tdt,double& jacobien_0,double& jacobien_tdt,
Vecteur& d_jacobien_tdt,TenseurHH& sigHH,Tableau <TenseurHH *>& d_sigHH
,EnergieMeca & energ,const EnergieMeca & energ_t,double& module_compressibilite,double& module_cisaillement
,const Met_abstraite::Impli& )
{
#ifdef MISE_AU_POINT
if (tab_ddl.NbDdl() != d_gijBB_tdt.Taille())
{ cout << "\nErreur : le nb de ddl est != de la taille de d_gijBB_t !\n";
cout << " HyperDN::Calcul_DsigmaHH_tdt\n";
Sortie(1);
};
#endif
int nbddl = tab_ddl.NbDdl();
TensBH epsBH; // la déformation en mixte
Tableau<TensBH> depsBH(nbddl);
TensBH * IdGBH; // l'identitée
// def des alphas
double alpha_0,alpha_1,alpha_2;
Tableau<double> dalpha_0(nbddl),dalpha_1(nbddl),dalpha_2(nbddl);
// calcul des trois invariants et de leurs variations, de epsBH,
// et de sa variation, puis retour de IdGBH qui pointe sur le bon tenseur
double Ieps,V,bIIb;
Tableau<double> dIeps(nbddl),dV(nbddl),dbIIb(nbddl);
IdGBH = Invariants_et_var (epsBB_tdt,gijBB_tdt,d_gijBB_tdt,
gijHH_tdt,d_gijHH_tdt,jacobien_0,jacobien_tdt,d_jacobien_tdt,
Ieps,dIeps,V,dV,bIIb,dbIIb,epsBH,depsBH);
// calcul du potentiel et de ses dérivées premières et secondes
double E,EV,EbIIb,EIeps,EVV,EbIIb2,EIeps2,EVbIIb,EVIeps,EbIIbIeps; // potentiel et variations
Potentiel_et_var(jacobien_0,Ieps,V,bIIb,E,EV,EbIIb,EIeps,EVV,EbIIb2,EIeps2,
EVbIIb,EVIeps,EbIIbIeps );
// calcul des coefficients alpha et de leurs variations par rapport aux degrés de libertés
Tableau<double> dE(nbddl),dEV(nbddl),dEbIIb(nbddl),dEVV(nbddl),dEQQ(nbddl),dEVQ(nbddl);
Alpha_var ( E, EV, EbIIb, EIeps,EVV, EbIIb2, EIeps2, EVbIIb, EVIeps, EbIIbIeps,
jacobien_0,Ieps,dIeps,V,dV,bIIb,dbIIb,alpha_0,dalpha_0,alpha_1,dalpha_1,alpha_2,dalpha_2);
// calcul des contraintes finales
TensBH epsepsBH = epsBH * epsBH; // epsBH:epsBH
TensBH sigBH = alpha_0 * (*IdGBH) + alpha_1 * epsBH
+ alpha_2 * epsepsBH;
sigHH = gijHH_tdt * sigBH; // en deux fois contravariants
// cas le la variation du tenseur des contraintes
for (int i = 1; i<= nbddl; i++)
{ const TensHH & dgijHH = (*(d_gijHH_tdt(i))) ; // pour simplifier l'ecriture
TensHH& dsigHH = *((TensHH*) (d_sigHH(i))); //
TensBB depsBB = 0.5 * (*(d_gijBB_tdt(i)));
TensBH depsBH = epsBB_tdt * dgijHH + depsBB * gijHH_tdt ;
TensBH dsigBH = dalpha_0(i) * (*IdGBH) + dalpha_1(i) * epsBH
+ dalpha_2(i) * epsepsBH + alpha_1 * depsBH + alpha_2 * (depsBH*epsBH + epsBH*depsBH);
dsigHH = dgijHH * sigBH + gijHH_tdt * dsigBH;
}
}
// -------------------------------------------------------------------------------------------------------------
// METHODES PROTEGEES :
// -------------------------------------------------------------------------------------------------------------
// calcul des coefficients alpha
template <class TensHH,class TensBB,class TensBH,class TensHB,
class Tens_nHH,class Tens_nBB>
inline void HyperDN<TensHH,TensBB,TensBH,TensHB,Tens_nHH,Tens_nBB>
::Alpha
(double& ,double& EV,double& EbIIb,double& EIeps,
double& jacobien_0,double & Ieps,double & V,double& ,
double & alpha_0,double & alpha_1,double & alpha_2)
{
double unsurrgOV = 1./(jacobien_0*V);
alpha_0 = unsurrgOV * ( EIeps + V * EV - Ieps/3. * EbIIb);
alpha_1 = unsurrgOV * ( -2. * EIeps + (2 * Ieps/3. + 1.) * EbIIb);
alpha_2 = unsurrgOV * ( - 2. * EbIIb);
}
// calcul des coefficients alpha et de leurs variations par rapport aux degrés de libertés
template <class TensHH,class TensBB,class TensBH,class TensHB,
class Tens_nHH,class Tens_nBB>
inline void HyperDN<TensHH,TensBB,TensBH,TensHB,Tens_nHH,Tens_nBB>
:: Alpha_var (
double& ,double& EV,double& EbIIb,double& EIeps,
double& EVV,double& EbIIb2,double& EIeps2,
double& EVbIIb,double& EVIeps,double& EbIIbIeps,
double& jacobien_0,double& Ieps,Tableau<double> & dIeps,double& V,Tableau<double> & dV,
double& ,Tableau<double> & dbIIb,
double& alpha_0,Tableau<double> & dalpha_0,
double& alpha_1,Tableau<double> & dalpha_1,
double& alpha_2,Tableau<double> & dalpha_2)
{ // calcul des coefficients alpha
double unsurrgOV = 1./(jacobien_0*V);
double al1 = ( EIeps + V * EV - Ieps/3. * EbIIb); // pour optimiser les calculs de variations
alpha_0 = unsurrgOV * al1;
double al2 = ( -2. * EIeps + (2 * Ieps/3. + 1.) * EbIIb); // pour optimiser les calculs de variations
alpha_1 = unsurrgOV * al2;
double al3 = ( - 2. * EbIIb); // pour optimiser les calculs de variations
alpha_2 = unsurrgOV * al3;
// cas le la variation des coefficients alpha
int nbddl = dIeps.Taille();
for (int i=1; i<= nbddl; i++)
{ // variation de unsurrgOV
double d_unsurrgOV = - unsurrgOV* unsurrgOV * (jacobien_0 * dV(i));
// variations des alphas
dalpha_0(i) = d_unsurrgOV * al1 + unsurrgOV *
( dIeps(i) * EIeps2 + dbIIb(i) * EbIIbIeps + dV(i) * EVIeps);
dalpha_1(i) = d_unsurrgOV * al2 + unsurrgOV *
( dIeps(i) * EbIIbIeps + dbIIb(i) * EbIIb2 + dV(i) * EVbIIb);
dalpha_2(i) = d_unsurrgOV * al3 + unsurrgOV *
( dIeps(i) * EVIeps + dbIIb(i) * EVbIIb + dV(i) * EVV);
};
}
Reference in a new issue View git blame Copy permalink