1547 lines
75 KiB
C++
1547 lines
75 KiB
C++
// FICHIER : BielletteQ.cp
|
|
// CLASSE : BielletteQ
|
|
|
|
// 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-2022 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 "Debug.h"
|
|
|
|
|
|
# include <iostream>
|
|
using namespace std; //introduces namespace std
|
|
#include <stdlib.h>
|
|
#include "Sortie.h"
|
|
|
|
|
|
#include "BielletteQ.h"
|
|
#include "FrontPointF.h"
|
|
#include "FrontSegQuad.h"
|
|
#include "TypeConsTens.h"
|
|
#include "TypeQuelconqueParticulier.h"
|
|
|
|
//----------------------------------------------------------------
|
|
// def des donnees commune a tous les elements
|
|
// la taille n'est pas defini ici car elle depend de la lecture
|
|
//----------------------------------------------------------------
|
|
BielletteQ::DonneeCommune * BielletteQ::doCo = NULL;
|
|
BielletteQ::UneFois BielletteQ::unefois;
|
|
BielletteQ::NombresConstruire BielletteQ::nombre_V;
|
|
BielletteQ::ConstrucElementbielQ BielletteQ::construcElementbielQ;
|
|
|
|
// constructeur définissant les nombres (de noeud, de point d'integ ..)
|
|
// utilisé dans la construction des éléments
|
|
BielletteQ::NombresConstruire::NombresConstruire()
|
|
{ nbne = 3; // le nombre de noeud de l'élément
|
|
nbneA = 3;// le nombre de noeud des aretes
|
|
nbi = 2; // le nombre de point d'intégration pour le calcul mécanique
|
|
nbiEr = 3;// le nombre de point d'intégration pour le calcul d'erreur
|
|
nbiA = 2; // le nombre de point d'intégration pour le calcul de second membre linéique
|
|
nbiMas = 3; // le nombre de points d'intégration pour le calcul de la matrice masse
|
|
nbiHour = 0; // le nombre de point d'intégration un blocage d'hourglass
|
|
};
|
|
|
|
// ---- definition du constructeur de la classe conteneur de donnees communes ------------
|
|
BielletteQ::DonneeCommune::DonneeCommune (GeomSeg& seg,DdlElement& tab,DdlElement& tabErr,DdlElement& tab_Err1Sig,
|
|
Met_biellette & met_bie,
|
|
Tableau <Vecteur *> & resEr,Mat_pleine& raidEr,
|
|
GeomSeg& seEr,Vecteur& residu_int,Mat_pleine& raideur_int,
|
|
Tableau <Vecteur* > & residus_extN,Tableau <Mat_pleine* >& raideurs_extN,
|
|
Tableau <Vecteur* > & residus_extA,Tableau <Mat_pleine* >& raideurs_extA,
|
|
Mat_pleine& mat_masse,GeomSeg& seMa,int nbi,GeomSeg* segHourg) :
|
|
segment(seg),tab_ddl(tab),tab_ddlErr(tabErr),tab_Err1Sig11(tab_Err1Sig)
|
|
,met_BielletteQ(met_bie)
|
|
,matGeom(tab.NbDdl(),tab.NbDdl())
|
|
,matInit(tab.NbDdl(),tab.NbDdl())
|
|
,d_epsBB(tab.NbDdl()),d_sigHH(tab.NbDdl()),d2_epsBB(nbi)
|
|
,resErr(resEr),raidErr(raidEr)
|
|
,segmentEr(seEr)
|
|
,residu_interne(residu_int),raideur_interne(raideur_int)
|
|
,residus_externeN(residus_extN),raideurs_externeN(raideurs_extN)
|
|
,residus_externeA(residus_extA),raideurs_externeA(raideurs_extA)
|
|
,matrice_masse(mat_masse)
|
|
,segmentMas(seMa),segmentHourg(segHourg)
|
|
{
|
|
int nbddl = tab.NbDdl();
|
|
for (int ni=1;ni<=nbi;ni++)
|
|
{d2_epsBB(ni).Change_taille(nbddl);
|
|
for (int i1=1; i1<= nbddl; i1++)
|
|
for (int i2=1; i2<= nbddl; i2++)
|
|
d2_epsBB(ni)(i1,i2) = NevezTenseurBB (1);
|
|
};
|
|
int tailledeps = d_epsBB.Taille();
|
|
for (int i=1;i<= tailledeps; i++)
|
|
d_epsBB(i) = NevezTenseurBB (1);
|
|
int tailledsig = d_sigHH.Taille();
|
|
for (int j=1;j<= tailledsig; j++)
|
|
d_sigHH(j) = NevezTenseurHH (1);
|
|
};
|
|
|
|
BielletteQ::DonneeCommune::DonneeCommune(DonneeCommune& a) :
|
|
segment(a.segment),tab_ddl(a.tab_ddl),tab_ddlErr(a.tab_ddlErr),tab_Err1Sig11(a.tab_Err1Sig11)
|
|
,met_BielletteQ(a.met_BielletteQ),matGeom(a.matGeom),matInit(a.matInit)
|
|
,d2_epsBB(a.d2_epsBB),resErr(a.resErr),raidErr(a.raidErr),segmentEr(a.segmentEr)
|
|
,d_epsBB(a.d_epsBB),d_sigHH(a.d_sigHH)
|
|
,residu_interne(a.residu_interne),raideur_interne(a.raideur_interne)
|
|
,residus_externeN(a.residus_externeN),raideurs_externeN(a.raideurs_externeN)
|
|
,residus_externeA(a.residus_externeA),raideurs_externeA(a.raideurs_externeA)
|
|
,matrice_masse(a.matrice_masse),segmentMas(a.segmentMas),segmentHourg(a.segmentHourg)
|
|
{ int nbddl = d_sigHH.Taille();
|
|
int nbi=d2_epsBB.Taille();
|
|
for (int ni=1;ni<=nbi;ni++)
|
|
for (int i1=1; i1<= nbddl; i1++)
|
|
for (int i2=1; i2<= nbddl; i2++)
|
|
d2_epsBB(ni)(i1,i2) = NevezTenseurBB (*(a.d2_epsBB(ni)(i1,i2)));
|
|
int tailledeps = d_epsBB.Taille();
|
|
for (int i=1;i<= tailledeps; i++)
|
|
d_epsBB(i) = NevezTenseurBB (*(a.d_epsBB(i)));
|
|
int tailledsig = d_sigHH.Taille();
|
|
for (int j=1;j<= tailledsig; j++)
|
|
d_sigHH(j) = NevezTenseurHH (*(d_sigHH(j)));
|
|
};
|
|
|
|
BielletteQ::DonneeCommune::~DonneeCommune()
|
|
{ int nbddl = tab_ddl.NbDdl();
|
|
int nbi=d2_epsBB.Taille();
|
|
for (int ni=1;ni<=nbi;ni++)
|
|
for (int i1=1; i1<= nbddl; i1++)
|
|
for (int i2=1; i2<= nbddl; i2++)
|
|
delete d2_epsBB(ni)(i1,i2);
|
|
int tailledeps = d_epsBB.Taille();
|
|
for (int i=1;i<= tailledeps; i++)
|
|
delete d_epsBB(i);
|
|
int tailledsig = d_sigHH.Taille();
|
|
for (int j=1;j<= tailledsig; j++)
|
|
delete d_sigHH(j);
|
|
};
|
|
|
|
// ---------- fin definition de la classe conteneur de donnees communes ------------
|
|
// -+-+ definition de la classe contenant tous les indicateurs qui sont modifiés une seule fois -+-+-+
|
|
BielletteQ::UneFois::UneFois () : // constructeur par défaut
|
|
doCoMemb(NULL),CalResPrem_t(0),CalResPrem_tdt(0),CalimpPrem(0),dualSortbiel(0)
|
|
,CalSMlin_t(0),CalSMlin_tdt(0),CalSMRlin(0)
|
|
,CalSMvol_t(0),CalSMvol_tdt(0),CalSMvol(0)
|
|
,CalDynamique(0),CalPt_0_t_tdt(0)
|
|
,nbelem_in_Prog(0)
|
|
{};
|
|
BielletteQ::UneFois::~UneFois ()
|
|
{ delete doCoMemb;
|
|
};
|
|
|
|
// -+-+ fin definition de la classe contenant tous les indicateurs qui sont modifiés une seule fois -+-+-+
|
|
|
|
|
|
// ---------- fin definition de la classe conteneur de donnees communes ------------
|
|
|
|
BielletteQ::BielletteQ () :
|
|
// Constructeur par defaut
|
|
ElemMeca(),lesPtMecaInt(),donnee_specif()
|
|
|
|
{ lesPtIntegMecaInterne = &lesPtMecaInt; // association avec le pointeur d'ElemMeca
|
|
id_interpol=BIE2; // donnees de la classe mere
|
|
id_geom=POUT; //
|
|
tab_noeud.Change_taille(nombre_V.nbne);
|
|
// initialisation par défaut
|
|
doCo = BielletteQ::Init ();
|
|
unefois.nbelem_in_Prog++;
|
|
};
|
|
|
|
BielletteQ::BielletteQ (double sect,int num_maill,int num_id):
|
|
// Constructeur utile si la section de l'element et
|
|
// le numero de l'element sont connus
|
|
ElemMeca(num_maill,num_id,BIE2,POUT),lesPtMecaInt(),donnee_specif(sect)
|
|
{ lesPtIntegMecaInterne = &lesPtMecaInt; // association avec le pointeur d'ElemMeca
|
|
tab_noeud.Change_taille(nombre_V.nbne);
|
|
// initialisation
|
|
doCo = BielletteQ::Init (donnee_specif);
|
|
unefois.nbelem_in_Prog++;
|
|
};
|
|
|
|
// Constructeur fonction d'un numero de maillage et d'identification
|
|
BielletteQ::BielletteQ (int num_maill,int num_id) :
|
|
ElemMeca(num_maill,num_id,BIE2,POUT),lesPtMecaInt(),donnee_specif()
|
|
{ lesPtIntegMecaInterne = &lesPtMecaInt; // association avec le pointeur d'ElemMeca
|
|
tab_noeud.Change_taille(nombre_V.nbne);
|
|
// initialisation par défaut
|
|
doCo = BielletteQ::Init ();
|
|
unefois.nbelem_in_Prog++;
|
|
};
|
|
|
|
BielletteQ::BielletteQ (double sect,int num_maill,int num_id,const Tableau<Noeud *>& tab):
|
|
// Constructeur utile si la section de l'element, le numero de l'element et
|
|
// le tableau des noeuds sont connus
|
|
ElemMeca(num_maill,num_id,tab,BIE2,POUT),lesPtMecaInt(),donnee_specif(sect)
|
|
{ lesPtIntegMecaInterne = &lesPtMecaInt; // association avec le pointeur d'ElemMeca
|
|
if (tab_noeud.Taille() != 3)
|
|
{ cout << "\n erreur de dimensionnement du tableau de noeud, on doit avoir 3 noeuds \n";
|
|
cout << " BielletteQ::BielletteQ (double sect,int num_id,const Tableau<Noeud *>& tab)\n";
|
|
Sortie (1);
|
|
};
|
|
// construction du tableau de ddl spécifique à l'élément pour ses
|
|
ConstTabDdl();
|
|
// initialisation
|
|
bool sans_init_noeud = true;
|
|
doCo = BielletteQ::Init (donnee_specif,sans_init_noeud);
|
|
unefois.nbelem_in_Prog++;
|
|
};
|
|
|
|
BielletteQ::BielletteQ (const BielletteQ& biel) :
|
|
ElemMeca (biel),lesPtMecaInt(biel.lesPtMecaInt),donnee_specif(biel.donnee_specif)
|
|
// Constructeur de copie
|
|
{ lesPtIntegMecaInterne = &lesPtMecaInt; // association avec le pointeur d'ElemMeca
|
|
// initialisation
|
|
unefois.nbelem_in_Prog++;
|
|
};
|
|
|
|
// Destructeur
|
|
BielletteQ::~BielletteQ ()
|
|
{ LibereTenseur();
|
|
unefois.nbelem_in_Prog--;
|
|
// defArete pointe sur la même grandeur que def donc pour éviter
|
|
// une destruction lors du destructeur générale dans elemméca on le met à null
|
|
defArete(1) = NULL;
|
|
Destruction();
|
|
};
|
|
|
|
// Lecture des donnees de la classe sur fichier
|
|
void
|
|
BielletteQ::LectureDonneesParticulieres
|
|
(UtilLecture * entreePrinc,Tableau<Noeud *> * tabMaillageNoeud)
|
|
{ int nb;
|
|
tab_noeud.Change_taille(nombre_V.nbne);
|
|
for (int i=1; i<= nombre_V.nbne; i++)
|
|
{ *(entreePrinc->entree) >> nb;
|
|
if ((entreePrinc->entree)->rdstate() == 0)
|
|
// pour mémoire ici on a
|
|
/* enum io_state
|
|
{ badbit = 1<<0, // -> 1 dans rdstate()
|
|
eofbit = 1<<1, // -> 2
|
|
failbit = 1<<2, // -> 4
|
|
goodbit = 0 // -> O
|
|
};*/
|
|
tab_noeud(i) = (*tabMaillageNoeud)(nb); // lecture normale
|
|
#ifdef ENLINUX
|
|
else if ((entreePrinc->entree)->fail())
|
|
// on a atteind la fin de la ligne et on appelle un nouvel enregistrement
|
|
{ entreePrinc->NouvelleDonnee(); // lecture d'un nouvelle enregistrement
|
|
*(entreePrinc->entree) >> nb;
|
|
tab_noeud(i) = (*tabMaillageNoeud)(nb); // lecture normale
|
|
}
|
|
#else
|
|
/* #ifdef SYSTEM_MAC_OS_X_unix
|
|
else if ((entreePrinc->entree)->fail())
|
|
// on a atteind la fin de la ligne et on appelle un nouvel enregistrement
|
|
{ entreePrinc->NouvelleDonnee(); // lecture d'un nouvelle enregistrement
|
|
*(entreePrinc->entree) >> nb;
|
|
tab_noeud(i) = (*tabMaillageNoeud)(nb); // lecture normale
|
|
}
|
|
#else*/
|
|
else if ((entreePrinc->entree)->eof())
|
|
// la lecture est bonne mais on a atteind la fin de la ligne
|
|
{ tab_noeud(i) = (*tabMaillageNoeud)(nb); // lecture
|
|
// si ce n'est pas la fin de la lecture on appelle un nouvel enregistrement
|
|
if (i != nombre_V.nbne) entreePrinc->NouvelleDonnee(); // lecture d'un nouvelle enregistrement
|
|
}
|
|
// #endif
|
|
#endif
|
|
else // cas d'une erreur de lecture
|
|
{ cout << "\n erreur de lecture inconnue ";
|
|
entreePrinc->MessageBuffer("** lecture des données particulières **");
|
|
cout << "BielletteQ::LectureDonneesParticulieres";
|
|
Affiche();
|
|
Sortie (1);
|
|
}
|
|
}
|
|
// construction du tableau de ddl des noeuds de BielletteQ
|
|
ConstTabDdl();
|
|
};
|
|
|
|
// calcul d'un point dans l'élément réel en fonction des coordonnées dans l'élément de référence associé
|
|
// temps: indique si l'on veut les coordonnées à t = 0, ou t ou tdt
|
|
// 1) cas où l'on utilise la place passée en argument
|
|
Coordonnee & BielletteQ::Point_physique(const Coordonnee& c_int,Coordonnee & co,Enum_dure temps)
|
|
{ BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// a) on commence par définir les bonnes grandeurs dans la métrique
|
|
if( unefois.CalPt_0_t_tdt == 0)
|
|
{ unefois.CalPt_0_t_tdt = 1;
|
|
Tableau<Enum_variable_metrique> tab(3);
|
|
tab(1)=iM0;tab(2)=iMt;tab(3)=iMtdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
// b) calcul de l'interpolation
|
|
const Vecteur& phi = doCo->segment.Phi_point(c_int);
|
|
// c) calcul du point
|
|
switch (temps)
|
|
{ case TEMPS_0 : co = doCo->met_BielletteQ.PointM_0(tab_noeud,phi); break;
|
|
case TEMPS_t : co = doCo->met_BielletteQ.PointM_t(tab_noeud,phi); break;
|
|
case TEMPS_tdt : co = doCo->met_BielletteQ.PointM_tdt(tab_noeud,phi); break;
|
|
}
|
|
// d) retour
|
|
return co;
|
|
};
|
|
|
|
// 3) cas où l'on veut les coordonnées aux 1, 2 ou trois temps selon la taille du tableau t_co
|
|
void BielletteQ::Point_physique(const Coordonnee& c_int,Tableau <Coordonnee> & t_co)
|
|
{ BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// a) on commence par définir les bonnes grandeurs dans la métrique
|
|
if( unefois.CalPt_0_t_tdt == 0)
|
|
{ unefois.CalPt_0_t_tdt = 1;
|
|
Tableau<Enum_variable_metrique> tab(3);
|
|
tab(1)=iM0;tab(2)=iMt;tab(3)=iMtdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
// b) calcul de l'interpolation
|
|
const Vecteur& phi = doCo->segment.Phi_point(c_int);
|
|
// c) calcul des point
|
|
switch (t_co.Taille())
|
|
{ case 3 : t_co(3) = doCo->met_BielletteQ.PointM_tdt(tab_noeud,phi);
|
|
case 2 : t_co(2) = doCo->met_BielletteQ.PointM_t(tab_noeud,phi);
|
|
case 1 : t_co(1) = doCo->met_BielletteQ.PointM_0(tab_noeud,phi);
|
|
}
|
|
};
|
|
|
|
// Calcul du residu local à t ou tdt en fonction du booleen atdt
|
|
Vecteur* BielletteQ::CalculResidu (bool atdt,const ParaAlgoControle & pa)
|
|
{ BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
Tableau <TenseurBB *>& d_epsBB = (doCo->d_epsBB);// "
|
|
// dimensionnement de la metrique
|
|
if (!atdt)
|
|
{if( unefois.CalResPrem_t == 0)
|
|
{ unefois.CalResPrem_t = 1;
|
|
Tableau<Enum_variable_metrique> tab(10);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igijBB_0;tab(4) = igijBB_t;
|
|
tab(5) = igijHH_t; tab(6) = id_giB_t; tab(7) = id_gijBB_t ;
|
|
tab(8) = id_gijBB_t ;tab(9) = igradVBB_t; tab(10) = iVt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};}
|
|
else
|
|
{if( unefois.CalResPrem_tdt == 0)
|
|
{unefois.CalResPrem_tdt = 1;
|
|
Tableau<Enum_variable_metrique> tab(11);
|
|
tab(1) = igiB_0; tab(2) = igiB_tdt; tab(3) = igijBB_0;tab(4) = igijBB_tdt;
|
|
tab(5) = igijHH_tdt; tab(6) = id_giB_tdt; tab(7) = id_gijBB_tdt ;
|
|
tab(8) = id_gijBB_t ;tab(9) = id_gijBB_tdt ;
|
|
tab(10) = igradVBB_tdt;tab(11) = iVtdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};};
|
|
// initialisation du résidu
|
|
residu->Zero();
|
|
Vecteur poids =(doCo->segment).TaWi();
|
|
// poids(1) *= donnee_specif.section;
|
|
ElemMeca::Cal_explicit ( doCo->tab_ddl,d_epsBB,nombre_V.nbi,poids,pa,atdt);
|
|
// prise en compte de la section sauf dans le cas d'une loi rien
|
|
double section_moyenne = 0.; // où on met tout à 0
|
|
// et surtout on ne recalcule pas la section car on n'a pas de compressibilité avec une loi rien !!
|
|
if (!Loi_rien(loiComp->Id_comport()))
|
|
section_moyenne = CalSectionMoyenne_et_vol_pti(atdt);
|
|
(*residu) *= section_moyenne;
|
|
energie_totale *= section_moyenne;
|
|
E_Hourglass *= section_moyenne; // meme si l'énergie d'hourglass est nulle
|
|
E_elem_bulk_tdt*= section_moyenne; // idem
|
|
P_elem_bulk *= section_moyenne; // idem
|
|
volume *= section_moyenne;
|
|
return residu;
|
|
};
|
|
|
|
// Calcul du residu local et de la raideur locale,
|
|
// pour le schema implicite
|
|
Element::ResRaid BielletteQ::Calcul_implicit (const ParaAlgoControle & pa)
|
|
{ BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
Tableau <TenseurBB *>& d_epsBB = (doCo->d_epsBB);// "
|
|
Tableau <TenseurHH *>& d_sigHH = (doCo->d_sigHH);// "
|
|
bool cald_Dvirtuelle = false;
|
|
if (unefois.CalimpPrem == 0)
|
|
{ unefois.CalimpPrem = 1;
|
|
Tableau<Enum_variable_metrique> tab(20);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igiB_tdt; tab(4) = igijBB_0;
|
|
tab(5) = igijBB_t;tab(6) = igijBB_tdt; tab(7) = igijHH_tdt; tab(8) = id_giB_tdt;
|
|
tab(9) = id_gijBB_tdt ;tab(10) = igiH_tdt;tab(11) = id_giH_tdt;
|
|
tab(12) = id_gijHH_tdt;tab(13) = id_jacobien_tdt;tab(14) = id2_gijBB_tdt;
|
|
tab(15) = id_gijBB_t ;tab(16) = id_gijBB_tdt ;tab(17) = idMtdt ;
|
|
tab(18) = igradVBB_tdt; tab(19) = iVtdt; tab(20) = idVtdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
// on ne calcul la dérivée de la déformation virtuelle qu'une fois
|
|
// car elle est constante dans le temps et indépendante des coordonnées
|
|
cald_Dvirtuelle=true;
|
|
};
|
|
|
|
// initialisation du résidu
|
|
residu->Zero();
|
|
// initialisation de la raideur
|
|
raideur->Zero();
|
|
|
|
Vecteur poids =(doCo->segment).TaWi();
|
|
// poids(1) *= donnee_specif.section;
|
|
ElemMeca::Cal_implicit(doCo->tab_ddl, d_epsBB,(doCo->d2_epsBB),d_sigHH,nombre_V.nbi
|
|
,poids,pa,cald_Dvirtuelle);
|
|
// prise en compte de la section sauf dans le cas d'une loi rien
|
|
double section_moyenne = 0.; // où on met tout à 0
|
|
// et surtout on ne recalcule pas la section car on n'a pas de compressibilité avec une loi rien !!
|
|
if (!Loi_rien(loiComp->Id_comport()))
|
|
{const bool atdt=true;
|
|
section_moyenne = CalSectionMoyenne_et_vol_pti(atdt);
|
|
};
|
|
|
|
(*residu) *= section_moyenne;
|
|
(*raideur) *= section_moyenne;
|
|
energie_totale *= section_moyenne;
|
|
E_Hourglass *= section_moyenne; // meme si l'énergie d'hourglass est nulle
|
|
E_elem_bulk_tdt*= section_moyenne; // idem
|
|
P_elem_bulk *= section_moyenne; // idem
|
|
volume *= section_moyenne;
|
|
Element::ResRaid el;
|
|
el.res = residu;
|
|
el.raid = raideur;
|
|
return el;
|
|
};
|
|
|
|
// Calcul de la matrice masse pour l'élément
|
|
Mat_pleine * BielletteQ::CalculMatriceMasse (Enum_calcul_masse type_calcul_masse)
|
|
{ BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// dimensionement de la métrique si nécessaire
|
|
if (unefois.CalDynamique == 0)
|
|
{ unefois.CalDynamique = 1;
|
|
Tableau<Enum_variable_metrique> tab(5);
|
|
tab(1) = igiB_0; tab(2) = igiB_tdt; tab(3) = igijBB_0;
|
|
tab(4) = igijBB_tdt;tab(5) = igradVmoyBB_t;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
// on vérifie le bon dimensionnement de la matrice
|
|
if (type_calcul_masse == MASSE_CONSISTANTE)
|
|
// dans le cas où la masse est consistante il faut la redimensionner
|
|
{ int nbddl = doCo->tab_ddl.NbDdl();
|
|
(doCo->matrice_masse).Initialise (nbddl,nbddl,0.);
|
|
}
|
|
};
|
|
Vecteur poids =(doCo->segmentMas).TaWi();
|
|
// poids *= donnee_specif.section; // prise en compte de la section
|
|
// appel de la routine générale
|
|
ElemMeca::Cal_Mat_masse (doCo->tab_ddl,type_calcul_masse,
|
|
nombre_V.nbiMas,(doCo->segmentMas).TaPhi(),nombre_V.nbne
|
|
,poids);
|
|
(*mat_masse) *= donnee_specif.secti.section0; // prise en compte de l'épaisseur
|
|
return mat_masse;
|
|
};
|
|
|
|
//------- calcul d'erreur, remontée des contraintes -------------------
|
|
// 1) calcul du résidu et de la matrice de raideur pour le calcul d'erreur
|
|
Element::Er_ResRaid BielletteQ::ContrainteAuNoeud_ResRaid()
|
|
{
|
|
if(( unefois.CalResPrem_t == 0)|| (unefois.CalimpPrem == 0))
|
|
{ unefois.CalResPrem_t = 1;
|
|
Tableau<Enum_variable_metrique> tab(10);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igijBB_0;tab(4) = igijBB_t;
|
|
tab(5) = igijHH_t; tab(6) = id_giB_t; tab(7) = id_gijBB_t ;
|
|
tab(8) = id_gijBB_t ;tab(9) = igradVBB_t; tab(10) = iVt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
// appel du programme général
|
|
int tabn_taille = tab_noeud.Taille();
|
|
ElemMeca:: SigmaAuNoeud_ResRaid(tab_noeud.Taille()
|
|
,(doCo->segment).TaPhi()
|
|
,(doCo->segment).TaWi()
|
|
,doCo-> resErr,doCo->raidErr
|
|
,(doCo->segmentEr).TaPhi()
|
|
,(doCo->segmentEr).TaWi());
|
|
// on tient compte de la section à t, supposée déjà calculée ou récupérée
|
|
for (int i=1;i<= tabn_taille;i++)
|
|
(*doCo-> resErr(i)) *= donnee_specif.secti.section_t;
|
|
doCo->raidErr *= donnee_specif.secti.section_t;
|
|
return (Element::Er_ResRaid( &(doCo-> resErr),&(doCo->raidErr)));
|
|
};
|
|
|
|
// 2) remontée aux erreurs aux noeuds
|
|
Element::Er_ResRaid BielletteQ::ErreurAuNoeud_ResRaid()
|
|
{ if(( unefois.CalResPrem_t == 0)|| (unefois.CalimpPrem == 0))
|
|
{ unefois.CalResPrem_t = 1;
|
|
Tableau<Enum_variable_metrique> tab(10);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igijBB_0;tab(4) = igijBB_t;
|
|
tab(5) = igijHH_t; tab(6) = id_giB_t; tab(7) = id_gijBB_t ;
|
|
tab(8) = id_gijBB_t ;tab(9) = igradVBB_t; tab(10) = iVt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
// appel du programme général
|
|
int tabn_taille = tab_noeud.Taille();
|
|
ElemMeca::Cal_ErrAuxNoeuds(tab_noeud.Taille(), (doCo->segment).TaPhi(),
|
|
(doCo->segment).TaWi(),doCo-> resErr );
|
|
// on tient compte de la section à t, supposée déjà calculée ou récupérée
|
|
for (int i=1;i<= tabn_taille;i++)
|
|
(*doCo-> resErr(i)) *= donnee_specif.secti.section_t;
|
|
doCo->raidErr *= donnee_specif.secti.section_t;
|
|
return (Element::Er_ResRaid( &(doCo-> resErr),&(doCo->raidErr)));
|
|
};
|
|
|
|
// actualisation des ddl et des grandeurs actives de t+dt vers t
|
|
void BielletteQ::TdtversT()
|
|
{ lesPtMecaInt.TdtversT(); // contrainte
|
|
for (int ni=1;ni<= nombre_V.nbi; ni++)
|
|
{ if (tabSaveDon(ni) != NULL) tabSaveDon(ni)->TdtversT();
|
|
if (tabSaveTP(ni) != NULL) tabSaveTP(ni)->TdtversT();
|
|
if (tabSaveDefDon(ni) != NULL) tabSaveDefDon(ni)->TdtversT();
|
|
}
|
|
donnee_specif.secti.section_t = donnee_specif.secti.section_tdt;
|
|
ElemMeca::TdtversT_(); // appel de la procédure mère
|
|
};
|
|
// actualisation des ddl et des grandeurs actives de t vers tdt
|
|
void BielletteQ::TversTdt()
|
|
{ lesPtMecaInt.TversTdt(); // contrainte
|
|
for (int ni=1;ni<= nombre_V.nbi; ni++)
|
|
{ if (tabSaveDon(ni) != NULL) tabSaveDon(ni)->TversTdt();
|
|
if (tabSaveTP(ni) != NULL) tabSaveTP(ni)->TversTdt();
|
|
if (tabSaveDefDon(ni) != NULL) tabSaveDefDon(ni)->TversTdt();
|
|
};
|
|
donnee_specif.secti.section_tdt = donnee_specif.secti.section_t;
|
|
ElemMeca::TversTdt_(); // appel de la procédure mère
|
|
};
|
|
|
|
// calcul de l'erreur sur l'élément. Ce calcul n'est disponible
|
|
// qu'une fois la remontée aux contraintes effectuées sinon aucune
|
|
// action. En retour la valeur de l'erreur sur l'élément
|
|
// type indique le type de calcul d'erreur :
|
|
void BielletteQ::ErreurElement(int type,double& errElemRelative
|
|
,double& numerateur, double& denominateur)
|
|
{ if(( unefois.CalResPrem_t == 0)|| (unefois.CalimpPrem == 0))
|
|
{ unefois.CalResPrem_t = 1;
|
|
Tableau<Enum_variable_metrique> tab(10);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igijBB_0;tab(4) = igijBB_t;
|
|
tab(5) = igijHH_t; tab(6) = id_giB_t; tab(7) = id_gijBB_t ;
|
|
tab(8) = id_gijBB_t ;tab(9) = igradVBB_t; tab(10) = iVt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
// appel du programme général
|
|
ElemMeca::Cal_ErrElem(type,errElemRelative,numerateur,denominateur,
|
|
tab_noeud.Taille(),(doCo->segment).TaPhi(),
|
|
(doCo->segment).TaWi(),
|
|
(doCo->segmentEr).TaPhi(),(doCo->segmentEr).TaWi());
|
|
};
|
|
|
|
// mise à jour de la boite d'encombrement de l'élément, suivant les axes I_a globales
|
|
// en retour coordonnées du point mini dans retour.Premier() et du point maxi dans .Second()
|
|
// la méthode est différente de la méthode générale car il faut prendre en compte la section de l'élément
|
|
const DeuxCoordonnees& BielletteQ::Boite_encombre_element(Enum_dure temps)
|
|
{ // on commence par calculer la boite d'encombrement pour l'élément médian
|
|
Element::Boite_encombre_element( temps);
|
|
// ensuite on augmente sytématiquement dans toutes directions d'une valeur sqrt(s)/2 majorée
|
|
double sSur2maj = sqrt(donnee_specif.secti.section_tdt) * 0.5
|
|
* ParaGlob::param->ParaAlgoControleActifs().Extra_boite_prelocalisation();
|
|
// ajout d'un extra dans toutes les directions
|
|
sSur2maj += ParaGlob::param->ParaAlgoControleActifs().Ajout_extra_boite_prelocalisation();
|
|
// mise à jour
|
|
boite_encombre.Premier().Ajout_meme_valeur(-sSur2maj); // le min
|
|
boite_encombre.Second().Ajout_meme_valeur(sSur2maj); // le max
|
|
// retour
|
|
return boite_encombre;
|
|
};
|
|
|
|
//============= lecture écriture dans base info ==========
|
|
|
|
// cas donne le niveau de la récupération
|
|
// = 1 : on récupère tout
|
|
// = 2 : on récupère uniquement les données variables (supposées comme telles)
|
|
void BielletteQ::Lecture_base_info
|
|
(ifstream& ent,const Tableau<Noeud *> * tabMaillageNoeud,const int cas)
|
|
{// tout d'abord appel de la lecture de la classe elem_meca
|
|
ElemMeca::Lecture_bas_inf(ent,tabMaillageNoeud,cas);
|
|
// traitement du cas particulier de la BielletteQ
|
|
switch (cas)
|
|
{ case 1 : // ------- on récupère tout -------------------------
|
|
{ // construction du tableau de ddl des noeuds du triangle
|
|
ConstTabDdl();
|
|
// récup contraintes et déformation
|
|
lesPtMecaInt.Lecture_base_info(ent,cas);
|
|
// les données spécifiques
|
|
string nom;
|
|
ent >> nom;
|
|
if (nom == "epaisStockeDansElement")
|
|
ent >> nom >> donnee_specif.secti.section0
|
|
>> nom >> donnee_specif.secti.section_t
|
|
>> nom >> donnee_specif.secti.section_tdt
|
|
>> nom >> donnee_specif.variation_section;
|
|
break;
|
|
}
|
|
case 2 : // ----------- lecture uniquement de se qui varie --------------------
|
|
{ // récup contraintes et déformation
|
|
lesPtMecaInt.Lecture_base_info(ent,cas);
|
|
// les données spécifiques
|
|
string nom;
|
|
ent >> nom >> donnee_specif.secti.section_tdt;
|
|
donnee_specif.secti.section_t = donnee_specif.secti.section_tdt;
|
|
break;
|
|
}
|
|
default :
|
|
{ cout << "\nErreur : valeur incorrecte du type de lecture !\n";
|
|
cout << "BielletteQ::Lecture_base_info(ofstream& sort,int cas)"
|
|
<< " cas= " << cas << endl;
|
|
Sortie(1);
|
|
}
|
|
}
|
|
};
|
|
|
|
// cas donne le niveau de sauvegarde
|
|
// = 1 : on sauvegarde tout
|
|
// = 2 : on sauvegarde uniquement les données variables (supposées comme telles)
|
|
void BielletteQ::Ecriture_base_info(ofstream& sort,const int cas)
|
|
{// tout d'abord appel de l'écriture de la classe elem_meca
|
|
ElemMeca::Ecriture_bas_inf(sort,cas);
|
|
// traitement du cas particulier de la BielletteQ
|
|
switch (cas)
|
|
{ case 1 : // ------- on sauvegarde tout -------------------------
|
|
{
|
|
// des tenseurs déformation et contrainte,
|
|
lesPtMecaInt.Ecriture_base_info(sort,cas);
|
|
// les données spécifiques
|
|
sort << "\n epaisStockeDansElement " << " section0= " << donnee_specif.secti.section0
|
|
<< " section_t= " << donnee_specif.secti.section_t
|
|
<< " section_tdt= " << donnee_specif.secti.section_tdt
|
|
<< " variation= " << donnee_specif.variation_section;
|
|
break;
|
|
}
|
|
case 2 : // ----------- sauvegarde uniquement de se qui varie --------------------
|
|
{ // des tenseurs déformation et contrainte,
|
|
lesPtMecaInt.Ecriture_base_info(sort,cas);
|
|
sort << "\n section_tdt= " << donnee_specif.secti.section_tdt;
|
|
break;
|
|
}
|
|
default :
|
|
{ cout << "\nErreur : valeur incorrecte du type d'écriture !\n";
|
|
cout << "BielletteQ::Ecriture_base_info(ofstream& sort,int cas)"
|
|
<< " cas= " << cas << endl;
|
|
Sortie(1);
|
|
}
|
|
}
|
|
};
|
|
|
|
|
|
// Calcul de la matrice géométrique et initiale
|
|
ElemMeca::MatGeomInit BielletteQ::MatricesGeometrique_Et_Initiale (const ParaAlgoControle & pa)
|
|
{ BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
Tableau <TenseurBB *>& d_epsBB = (doCo->d_epsBB);// "
|
|
Tableau <TenseurHH *>& d_sigHH = (doCo->d_sigHH);// "
|
|
bool cald_Dvirtuelle = false;
|
|
if (unefois.CalimpPrem == 0)
|
|
{ unefois.CalimpPrem = 1;
|
|
Tableau<Enum_variable_metrique> tab(20);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igiB_tdt; tab(4) = igijBB_0;
|
|
tab(5) = igijBB_t;tab(6) = igijBB_tdt; tab(7) = igijHH_tdt; tab(8) = id_giB_tdt;
|
|
tab(9) = id_gijBB_tdt ;tab(10) = igiH_tdt;tab(11) = id_giH_tdt;
|
|
tab(12) = id_gijHH_tdt;tab(13) = id_jacobien_tdt;tab(14) = id2_gijBB_tdt;
|
|
tab(15) = id_gijBB_t ;tab(16) = id_gijBB_tdt ;tab(17) = idMtdt ;
|
|
tab(18) = igradVBB_tdt; tab(19) = iVtdt; tab(20) = idVtdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
// on ne calcul la dérivée de la déformation virtuelle qu'une fois
|
|
// car elle est constante dans le temps et indépendante des coordonnées
|
|
cald_Dvirtuelle=true;
|
|
};
|
|
// Par simplicité
|
|
Mat_pleine & matGeom = doCo->matGeom;
|
|
Mat_pleine & matInit = doCo->matInit;
|
|
// mise à zéro de la matrice géométrique
|
|
matGeom.Initialise();
|
|
Vecteur poids =(doCo->segment).TaWi();
|
|
// poids(1) *= donnee_specif.section;
|
|
|
|
// prise en compte de la section sauf dans le cas d'une loi rien
|
|
double section_moyenne = 0.; // où on met tout à 0
|
|
// et surtout on ne recalcule pas la section car on n'a pas de compressibilité avec une loi rien !!
|
|
if (!Loi_rien(loiComp->Id_comport()))
|
|
{const bool atdt=true;
|
|
section_moyenne = CalSectionMoyenne_et_vol_pti(atdt);
|
|
};
|
|
|
|
poids(1) *= section_moyenne;
|
|
ElemMeca::Cal_matGeom_Init
|
|
(matGeom,matInit,doCo->tab_ddl, d_epsBB,
|
|
doCo->d2_epsBB,d_sigHH,1,poids,pa,cald_Dvirtuelle);
|
|
return MatGeomInit(&matGeom,&matInit);
|
|
} ;
|
|
|
|
// retourne les tableaux de ddl associés aux noeuds, gere par l'element
|
|
// ce tableau et specifique a l'element
|
|
const DdlElement & BielletteQ::TableauDdl() const
|
|
{ return doCo->tab_ddl; };
|
|
// liberation de la place pointee
|
|
void BielletteQ::Libere ()
|
|
{Element::Libere (); // liberation de residu et raideur
|
|
LibereTenseur() ; // liberation des tenseur intermediaires
|
|
};
|
|
|
|
// acquisition ou modification d'une loi de comportement
|
|
void BielletteQ::DefLoi (LoiAbstraiteGeneral * NouvelleLoi)
|
|
{ // verification du type de loi
|
|
if ((NouvelleLoi->Dimension_loi() != 1) && (NouvelleLoi->Dimension_loi() != 4))
|
|
{ cout << "\n Erreur, la loi de comportement a utiliser avec des BielletteQs";
|
|
cout << " doit etre de type 1D, \n ici est de type = "
|
|
<< (NouvelleLoi->Dimension_loi()) << " !!! " << endl;
|
|
Sortie(1);
|
|
}
|
|
// cas d'une loi mécanique
|
|
if (GroupeMecanique(NouvelleLoi->Id_categorie()))
|
|
{loiComp = (Loi_comp_abstraite *) NouvelleLoi;
|
|
// initialisation du stockage particulier, ici en fonction du nb de pt d'integ
|
|
int imax = tabSaveDon.Taille();
|
|
for (int i=1;i<= imax;i++) tabSaveDon(i) = loiComp->New_et_Initialise();
|
|
// idem pour le type de déformation mécanique associé
|
|
int iDefmax = tabSaveDefDon.Taille();
|
|
for (int i=1;i<= iDefmax;i++) tabSaveDefDon(i) = def->New_et_Initialise();
|
|
// définition du type de déformation associé à la loi
|
|
loiComp->Def_type_deformation(*def);
|
|
// on active les données particulières nécessaires au fonctionnement de la loi de comp
|
|
loiComp->Activation_donnees(tab_noeud,dilatation,lesPtMecaInt);
|
|
};
|
|
// cas d'une loi thermo physique
|
|
if (GroupeThermique(NouvelleLoi->Id_categorie()))
|
|
{loiTP = (CompThermoPhysiqueAbstraite *) NouvelleLoi;
|
|
// initialisation du stockage particulier,
|
|
int imax = tabSaveTP.Taille();
|
|
for (int i=1;i<= imax;i++) tabSaveTP(i) = loiTP->New_et_Initialise();
|
|
// on active les données particulières nécessaires au fonctionnement de la loi de comp
|
|
loiTP->Activation_donnees(tab_noeud);
|
|
};
|
|
// cas d'une loi de frottement
|
|
if (GroupeFrottement(NouvelleLoi->Id_categorie()))
|
|
loiFrot = (CompFrotAbstraite *) NouvelleLoi;
|
|
};
|
|
|
|
// test si l'element est complet
|
|
int BielletteQ::TestComplet()
|
|
{ int res = ElemMeca::TestComplet(); // test dans la fonction mere
|
|
if ( donnee_specif.secti.section0 == section_defaut)
|
|
{ cout << "\n la section de la BielletteQ n'est pas defini \n";
|
|
res = 0; }
|
|
if ( tab_noeud(1) == NULL)
|
|
{ cout << "\n les noeuds de la BielletteQ ne sont pas defini \n";
|
|
res = 0; }
|
|
else
|
|
{ int testi =1;
|
|
int posi = Id_nom_ddl("X1") -1;
|
|
int dim = ParaGlob::Dimension();
|
|
for (int i =1; i<= dim; i++)
|
|
for (int j=1;j<=nombre_V.nbne;j++)
|
|
if(!(tab_noeud(j)->Existe_ici(Enum_ddl(posi+i))))
|
|
testi = 0;
|
|
if(testi == 0)
|
|
{ cout << "\n les ddls X1,X2 etc des noeuds de la BielletteQ ne sont pas defini \n";
|
|
cout << " \n utilisez BielletteQ::ConstTabDdl() pour completer " ;
|
|
res = 0; }
|
|
}
|
|
return res;
|
|
};
|
|
|
|
// ajout du tableau de ddl des noeuds de BielletteQ
|
|
void BielletteQ::ConstTabDdl()
|
|
{
|
|
Tableau <Ddl> ta(ParaGlob::Dimension());
|
|
int posi = Id_nom_ddl("X1") -1;
|
|
int dim = ParaGlob::Dimension();
|
|
for (int i =1; i<= dim; i++)
|
|
{Ddl inter((Enum_ddl(i+posi)),0.,LIBRE);
|
|
ta(i) = inter;
|
|
}
|
|
// attribution des ddls aux noeuds
|
|
for (int ne=1; ne<= nombre_V.nbne; ne++)
|
|
tab_noeud(ne)->PlusTabDdl(ta);
|
|
};
|
|
|
|
// procesure permettant de completer l'element apres
|
|
// sa creation avec les donnees du bloc transmis
|
|
// peut etre appeler plusieurs fois
|
|
Element* BielletteQ::Complete(BlocGen & bloc,LesFonctions_nD* lesFonctionsnD)
|
|
{ // complétion avec bloc
|
|
if (bloc.Nom(1) == "sections")
|
|
{ // 2 cas: soit la section est fixe, soit elle dépend d'une fct nD
|
|
if (bloc.Nom(2) == "_") // cas fixe
|
|
{donnee_specif.secti.section0 = bloc.Val(1);
|
|
// on initialise aussi les grandeurs à t et tdt
|
|
donnee_specif.secti.section_tdt = donnee_specif.secti.section_t = bloc.Val(1);
|
|
}
|
|
else // cas fct nD, on met à 0 la section pour être sûr qu'il faut pas s'en servir tel quel
|
|
{donnee_specif.secti.section0 = 0.;
|
|
// on initialise aussi les grandeurs à t et tdt
|
|
donnee_specif.secti.section_tdt = donnee_specif.secti.section_t = 0.;
|
|
// on s'occupe de la fct nD
|
|
if (lesFonctionsnD->Existe(bloc.Nom(2)))
|
|
{donnee_specif.fctnD_section = lesFonctionsnD->Trouve(bloc.Nom(2));}
|
|
else // sinon c'est un problème
|
|
{ cout << "\n *** erreur dans la definition de la section via la fonction nD "
|
|
<< bloc.Nom(2) << ", a priori cette fonction n'est pas disponible!! ";
|
|
Sortie(1);
|
|
};
|
|
};
|
|
return this;
|
|
}
|
|
else if (bloc.Nom(1) == "variation_section")
|
|
{ donnee_specif.variation_section = bloc.Val(1);
|
|
// on initialise aussi les grandeurs à t et tdt
|
|
return this;
|
|
}
|
|
else
|
|
return ElemMeca::Complete_ElemMeca(bloc,lesFonctionsnD);
|
|
};
|
|
|
|
// Compléter pour la mise en place de la gestion de l'hourglass
|
|
Element* BielletteQ::Complet_Hourglass(LoiAbstraiteGeneral * loiHourglass, const BlocGen & bloc)
|
|
{ // on initialise le traitement de l'hourglass
|
|
string str_precision; // string vide indique que l'on veut utiliser un élément normal
|
|
BielletteQ::DonneeCommune* doCo = unefois.doCoMemb;
|
|
ElemMeca::Init_hourglass_comp(*(doCo->segmentHourg),str_precision,loiHourglass,bloc);
|
|
// dans le cas où l'hourglass a été activé mais que l'élément n'a pas
|
|
// de traitement particulier associé, alors on désactive l'hourglass
|
|
if ( ((type_stabHourglass == STABHOURGLASS_PAR_COMPORTEMENT) || (type_stabHourglass == STABHOURGLASS_PAR_COMPORTEMENT_REDUIT))
|
|
&&(doCo->segmentHourg == NULL))
|
|
type_stabHourglass = STABHOURGLASS_NON_DEFINIE;
|
|
return this;
|
|
};
|
|
|
|
|
|
// affichage dans la sortie transmise, des variables duales "nom"
|
|
// dans le cas ou nom est vide, affichage de "toute" les variables
|
|
void BielletteQ::AfficheVarDual(ofstream& sort, Tableau<string>& nom)
|
|
{ // affichage de l'entête de l'element
|
|
sort << "\n******************************************************************";
|
|
sort << "\n Element bieletteQ (3 noeuds 2 point d'integration) ";
|
|
sort << "\n******************************************************************";
|
|
// appel de la procedure de elem meca
|
|
if ((unefois.dualSortbiel == 0) && (unefois.CalimpPrem != 0))
|
|
{ VarDualSort(sort,nom,1,1);
|
|
unefois.dualSortbiel = 1;
|
|
}
|
|
else if ((unefois.dualSortbiel == 1) && (unefois.CalimpPrem != 0))
|
|
VarDualSort(sort,nom,1,11);
|
|
else if ((unefois.dualSortbiel == 0) && (unefois.CalResPrem_tdt != 0))
|
|
{ VarDualSort(sort,nom,1,2);
|
|
unefois.dualSortbiel = 1;
|
|
}
|
|
else if ((unefois.dualSortbiel == 1) && (unefois.CalResPrem_tdt != 0))
|
|
VarDualSort(sort,nom,1,12);
|
|
// sinon on ne fait rien
|
|
};
|
|
|
|
// récupération des valeurs au numéro d'ordre = iteg pour
|
|
// les grandeur enu
|
|
// absolue: indique si oui ou non on sort les tenseurs dans la base absolue ou une base particulière
|
|
Tableau <double> BielletteQ::Valeur_a_diff_temps(bool absolue,Enum_dure enu_t,const List_io<Ddl_enum_etendu>& enu,int iteg)
|
|
{ // appel de la procedure de elem meca
|
|
int cas;
|
|
if ((unefois.dualSortbiel == 0) && (unefois.CalimpPrem != 0))
|
|
{ cas=1;unefois.dualSortbiel = 1;
|
|
}
|
|
else if ((unefois.dualSortbiel == 1) && (unefois.CalimpPrem != 0))
|
|
{ cas = 11;}
|
|
else if ((unefois.dualSortbiel == 0) && (unefois.CalResPrem_tdt != 0))
|
|
{ cas=2;unefois.dualSortbiel = 1;
|
|
}
|
|
else if ((unefois.dualSortbiel == 1) && (unefois.CalResPrem_tdt != 0))
|
|
{ cas = 12;}
|
|
// sinon pour l'instant pb, car il faut définir des variable dans la métrique
|
|
else
|
|
{ cout << "\n warning: les grandeurs ne sont pas calculees : il faudrait au moins un pas de calcul"
|
|
<< " pour inialiser les conteneurs des tenseurs resultats ";
|
|
if (ParaGlob::NiveauImpression() >= 4)
|
|
cout << "\n cas non prévu, unefois.dualSortbiel= " << unefois.dualSortbiel
|
|
<< " unefois.CalimpPrem= " << unefois.CalimpPrem
|
|
<< "\n BielletteQ::Valeur_a_diff_temps(Enum_dure enu_t...";
|
|
Sortie(1);
|
|
};
|
|
return ElemMeca::Valeur_multi(absolue,enu_t,enu,iteg,cas);
|
|
};
|
|
|
|
// récupération des valeurs au numéro d'ordre = iteg pour les grandeurs enu
|
|
// ici il s'agit de grandeurs tensorielles, le retour s'effectue dans la liste
|
|
// de conteneurs quelconque associée
|
|
void BielletteQ::ValTensorielle_a_diff_temps(bool absolue,Enum_dure enu_t,List_io<TypeQuelconque>& enu,int iteg)
|
|
{ // appel de la procedure de elem meca
|
|
int cas;
|
|
if ((unefois.dualSortbiel == 0) && (unefois.CalimpPrem != 0))
|
|
{ cas=1;unefois.dualSortbiel = 1;
|
|
}
|
|
else if ((unefois.dualSortbiel == 1) && (unefois.CalimpPrem != 0))
|
|
{ cas = 11;}
|
|
else if ((unefois.dualSortbiel == 0) && (unefois.CalResPrem_tdt != 0))
|
|
{ cas=2;unefois.dualSortbiel = 1;
|
|
}
|
|
else if ((unefois.dualSortbiel == 1) && (unefois.CalResPrem_tdt != 0))
|
|
{ cas = 12;}
|
|
// sinon pour l'instant pb, car il faut définir des variable dans la métrique
|
|
else
|
|
{ cout << "\n warning: les grandeurs ne sont pas calculees : il faudrait au moins un pas de calcul"
|
|
<< " pour inialiser les conteneurs des tenseurs resultats ";
|
|
if (ParaGlob::NiveauImpression() >= 4)
|
|
cout << "\n cas non prévu, unefois.dualSortbiel= " << unefois.dualSortbiel
|
|
<< " unefois.CalimpPrem= " << unefois.CalimpPrem
|
|
<< "\n BielletteQ::ValTensorielle_a_diff_temps(Enum_dure enu_t...";
|
|
Sortie(1);
|
|
};
|
|
ElemMeca::Valeurs_Tensorielles(absolue,enu_t,enu,iteg,cas);
|
|
};
|
|
|
|
// cas d'un chargement volumique,
|
|
// force indique la force volumique appliquée
|
|
// retourne le second membre résultant
|
|
// ici on considère la section de la BielletteQ pour constituer le volume
|
|
// -> explicite à t ou tdt en fonction de la variable booléenne atdt
|
|
Vecteur BielletteQ::SM_charge_volumique_E
|
|
(const Coordonnee& force,Fonction_nD* pt_fonct,bool atdt,const ParaAlgoControle & pa,bool sur_volume_finale_)
|
|
{ BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// dimensionnement de la metrique
|
|
if (!atdt)
|
|
{if(( unefois.CalResPrem_t == 0) && (unefois.CalimpPrem == 0) && (unefois.CalSMvol_t == 0))
|
|
{ unefois.CalSMvol_t = 1;
|
|
Tableau<Enum_variable_metrique> tab(10);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igijBB_0;tab(4) = igijBB_t;
|
|
tab(5) = igijHH_t; tab(6) = id_giB_t; tab(7) = id_gijBB_t ;
|
|
tab(8) = id_gijBB_t ;tab(9) = igradVBB_t; tab(10) = iVt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
}
|
|
else
|
|
{if(( unefois.CalResPrem_tdt == 0) && (unefois.CalimpPrem == 0) && (unefois.CalSMvol_tdt == 0))
|
|
{ unefois.CalSMvol_tdt = 1;
|
|
Tableau<Enum_variable_metrique> tab(11);
|
|
tab(1) = igiB_0; tab(2) = igiB_tdt; tab(3) = igijBB_0;tab(4) = igijBB_tdt;
|
|
tab(5) = igijHH_tdt; tab(6) = id_giB_tdt; tab(7) = id_gijBB_tdt ;
|
|
tab(8) = id_gijBB_t ;tab(9) = id_gijBB_tdt ;
|
|
tab(10) = igradVBB_tdt;tab(11) = iVtdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
};
|
|
// initialisation du résidu
|
|
residu->Zero();
|
|
// appel du programme général d'elemmeca et retour du vecteur second membre
|
|
// multiplié par la section pour avoir le volume, on considère que la section est à jour
|
|
double section = donnee_specif.secti.section_tdt;
|
|
if (!atdt) section = donnee_specif.secti.section_t;
|
|
return (section * ElemMeca::SM_charge_vol_E (doCo->tab_ddl,(doCo->segment).TaPhi()
|
|
,tab_noeud.Taille(),(doCo->segment).TaWi(),force,pt_fonct,pa,sur_volume_finale_));
|
|
};
|
|
|
|
// calcul des seconds membres suivant les chargements
|
|
// cas d'un chargement volumique,
|
|
// force indique la force volumique appliquée
|
|
// retourne le second membre résultant
|
|
// ici on l'épaisseur de l'élément pour constituer le volume -> implicite,
|
|
// pa: permet de déterminer si oui ou non on calcul la contribution à la raideur
|
|
// retourne le second membre et la matrice de raideur correspondant
|
|
Element::ResRaid BielletteQ::SMR_charge_volumique_I
|
|
(const Coordonnee& force,Fonction_nD* pt_fonct,const ParaAlgoControle & pa,bool sur_volume_finale_)
|
|
{ BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
|
|
// initialisation du résidu
|
|
residu->Zero();
|
|
// initialisation de la raideur
|
|
raideur->Zero();
|
|
|
|
// -- définition des constantes de la métrique si nécessaire
|
|
// en fait on fait appel aux même éléments que pour le calcul implicite
|
|
if ((unefois.CalimpPrem == 0) && (unefois.CalSMvol_tdt == 0))
|
|
{ unefois.CalSMvol_tdt = 1;
|
|
Tableau<Enum_variable_metrique> tab(20);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igiB_tdt; tab(4) = igijBB_0;
|
|
tab(5) = igijBB_t;tab(6) = igijBB_tdt; tab(7) = igijHH_tdt; tab(8) = id_giB_tdt;
|
|
tab(9) = id_gijBB_tdt ;tab(10) = igiH_tdt;tab(11) = id_giH_tdt;
|
|
tab(12) = id_gijHH_tdt;tab(13) = id_jacobien_tdt;tab(14) = id2_gijBB_tdt;
|
|
tab(15) = id_gijBB_t ;tab(16) = id_gijBB_tdt ;tab(17) = idMtdt ;
|
|
tab(18) = igradVBB_tdt; tab(19) = iVtdt; tab(20) = idVtdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
|
|
// appel du programme général d'elemmeca
|
|
ElemMeca::SMR_charge_vol_I (doCo->tab_ddl,(doCo->segment).TaPhi()
|
|
,tab_noeud.Taille(),(doCo->segment).TaWi(),force,pt_fonct,pa,sur_volume_finale_);
|
|
// prise en compte de la section à t, pour éviter des oscillations de convergence
|
|
// peut-être à changer par la suite ??
|
|
(*residu) *= donnee_specif.secti.section_t;
|
|
(*raideur) *= donnee_specif.secti.section_t;
|
|
Element::ResRaid el;
|
|
el.res = residu;
|
|
el.raid = raideur;
|
|
return el;
|
|
};
|
|
|
|
// cas d'un chargement lineique, sur l'aretes de la BielletteQ
|
|
// force indique la force lineique appliquée
|
|
// numarete indique le numéro de l'arete chargée ici 1
|
|
// retourne le second membre résultant
|
|
// -> explicite à t ou tdt en fonction de la variable booléenne atdt
|
|
Vecteur BielletteQ::SM_charge_lineique_E(const Coordonnee& force,Fonction_nD* pt_fonct,int ,bool atdt,const ParaAlgoControle & pa)
|
|
{ // initialisation du vecteur résidu
|
|
((*res_extA)(1))->Zero();
|
|
// on récupère ou on crée la frontière arrête
|
|
Frontiere_lineique(1,true);
|
|
BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// dimensionnement de la metrique
|
|
if (!atdt)
|
|
{if( unefois.CalSMlin_t == 0)
|
|
{ unefois.CalSMlin_t = 1;
|
|
Tableau<Enum_variable_metrique> tab(8);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igijBB_0;tab(4) = igijBB_t;
|
|
tab(5) = igijHH_t; tab(6) = id_giB_t; tab(7) = id_gijBB_t ;
|
|
tab(8) = igradVBB_t;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
}
|
|
else
|
|
{if( unefois.CalSMlin_tdt == 0)
|
|
{ unefois.CalSMlin_tdt = 1;
|
|
Tableau<Enum_variable_metrique> tab(8);
|
|
tab(1) = igiB_0; tab(2) = igiB_tdt; tab(3) = igijBB_0;tab(4) = igijBB_tdt;
|
|
tab(5) = igijHH_tdt; tab(6) = id_giB_tdt; tab(7) = id_gijBB_tdt ;
|
|
tab(8) = igradVBB_tdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
};
|
|
// on définit la déformation a doc si elle n'existe pas déjà
|
|
if (defArete(1) == NULL)
|
|
defArete(1) = def; // a priori idem que la BielletteQ
|
|
|
|
// appel du programme général d'elemmeca et retour du vecteur second membre
|
|
return ElemMeca::SM_charge_line_E (doCo->tab_ddl,1,(doCo->segment).TaPhi()
|
|
,tab_noeud.Taille(),(doCo->segment).TaWi(),force,pt_fonct,pa);
|
|
};
|
|
|
|
// cas d'un chargement lineique, sur les aretes frontières des éléments
|
|
// force indique la force lineique appliquée
|
|
// numarete indique le numéro de l'arete chargée ici 1 par défaut
|
|
// retourne le second membre résultant
|
|
// -> implicite,
|
|
// pa: permet de déterminer si oui ou non on calcul la contribution à la raideur
|
|
// retourne le second membre et la matrice de raideur correspondant
|
|
Element::ResRaid BielletteQ::SMR_charge_lineique_I
|
|
(const Coordonnee& force,Fonction_nD* pt_fonct,int ,const ParaAlgoControle & pa)
|
|
{ // initialisation du vecteur résidu et de la raideur
|
|
((*res_extA)(1))->Zero();
|
|
((*raid_extA)(1))->Zero();
|
|
// on récupère ou on crée la frontière arrête
|
|
Frontiere_lineique(1,true);
|
|
|
|
BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// dimensionnement de la metrique
|
|
if( unefois.CalSMRlin == 0)
|
|
{ unefois.CalSMRlin = 1;
|
|
Tableau<Enum_variable_metrique> tab(15);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igiB_tdt; tab(4) = igijBB_0;
|
|
tab(5) = igijBB_t;tab(6) = igijBB_tdt; tab(7) = igijHH_tdt; tab(8) = id_giB_tdt;
|
|
tab(9) = id_gijBB_tdt ;tab(10) = igiH_tdt;tab(11) = id_giH_tdt;
|
|
tab(12) = id_gijHH_tdt;tab(13) = id_jacobien_tdt;tab(14) = id2_gijBB_tdt;
|
|
tab(15) = igradVBB_tdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
// on définit la déformation a doc si elle n'existe pas déjà
|
|
if (defArete(1) == NULL)
|
|
defArete(1) = def; // a priori idem que la BielletteQ
|
|
|
|
// appel du programme général d'elemmeca et retour du vecteur second membre
|
|
return ElemMeca::SMR_charge_line_I (doCo->tab_ddl,1
|
|
,(doCo->segment).TaPhi(),tab_noeud.Taille(),(doCo->segment).TaWi(),force,pt_fonct,pa);
|
|
};
|
|
|
|
|
|
// cas d'un chargement lineique suiveuse, sur l'aretes frontière de la BielletteQ (2D uniquement)
|
|
// force indique la force lineique appliquée
|
|
// numarete indique le numéro de l'arete chargée
|
|
// retourne le second membre résultant
|
|
// -> explicite à t ou tdt en fonction de la variable booléenne atdt
|
|
Vecteur BielletteQ::SM_charge_lineique_Suiv_E(const Coordonnee& force,Fonction_nD* pt_fonct,int ,bool atdt,const ParaAlgoControle & pa)
|
|
{ // initialisation du vecteur résidu
|
|
((*res_extA)(1))->Zero();
|
|
// on récupère ou on crée la frontière arrête
|
|
Frontiere_lineique(1,true);
|
|
BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// dimensionnement de la metrique
|
|
if (!atdt)
|
|
{if( unefois.CalSMlin_t == 0)
|
|
{ unefois.CalSMlin_t = 1;
|
|
Tableau<Enum_variable_metrique> tab(8);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igijBB_0;tab(4) = igijBB_t;
|
|
tab(5) = igijHH_t; tab(6) = id_giB_t; tab(7) = id_gijBB_t ;
|
|
tab(8) = igradVBB_t;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
}
|
|
else
|
|
{if( unefois.CalSMlin_tdt == 0)
|
|
{ unefois.CalSMlin_tdt = 1;
|
|
Tableau<Enum_variable_metrique> tab(8);
|
|
tab(1) = igiB_0; tab(2) = igiB_tdt; tab(3) = igijBB_0;tab(4) = igijBB_tdt;
|
|
tab(5) = igijHH_tdt; tab(6) = id_giB_tdt; tab(7) = id_gijBB_tdt ;
|
|
tab(8) = igradVBB_tdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
};
|
|
// on définit la déformation a doc si elle n'existe pas déjà
|
|
if (defArete(1) == NULL)
|
|
defArete(1) = def; // a priori idem que la BielletteQ
|
|
|
|
// appel du programme général d'elemmeca et retour du vecteur second membre
|
|
return ElemMeca::SM_charge_line_Suiv_E (doCo->tab_ddl,1,(doCo->segment).TaPhi()
|
|
,tab_noeud.Taille(),(doCo->segment).TaWi(),force,pt_fonct,pa);
|
|
};
|
|
|
|
// cas d'un chargement lineique suiveuse, sur l'aretes frontière de la BielletteQ (2D uniquement)
|
|
// force indique la force lineique appliquée
|
|
// numarete indique le numéro de l'arete chargée
|
|
// retourne le second membre résultant
|
|
// -> implicite,
|
|
// pa: permet de déterminer si oui ou non on calcul la contribution à la raideur
|
|
// retourne le second membre et la matrice de raideur correspondant
|
|
Element::ResRaid BielletteQ::SMR_charge_lineique_Suiv_I
|
|
(const Coordonnee& force,Fonction_nD* pt_fonct,int ,const ParaAlgoControle & pa)
|
|
{ // initialisation du vecteur résidu et de la raideur
|
|
((*res_extA)(1))->Zero();
|
|
((*raid_extA)(1))->Zero();
|
|
// on récupère ou on crée la frontière arrête
|
|
Frontiere_lineique(1,true);
|
|
|
|
BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// dimensionnement de la metrique
|
|
if( unefois.CalSMRlin == 0)
|
|
{ unefois.CalSMRlin = 1;
|
|
Tableau<Enum_variable_metrique> tab(15);
|
|
tab(1) = igiB_0; tab(2) = igiB_t; tab(3) = igiB_tdt; tab(4) = igijBB_0;
|
|
tab(5) = igijBB_t;tab(6) = igijBB_tdt; tab(7) = igijHH_tdt; tab(8) = id_giB_tdt;
|
|
tab(9) = id_gijBB_tdt ;tab(10) = igiH_tdt;tab(11) = id_giH_tdt;
|
|
tab(12) = id_gijHH_tdt;tab(13) = id_jacobien_tdt;tab(14) = id2_gijBB_tdt;
|
|
tab(15) = igradVBB_tdt;
|
|
doCo->met_BielletteQ.PlusInitVariables(tab) ;
|
|
};
|
|
// on définit la déformation a doc si elle n'existe pas déjà
|
|
if (defArete(1) == NULL)
|
|
defArete(1) = def; // a priori idem que la BielletteQ
|
|
|
|
// appel du programme général d'elemmeca et retour du vecteur second membre
|
|
return ElemMeca::SMR_charge_line_Suiv_I (doCo->tab_ddl,1
|
|
,(doCo->segment).TaPhi(),tab_noeud.Taille(),(doCo->segment).TaWi(),force,pt_fonct,pa);
|
|
};
|
|
|
|
// cas d'un chargement surfacique hydro-dynamique,
|
|
// Il y a trois forces: une suivant la direction de la vitesse: de type traînée aerodynamique
|
|
// Fn = poids_volu * fn(V) * S * (normale*u) * u, u étant le vecteur directeur de V (donc unitaire)
|
|
// une suivant la direction normale à la vitesse de type portance
|
|
// Ft = poids_volu * ft(V) * S * (normale*u) * w, w unitaire, normal à V, et dans le plan n et V
|
|
// une suivant la vitesse tangente de type frottement visqueux
|
|
// T = to(Vt) * S * ut, Vt étant la vitesse tangentielle et ut étant le vecteur directeur de Vt
|
|
// coef_mul: est un coefficient multiplicateur global (de tout)
|
|
// retourne le second membre résultant
|
|
// // -> explicite à t ou à tdt en fonction de la variable booléenne atdt
|
|
Vecteur BielletteQ::SM_charge_hydrodynamique_E( Courbe1D* frot_fluid,const double& poidvol
|
|
, Courbe1D* coef_aero_n,int num,const double& coef_mul
|
|
, Courbe1D* coef_aero_t,bool atdt,const ParaAlgoControle & pa)
|
|
{ // initialisation du vecteur résidu
|
|
((*res_extN)(num))->Zero(); // ici les frontières sont des points
|
|
BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// on récupère ou on crée la frontière points
|
|
Frontiere_points(num,true);
|
|
// ici on fait l'hypothèse d'une section circulaire, n'ayant pas d'autres infos pour l'instant sur les
|
|
// dimensions de la section, d'où la surface frontale approximativement s'en déduie
|
|
double section = donnee_specif.secti.section_tdt;
|
|
if (!atdt) section = donnee_specif.secti.section_t;
|
|
double dimension_transversale = sqrt(section/ConstMath::Pi);
|
|
// on utilise l'élément géométique de l'élément, pour TaPhi et TaWi pour les passages de paramètres,
|
|
// mais ces infos ne sont pas utilisées car pour les frontières points, pas de points d'integ
|
|
// appel du programme général d'ElemMeca et retour du vecteur second membre
|
|
return (dimension_transversale*ElemMeca::SM_charge_hydrodyn_E (poidvol,doCo->segment.TaPhi(),1
|
|
,frot_fluid,doCo->segment.TaWi()
|
|
,coef_aero_n,num,coef_mul,coef_aero_t,pa,atdt));
|
|
};
|
|
|
|
// -> implicite,
|
|
// pa: permet de déterminer si oui ou non on calcul la contribution à la raideur
|
|
// retourne le second membre et la matrice de raideur correspondant
|
|
Element::ResRaid BielletteQ::SMR_charge_hydrodynamique_I( Courbe1D* frot_fluid,const double& poidvol
|
|
, Courbe1D* coef_aero_n,int num,const double& coef_mul
|
|
, Courbe1D* coef_aero_t,const ParaAlgoControle & pa)
|
|
{ // initialisation du vecteur résidu et de la raideur
|
|
((*res_extN)(num))->Zero();
|
|
((*raid_extN)(num))->Zero();
|
|
BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
// on récupère ou on crée la frontière points
|
|
Frontiere_points(num,true);
|
|
// ici on fait l'hypothèse d'une section circulaire, n'ayant pas d'autres infos pour l'instant sur les
|
|
// dimensions de la section, d'où la surface frontale approximativement s'en déduie
|
|
double section = donnee_specif.secti.section_tdt;
|
|
double dimension_transversale = sqrt(section/ConstMath::Pi);
|
|
// on utilise l'élément géométique de l'élément, pour TaPhi et TaWi pour les passages de paramètres,
|
|
// mais ces infos ne sont pas utilisées car pour les frontières points, pas de points d'integ
|
|
// appel du programme général d'ElemMeca et retour du vecteur second membre
|
|
Element::ResRaid el(ElemMeca::SM_charge_hydrodyn_I (poidvol,doCo->segment.TaPhi(),1
|
|
,frot_fluid,doCo->segment.TaWi(),tabb(posi_tab_front_point+num)->DdlElem()
|
|
,coef_aero_n,num,coef_mul,coef_aero_t,pa));
|
|
(*el.res) *= dimension_transversale;
|
|
(*el.raid) *= dimension_transversale;
|
|
return el;
|
|
};
|
|
|
|
// calcul de la nouvelle section moyenne finale (sans raideur)
|
|
// ramène la section moyenne calculée à atdt
|
|
const double& BielletteQ::CalSectionMoyenne_et_vol_pti(const bool atdt)
|
|
{ // dépend si l'on veut une variation ou non
|
|
if (donnee_specif.variation_section)
|
|
{// S(t+\Delta t) = \frac{(S(t)~\sqrt{g(t)})}{\sqrt{g(t+\Delta t)}}\left ( \frac{K_{t}}{(K_{t} + (P(t+\Delta t) - P(t)))} \right )
|
|
// il faut avoir en tête que P = - trace(sig)/3 !!!
|
|
|
|
|
|
// .. en bouclant sur les pt d'integ enregistré ..
|
|
// -- on récupère et on calcule les jacobiens moyens à t(ou 0) et final
|
|
double jacobien_moy_ini = 0.; double jacobien_moy_fin = 0.; // init
|
|
// -- de même on récupère et on calcul la trace moyenne de la contrainte
|
|
double traceSig_moy = 0.; double traceSig_moy_ini = 0.;
|
|
// -- de même on récupère et on calcul le module de compressibilité moyen
|
|
double troisK_moy = 0.;
|
|
Sect& sect = donnee_specif.secti; // pour simplifier
|
|
for (int i=1;i<= nombre_V.nbi;i++)
|
|
{ // cas de la compressibilité
|
|
const double& troisK = 3. * (*lesPtIntegMecaInterne)(i).ModuleCompressibilite_const();
|
|
troisK_moy += troisK;
|
|
// cas des jacobiens
|
|
// const double& jacobien_0 = *(tabSaveDefDon(i)->Meti_00().jacobien_);
|
|
if (atdt)
|
|
{ const double& jacobien_ini = *(tabSaveDefDon(i)->Meti_t().jacobien_);
|
|
jacobien_moy_ini += jacobien_ini;
|
|
const double& jacobien_fin = *(tabSaveDefDon(i)->Meti_tdt().jacobien_);
|
|
jacobien_moy_fin += jacobien_fin;
|
|
// cas de la trace de sigma
|
|
const double traceSig = (*lesPtIntegMecaInterne)(i).SigHH_const() && (*(tabSaveDefDon(i)->Meti_tdt().gijBB_));
|
|
traceSig_moy += traceSig;
|
|
const double traceSig_ini = (*lesPtIntegMecaInterne)(i).SigHH_t_const() && (*(tabSaveDefDon(i)->Meti_t().gijBB_));
|
|
////---- debug
|
|
//cout << "\n BielletteQ::CalSectionMoyenne_et_vol_pti ";
|
|
//cout << "\n SigHH= "<<((*lesPtIntegMecaInterne)(i).SigHH_const())(1,1)
|
|
// <<" gijBB_tdt= " << (*(tabSaveDefDon(i)->Meti_tdt().gijBB_))(1,1)
|
|
// << " SigHH_t= "<< ((*lesPtIntegMecaInterne)(i).SigHH_t_const())(1,1)
|
|
// << " gijBB_t= "<<(*(tabSaveDefDon(i)->Meti_t().gijBB_))(1,1)
|
|
// << endl;
|
|
////---- fin debug
|
|
traceSig_moy_ini += traceSig_ini;
|
|
(*lesPtIntegMecaInterne)(i).Volume_pti() *= (sect.section_t * jacobien_ini) / jacobien_fin
|
|
* troisK / (troisK - traceSig+traceSig_ini);
|
|
}
|
|
else
|
|
{ const double& jacobien_ini = *(tabSaveDefDon(i)->Meti_00().jacobien_);
|
|
jacobien_moy_ini += jacobien_ini;
|
|
const double& jacobien_fin = *(tabSaveDefDon(i)->Meti_t().jacobien_);
|
|
jacobien_moy_fin += jacobien_fin;
|
|
// cas de la trace de sigma
|
|
const double traceSig = (*lesPtIntegMecaInterne)(i).SigHH_const() && (*(tabSaveDefDon(i)->Meti_t().gijBB_));
|
|
traceSig_moy += traceSig;
|
|
(*lesPtIntegMecaInterne)(i).Volume_pti() *= (sect.section0 * jacobien_ini) / jacobien_fin
|
|
* troisK / (troisK - traceSig);
|
|
};
|
|
};
|
|
jacobien_moy_ini /= nombre_V.nbi;
|
|
jacobien_moy_fin /= nombre_V.nbi;
|
|
traceSig_moy /= nombre_V.nbi;
|
|
troisK_moy /= nombre_V.nbi;
|
|
// d'où le calcul de la nouvelle section en utilisant la relation:
|
|
// (V-V_0)/V = trace(sigma)/3 /K_moy
|
|
if (atdt)
|
|
{ // là on utilise uniquement l'incrément: cf. partie théorique d'herezh
|
|
// S(t+\Delta t) = \frac{(S(t)~\sqrt{g(t)})}{\sqrt{g(t+\Delta t)}}\left ( \frac{K_{t}}{(K_{t} + (P(t+\Delta t) - P(t)))} \right )
|
|
sect.section_tdt = (sect.section_t * jacobien_moy_ini) / jacobien_moy_fin
|
|
* troisK_moy / (troisK_moy - traceSig_moy+traceSig_moy_ini);
|
|
////---- debug
|
|
//cout << "\n BielletteQ::CalSectionMoyenne_et_vol_pti ";
|
|
//cout << "\n jacobien_moy_ini= "<<jacobien_moy_ini<<" jacobien_moy_fin= " << jacobien_moy_fin
|
|
// << " traceSig_moy_ini= "<< traceSig_moy_ini << " traceSig_moy= "<<traceSig_moy
|
|
// << " sect.section_t= "<<sect.section_t << " sect.section_tdt= "<<sect.section_tdt
|
|
// << endl;
|
|
////---- fin debug
|
|
|
|
// ancienne nouvelle formule, via la def logarithmique
|
|
// double ratio = traceSig_moy/troisK_moy;
|
|
// sect.section_tdt = (sect.section0 * jacobien_moy_t) / jacobien_moy_fin
|
|
// * exp(ratio);
|
|
// ancienne formule * troisK_moy / (troisK_moy - traceSig_moy);
|
|
//--debug
|
|
//if (num_elt==1)
|
|
// { Noeud* noe = tab_noeud(1);
|
|
// double R_0 = noe->Coord0().Norme(); double R = noe->Coord2().Norme();
|
|
//cout << "\n e0= " << epais.section_tdt<< " troisK_moy=" << troisK_moy << " traceSig_moy=" << traceSig_moy
|
|
// << " J0= " << jacobien_moy_0 << " J= " << jacobien_moy_fin << " R_0 " << R_0 << " R= " << R;
|
|
// };
|
|
//-- fin debug
|
|
return sect.section_tdt;
|
|
}
|
|
else
|
|
{// dans le cas où on n'a pas d'incrément, on retiend la formule intégrée
|
|
// dans laquelle on suppose que le K ne varie pas ... donc une formule qui n'est
|
|
// pas vraiment valide pour des lois complexes !!
|
|
double ratio = traceSig_moy/troisK_moy;
|
|
sect.section_t = (sect.section0 * jacobien_moy_ini) / jacobien_moy_fin
|
|
* exp(ratio);
|
|
// * troisK_moy / (troisK_moy - traceSig_moy);
|
|
return sect.section_t;
|
|
};
|
|
}
|
|
else
|
|
// si on ne veut pas de variation
|
|
{if (atdt) { return donnee_specif.secti.section_tdt;}
|
|
else { return donnee_specif.secti.section_t;};
|
|
};
|
|
};
|
|
|
|
// retourne la liste abondée de tous les données particulières interne actuellement utilisées
|
|
// par l'élément (actif ou non), sont exclu de cette liste les données particulières des noeuds
|
|
// reliées à l'élément
|
|
// absolue: indique si oui ou non on sort les tenseurs dans la base absolue ou une base particulière
|
|
List_io <TypeQuelconque> BielletteQ::Les_types_particuliers_internes(bool absolue) const
|
|
{ // on commence par récupérer la liste général provenant d'ElemMeca
|
|
List_io <TypeQuelconque> ret = ElemMeca::Les_types_particuliers_internes(absolue);
|
|
// ensuite on va ajouter les données particulières aux sfe
|
|
Grandeur_scalaire_double grand_courant; // def d'une grandeur courante
|
|
// $$$ cas de la section initiale
|
|
TypeQuelconque typQ1(SECTION_MOY_INITIALE,SIG11,grand_courant);
|
|
ret.push_back(typQ1);
|
|
// $$$ cas de la section finale
|
|
TypeQuelconque typQ2(SECTION_MOY_FINALE,SIG11,grand_courant);
|
|
ret.push_back(typQ2);
|
|
return ret;
|
|
};
|
|
|
|
// récupération de grandeurs particulières au numéro d'ordre = iteg
|
|
// celles-ci peuvent être quelconques
|
|
// en retour liTQ est modifié et contiend les infos sur les grandeurs particulières
|
|
// absolue: indique si oui ou non on sort les tenseurs dans la base absolue ou une base particulière
|
|
void BielletteQ::Grandeur_particuliere (bool absolue,List_io<TypeQuelconque>& liTQ,int iteg)
|
|
{
|
|
// on balaie la liste transmise pour les grandeurs propres
|
|
List_io<TypeQuelconque>::iterator il,ilfin = liTQ.end();
|
|
// on commence par appeler la fonction de la classe m่re
|
|
// il n'y aura pas de calcul des grandeurs inactivées
|
|
ElemMeca::Grandeur_particuliere (absolue,liTQ,iteg);
|
|
|
|
// puis les grandeurs sp้cifiques
|
|
for (il=liTQ.begin();il!=ilfin;il++)
|
|
{TypeQuelconque& tipParticu = (*il); // pour simplifier
|
|
|
|
if (tipParticu.EnuTypeQuelconque().Nom_vide()) // veut dire que c'est un enum pur
|
|
switch (tipParticu.EnuTypeQuelconque().EnumTQ())
|
|
{
|
|
// 1) -----cas de la section moyenne initiale, ici elle ne d้pend pas du point d'int้gration
|
|
case SECTION_MOY_INITIALE:
|
|
{ *((Grandeur_scalaire_double*) tipParticu.Grandeur_pointee())=donnee_specif.secti.section0;
|
|
(*il).Active();
|
|
break;
|
|
}
|
|
// 2) -----cas de la section moyenne finale, ici elle ne d้pend pas du point d'int้gration
|
|
case SECTION_MOY_FINALE: // on inactive la grandeur quelconque
|
|
{ *((Grandeur_scalaire_double*) tipParticu.Grandeur_pointee())=donnee_specif.secti.section_tdt;
|
|
(*il).Active();
|
|
break;
|
|
}
|
|
default: // on ne fait rien
|
|
break;
|
|
};
|
|
};
|
|
};
|
|
|
|
|
|
// Calcul des frontieres de l'element
|
|
// creation des elements frontieres et retour du tableau de ces elements
|
|
// la création n'a lieu qu'au premier appel
|
|
// ou lorsque l'on force le paramètre force a true
|
|
// dans ce dernier cas seul les frontière effacées sont recréée
|
|
Tableau <ElFrontiere*> const & BielletteQ::Frontiere(bool force)
|
|
{ int cas = 6; // on veut des lignes et des points
|
|
return Frontiere_elemeca(cas,force);
|
|
};
|
|
|
|
|
|
// =====>>>> methodes privées appelees par les classes dérivees <<<<=====
|
|
|
|
// fonction d'initialisation servant au niveau du constructeur
|
|
BielletteQ::DonneeCommune * BielletteQ::Init
|
|
(Donnee_specif donnee_spec,bool sans_init_noeud)
|
|
{ // bien que la grandeur donnee_specif est défini dans la classe generique
|
|
// le fait de le passer en paramètre permet de tout initialiser dans Init
|
|
// et ceci soit avec les valeurs par défaut soit avec les bonnes valeurs
|
|
donnee_specif =donnee_spec;
|
|
// le fait de mettre les pointeurs a null permet
|
|
// de savoir que l'element n'est pas complet
|
|
tab_noeud.Change_taille(nombre_V.nbne);
|
|
// dans le cas d'un constructeur avec tableau de noeud, il ne faut pas mettre
|
|
// les pointeurs à nuls d'où le test
|
|
if (!sans_init_noeud)
|
|
for (int i =1;i<= nombre_V.nbne;i++) tab_noeud(i) = NULL;
|
|
// definition des donnees communes aux BielletteQxxx
|
|
// a la premiere definition d'une instance
|
|
if (unefois.doCoMemb == NULL)
|
|
BielletteQ::Def_DonneeCommune();
|
|
unefois.doCoMemb = doCo ;
|
|
met = &(doCo->met_BielletteQ); // met est defini dans elemeca
|
|
// def pointe sur la deformation specifique a l'element pour le calcul mecanique
|
|
def = new Deformation(*met,tab_noeud,(doCo->segment).TaDphi(),(doCo->segment).TaPhi());
|
|
// idem pour la remontee aux contraintes et le calcul d'erreur
|
|
defEr = new Deformation(*met,tab_noeud,(doCo->segmentEr).TaDphi(),(doCo->segmentEr).TaPhi());
|
|
// idem pour la remontee aux contraintes et le calcul d'erreur
|
|
defMas = new Deformation(*met,tab_noeud,(doCo->segmentMas).TaDphi(),(doCo->segmentMas).TaPhi());
|
|
// idem pour le calcul de second membre
|
|
defArete.Change_taille(1); // 1 arrête utilisée pour le second membre
|
|
// la déformation sera construite si nécessaire au moment du calcul de second membre
|
|
defArete(1) = NULL;
|
|
|
|
//dimensionnement des deformations et contraintes etc..
|
|
int dimtens = 1;
|
|
lesPtMecaInt.Change_taille_PtIntegMeca(nombre_V.nbi,dimtens);
|
|
// attribution des numéros de référencement dans le conteneur
|
|
for (int ni = 1; ni<= nombre_V.nbi; ni++)
|
|
{lesPtMecaInt(ni).Change_Nb_mail(this->num_maillage);
|
|
lesPtMecaInt(ni).Change_Nb_ele(this->num_elt);
|
|
lesPtMecaInt(ni).Change_Nb_pti(ni);
|
|
};
|
|
|
|
// stockage des donnees particulieres de la loi de comportement au point d'integ
|
|
tabSaveDon.Change_taille(nombre_V.nbi);
|
|
tabSaveTP.Change_taille(nombre_V.nbi);
|
|
tabSaveDefDon.Change_taille(nombre_V.nbi);
|
|
tab_energ.Change_taille(nombre_V.nbi);
|
|
tab_energ_t.Change_taille(nombre_V.nbi);
|
|
// initialisation des pointeurs définis dans la classe Element concernant les résidus et
|
|
// raideur
|
|
// --- cas de la puissance interne ---
|
|
residu = &(doCo->residu_interne); // residu local
|
|
raideur = &(doCo->raideur_interne); // raideur locale
|
|
// --- cas de la dynamique -----
|
|
mat_masse = &(doCo->matrice_masse);
|
|
// --- cas des efforts externes concernant les noeuds ------
|
|
res_extN = &(doCo->residus_externeN); // pour les résidus et second membres
|
|
raid_extN= &(doCo->raideurs_externeN);// pour les raideurs
|
|
// --- cas des efforts externes concernant les aretes ------
|
|
res_extA = &(doCo->residus_externeA); // pour les résidus et second membres
|
|
raid_extA= &(doCo->raideurs_externeA);// pour les raideurs
|
|
|
|
return doCo;
|
|
};
|
|
|
|
|
|
// fonction privee
|
|
// dans cette fonction il ne doit y avoir que les données communes !!
|
|
void BielletteQ::Def_DonneeCommune()
|
|
{ int nbn = nombre_V.nbne;
|
|
// interpollation : element geometrique correspondant aux pt integ et nbn
|
|
GeomSeg seg(nombre_V.nbi,nbn) ;
|
|
// degre de liberte
|
|
int dim = ParaGlob::Dimension();
|
|
DdlElement tab_ddl(nbn,dim);
|
|
int posi = Id_nom_ddl("X1") -1;
|
|
for (int i =1; i<= dim; i++)
|
|
for (int j=1; j<= nbn; j++)
|
|
// tab_ddl (j,i) = Enum_ddl(i+posi);
|
|
tab_ddl.Change_Enum(j,i,Enum_ddl(i+posi));
|
|
// cas des ddl éléments secondaires pour le calcul d'erreur
|
|
// def du nombre de composantes du tenseur de contrainte en absolu
|
|
int nbcomposante = ParaGlob::NbCompTens ();
|
|
DdlElement tab_ddlErr(nbn,nbcomposante);
|
|
posi = Id_nom_ddl("SIG11") -1;
|
|
for (int j=1; j<= nbn; j++)
|
|
{// on definit le nombre de composante de sigma en absolu
|
|
if (nbcomposante == 3)
|
|
// dans l'énumération les composantes ne sont pas a suivre
|
|
{//tab_ddlErr (j,1) = Enum_ddl(1+posi); // cas de SIG11
|
|
//tab_ddlErr (j,2) = Enum_ddl(2+posi); // cas de SIG22
|
|
//tab_ddlErr (j,3) = Enum_ddl(4+posi); // cas de SIG12
|
|
tab_ddlErr.Change_Enum(j,1,Enum_ddl(1+posi)); // cas de SIG11
|
|
tab_ddlErr.Change_Enum(j,2,Enum_ddl(2+posi)); // cas de SIG22
|
|
tab_ddlErr.Change_Enum(j,3,Enum_ddl(4+posi)); // cas de SIG12
|
|
}
|
|
else if (nbcomposante == 6)
|
|
// les composantes sont a suivre dans l'enumération
|
|
for (int i= 1;i<= nbcomposante; i++)
|
|
// tab_ddlErr (j,i) = Enum_ddl(i+posi);
|
|
tab_ddlErr.Change_Enum(j,i,Enum_ddl(i+posi));
|
|
else
|
|
// tab_ddlErr (j,1) = Enum_ddl(1+posi); // uniquement SIG11
|
|
tab_ddlErr.Change_Enum(j,1,Enum_ddl(1+posi)); // uniquement SIG11
|
|
}
|
|
// egalement pour tab_Err1Sig11, def d'un tableau de un ddl : enum SIG11
|
|
// par noeud
|
|
DdlElement tab_Err1Sig11(nbn,DdlNoeudElement(SIG11));
|
|
// toujours pour le calcul d'erreur definition des fonctions d'interpolation
|
|
// pour le calcul du hession de la fonctionnelle :
|
|
GeomSeg segEr(nombre_V.nbiEr,nbn) ;
|
|
// pour le calcul de la matrice masse definition des fonctions d'interpolation
|
|
// en particulier pour le calcul de la masse consistante
|
|
GeomSeg segMa(nombre_V.nbiMas,nbn) ;
|
|
// pour le calcul relatifs à la stabilisation d'hourglass
|
|
GeomSeg* segmentHourg = NULL;
|
|
if (nombre_V.nbiHour > 0)
|
|
segmentHourg = new GeomSeg(nombre_V.nbiHour,nombre_V.nbne);
|
|
// def metrique
|
|
// on definit les variables a priori toujours utiles
|
|
Tableau<Enum_variable_metrique> tab(24);
|
|
tab(1) = iM0; tab(2) = iMt; tab(3) = iMtdt ;
|
|
tab(4) = igiB_0; tab(5) = igiB_t; tab(6) = igiB_tdt;
|
|
tab(7) = igiH_0; tab(8) = igiH_t; tab(9) = igiH_tdt ;
|
|
tab(10)= igijBB_0; tab(11)= igijBB_t; tab(12)= igijBB_tdt;
|
|
tab(13)= igijHH_0; tab(14)= igijHH_t; tab(15)= igijHH_tdt ;
|
|
tab(16)= id_gijBB_tdt; tab(17)= id_giH_tdt; tab(18)= id_gijHH_tdt;
|
|
tab(19)= idMtdt ; tab(20)= id_jacobien_tdt;tab(21)= id2_gijBB_tdt;
|
|
tab(22)= igradVBB_tdt; tab(23) = iVtdt; tab(24)= idVtdt;
|
|
// dim du pb , nb de vecteur de la base , tableau de ddl et la def de variables
|
|
Met_biellette metri(ParaGlob::Dimension(),tab_ddl,tab,nbn) ;
|
|
// ---- cas du calcul d'erreur sur sigma ou epsilon
|
|
// les tenseurs sont exprimees en absolu donc nombre de composante fonction
|
|
// de la dimension absolue
|
|
Tableau <Vecteur *> resEr(nbcomposante);
|
|
for (int i = 1;i<= nbcomposante; i++)
|
|
resEr(i)=new Vecteur (nbn); // une composante par noeud
|
|
Mat_pleine raidEr(nbn,nbn); // la raideur pour l'erreur
|
|
// dimensionnement des différents résidus et raideurs pour le calcul mécanique
|
|
int nbddl = tab_ddl.NbDdl();
|
|
Vecteur residu_int(nbddl); Mat_pleine raideur_int(nbddl,nbddl);
|
|
// cas de la dynamique
|
|
Mat_pleine matmasse(1,nbddl); // a priori on dimensionne en diagonale
|
|
// il y a deux extrémités mais identiques
|
|
Tableau <Vecteur* > residus_extN(2); residus_extN(1) = new Vecteur(dim);
|
|
residus_extN(2) = residus_extN(1);
|
|
int nbddlA = nombre_V.nbneA * dim; int nbA = 1; // 1 arêtes
|
|
Tableau <Vecteur* > residus_extA(nbA); residus_extA(1) = new Vecteur(nbddlA);
|
|
Tableau <Mat_pleine* > raideurs_extA(nbA); raideurs_extA(1) = new Mat_pleine(nbddlA,nbddlA);
|
|
Tableau <Mat_pleine* > raideurs_extN(2);raideurs_extN(1) = new Mat_pleine(dim,dim);
|
|
raideurs_extN(2) = raideurs_extN(1);
|
|
// definition de la classe static contenant toute les variables communes aux BielletteQs
|
|
doCo = new DonneeCommune(seg,tab_ddl,tab_ddlErr,tab_Err1Sig11,metri,resEr,raidEr,segEr,
|
|
residu_int,raideur_int,residus_extN,raideurs_extN,residus_extA,raideurs_extA
|
|
,matmasse,segMa,nombre_V.nbi,segmentHourg);
|
|
};
|
|
|
|
// destructions de certaines grandeurs pointées, créées au niveau de l'initialisation
|
|
void BielletteQ::Destruction()
|
|
{ // tout d'abord l'idée est de détruire certaines grandeurs pointées que pour le dernier élément
|
|
if ((unefois.nbelem_in_Prog == 0)&& (unefois.doCoMemb != NULL))
|
|
// cas de la destruction du dernier élément
|
|
{ BielletteQ::DonneeCommune* doCo = unefois.doCoMemb; // pour simplifier l'écriture
|
|
int resErrTaille = doCo->resErr.Taille();
|
|
for (int i=1;i<= resErrTaille;i++)
|
|
delete doCo->resErr(i);
|
|
delete doCo->residus_externeN(1);
|
|
delete doCo->raideurs_externeN(1);
|
|
delete doCo->residus_externeA(1);
|
|
delete doCo->raideurs_externeA(1);
|
|
if (doCo->segmentHourg != NULL) delete doCo->segmentHourg;
|
|
};
|
|
};
|
|
|
|
|
|
|