// FICHIER : Hypo_hooke1D.cc // CLASSE : Hypo_hooke1D // 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) . // // 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 . // // For more information, please consult: . //#include "Debug.h" # include using namespace std; //introduces namespace std #include #include #include "Sortie.h" #include "TypeConsTens.h" #include "ParaGlob.h" #include "ConstMath.h" #include "TypeQuelconqueParticulier.h" #include "CharUtil.h" #include "Hypo_hooke1D.h" // constructeur par défaut à ne pas utiliser Hypo_hooke1D::SaveResulLoi_Hypo1D::SaveResulLoi_Hypo1D() : Kc(0.),Kc_t(0.),f(0.),f_t(0.) ,eps22(0.),eps22_t(0.),eps33(0.),eps33_t(0.) ,eps_cumulBB(),eps_cumulBB_t() { cout << "\n erreur, le constructeur par defaut ne doit pas etre utilise !" << "\n Hypo_hooke1D::SaveResulLoi_Hypo1D::SaveResulLoi_Hypo1D()"; Sortie(1); }; // le constructeur courant Hypo_hooke1D::SaveResulLoi_Hypo1D::SaveResulLoi_Hypo1D (SaveResul* ): Kc(0.),Kc_t(0.),f(0.),f_t(0.) ,eps22(0.),eps22_t(0.),eps33(0.),eps33_t(0.) ,eps_cumulBB(),eps_cumulBB_t() { }; // 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 Hypo_hooke1D::SaveResulLoi_Hypo1D::Lecture_base_info (ifstream& ent,const int cas) { // ici toutes les données sont toujours a priori variables // ou en tout cas pour les méthodes appelées, elles sont gérées par le paramètre: cas string toto; ent >> toto; #ifdef MISE_AU_POINT if (toto != "HYPO_E_ISO_1D") { cout << "\n erreur en lecture du conteneur pour la loi 1D, on a lue: " <> toto >> Kc >> toto >> f >> toto >> eps22 >> toto >> eps33 >> toto >> eps_cumulBB; Kc_t = Kc; f_t=f; eps22_t=eps22;eps33_t=eps33; eps_cumulBB_t=eps_cumulBB; }; // cas donne le niveau de sauvegarde // = 1 : on sauvegarde tout // = 2 : on sauvegarde uniquement les données variables (supposées comme telles) void Hypo_hooke1D::SaveResulLoi_Hypo1D::Ecriture_base_info (ofstream& sort,const int cas) { // ici toutes les données sont toujours a priori variables // ou en tout cas pour les méthodes appelées, elles sont gérées par le paramètre: cas sort << "\n HYPO_E_ISO_1D "; // données de la loi sort << "\n Kc= "<< Kc << " f= " << f << " eps22= "<< eps22 << " eps33= "<< eps33 << " eps_cumulBB= "<NomCourbe() == "_") f_temperature = Courbe1D::New_Courbe1D(*(loi.f_temperature)); if (Kc_temperature != NULL) if (Kc_temperature->NomCourbe() == "_") Kc_temperature = Courbe1D::New_Courbe1D(*(loi.Kc_temperature));; if (f_IIeps != NULL) if (f_IIeps->NomCourbe() == "_") f_IIeps = Courbe1D::New_Courbe1D(*(loi.f_IIeps));; if (Kc_IIeps != NULL) if (Kc_IIeps->NomCourbe() == "_") Kc_IIeps = Courbe1D::New_Courbe1D(*(loi.Kc_IIeps));; // idem pour les fonctions nD if (f_nD != NULL) if (f_nD->NomFonction() == "_") {// comme il s'agit d'une fonction locale on la redéfinie (sinon pb lors du destructeur de loi) string non_courbe("_"); f_nD = Fonction_nD::New_Fonction_nD(*loi.f_nD); }; if (Kc_nD != NULL) if (Kc_nD->NomFonction() == "_") {// comme il s'agit d'une fonction locale on la redéfinie (sinon pb lors du destructeur de loi) string non_courbe("_"); Kc_nD = Fonction_nD::New_Fonction_nD(*loi.Kc_nD); }; }; Hypo_hooke1D::~Hypo_hooke1D () // Destructeur { if (f_temperature != NULL) if (f_temperature->NomCourbe() == "_") delete f_temperature; if (Kc_temperature != NULL) if (Kc_temperature->NomCourbe() == "_") delete Kc_temperature; if (f_IIeps != NULL) if (f_IIeps->NomCourbe() == "_") delete f_IIeps; if (Kc_IIeps != NULL) if (Kc_IIeps->NomCourbe() == "_") delete Kc_IIeps; if (f_nD != NULL) if (f_nD->NomFonction() == "_") delete f_nD; if (Kc_nD != NULL) if (Kc_nD->NomFonction() == "_") delete Kc_nD; }; // Lecture des donnees de la classe sur fichier void Hypo_hooke1D::LectureDonneesParticulieres (UtilLecture * entreePrinc,LesCourbes1D& lesCourbes1D ,LesFonctions_nD& lesFonctionsnD) { string nom; string nom_class_methode("Hypo_hooke1D::LectureDonneesParticulieres"); // lecture de la compressibilité instantanée *(entreePrinc->entree) >> nom ; #ifdef MISE_AU_POINT if (nom != "Kc=") { cout << "\n erreur en lecture de la loi de HYPO_ELAS1D, on attendait Kc= suivi d'un scalaire " << " ou du mot cle Kc_thermo_dependant_ et ensuite une courbe " << " ou du mot cle Kc_fonction_nD: et ensuite une fonction nD "; entreePrinc->MessageBuffer("**erreur1 Hypo_hooke1D::LectureDonneesParticulieres (...**"); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; #endif // gestion d'erreur de syntaxe non permise if(((strstr(entreePrinc->tablcar,"Kc_thermo_dependant_")!=0) || (strstr(entreePrinc->tablcar,"Kc_fonction_nD:")!=0) ) && (strstr(entreePrinc->tablcar,"Kc_avec_compressibilite_externe")!=0) ) { cout << "\n erreur en lecture de la loi de HYPO_ELAS1D, il n'est pas possible d'avoir une" << " compressibilite externe " << " et une compressibilite thermo_dependante ou via une fonction nD, simultanement: " << " c'est incoherent !! "; entreePrinc->MessageBuffer("**erreur1 Hypo_hooke1D::LectureDonneesParticulieres (...**"); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; // on regarde si Kc est thermo dépendant compress_thermophysique=false; // init par défaut if(strstr(entreePrinc->tablcar,"Kc_thermo_dependant_")!=0) { thermo_dependant=true; *(entreePrinc->entree) >> nom; if (nom != "Kc_thermo_dependant_") { cout << "\n erreur en lecture de la thermodependance de Kc, on aurait du lire le mot" << " cle Kc_thermo_dependant_" << " suivi du nom d'une courbe de charge ou de la courbe elle meme "; entreePrinc->MessageBuffer("**erreur2 Hypo_hooke1D::LectureDonneesParticulieres (...**"); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; // lecture de la loi d'évolution de Kc en fonction de la température *(entreePrinc->entree) >> nom; // on regarde si la courbe existe, si oui on récupère la référence if (lesCourbes1D.Existe(nom)) { Kc_temperature = lesCourbes1D.Trouve(nom);} else { // sinon il faut la lire maintenant string non_courbe("_"); Kc_temperature = Courbe1D::New_Courbe1D(non_courbe,Id_Nom_Courbe1D (nom.c_str())); // lecture de la courbe Kc_temperature->LectDonnParticulieres_courbes (non_courbe,entreePrinc); }; // prepa du flot de lecture entreePrinc->NouvelleDonnee(); } // sinon on regarde si la compressibilité sera fourni en externe else if (strstr(entreePrinc->tablcar,"Kc_avec_compressibilite_externe")!=0) {*(entreePrinc->entree) >> nom; compress_thermophysique=true; } // sinon on regarde si le module dépend d'une fonction nD else if(strstr(entreePrinc->tablcar,"Kc_fonction_nD:")!=0) { string nom; string mot_cle2="Kc_fonction_nD:"; // on lit le nom de la fonction string nom_fonct; bool lec = entreePrinc->Lecture_mot_cle_et_string(nom_class_methode,mot_cle2,nom_fonct); if (!lec ) { entreePrinc->MessageBuffer("**erreur02 en lecture** "+mot_cle2); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; // maintenant on définit la fonction if (lesFonctionsnD.Existe(nom_fonct)) {Kc_nD = lesFonctionsnD.Trouve(nom_fonct); } else {// sinon il faut la lire maintenant string non("_"); Kc_nD = Fonction_nD::New_Fonction_nD(non, Id_Nom_Fonction_nD(nom_fonct)); // lecture de la courbe Kc_nD->LectDonnParticulieres_Fonction_nD (non,entreePrinc); // maintenant on vérifie que la fonction est utilisable if (Kc_nD->NbComposante() != 1 ) { cout << "\n erreur en lecture de Kc, la fonction " << nom_fonct << " est une fonction vectorielle a " << Kc_nD->NbComposante() << " composante alors qu'elle devrait etre scalaire ! " << " elle n'est donc pas utilisable !! "; string message("\n**erreur03** \n"+nom_class_methode+"(..."); entreePrinc->MessageBuffer(message); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; }; // on regarde si la fonction nD intègre la température const Tableau & tab_enu = Kc_nD->Tab_enu_etendu(); if (tab_enu.Contient(TEMP)) thermo_dependant=true; // prepa du flot de lecture entreePrinc->NouvelleDonnee(); } // sinon c'est directement le module que l'on lit else { // lecture de la valeur de Kc *(entreePrinc->entree) >> Kc ; }; // on regarde si Kc dépend multiplicativement directement du deuxième invariant if(strstr(entreePrinc->tablcar,"Kc_IIeps_")!=0) { // on vérifie que Kc_avec_compressibilite_externe n'est pas valide car dans ce cas // la dépendance au deuxième invariant de la déformation n'est pas valide if (compress_thermophysique) { cout << "\n erreur on ne peut pas avoir Kc_avec_compressibilite_externe " << " et Kc dependant du deuxieme invariant d'epsilon, "; entreePrinc->MessageBuffer("**erreur5 Hypo_hooke1D::LectureDonneesParticulieres (...**"); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; // si tout est ok on lit *(entreePrinc->entree) >> nom; if (nom != "Kc_IIeps_") { cout << "\n erreur en lecture de la dependance de Kc au deuxieme invariant d'epsilon, " << " on aurait du lire le mot cle Kc_IIeps_" << " suivi du nom d'une courbe de charge ou de la courbe elle meme "; entreePrinc->MessageBuffer("**erreur6 Hypo_hooke1D::LectureDonneesParticulieres (...**"); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; // lecture de la loi d'évolution en fonction du deuxième invariant *(entreePrinc->entree) >> nom; // on regarde si la courbe existe, si oui on récupère la référence if (lesCourbes1D.Existe(nom)) { Kc_IIeps = lesCourbes1D.Trouve(nom); } else { // sinon il faut la lire maintenant string non_courbe("_"); Kc_IIeps = Courbe1D::New_Courbe1D(non_courbe,Id_Nom_Courbe1D (nom.c_str())); // lecture de la courbe Kc_IIeps->LectDonnParticulieres_courbes (non_courbe,entreePrinc); } // prepa du flot de lecture if(strstr(entreePrinc->tablcar,"fin_loi_HYPO_ELAS1D")==0) entreePrinc->NouvelleDonnee(); }; // ---- lecture de la fonction f reliant \dot sig_1^{.1} et D_1^{.1} *(entreePrinc->entree) >> nom; if (nom != "f=") { cout << "\n erreur en lecture de la fonction f , on aurait du lire le mot f="; entreePrinc->MessageBuffer("**erreur3 Hypo_hooke1D::LectureDonneesParticulieres (...**"); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; // on regarde si f est thermo dépendante if(strstr(entreePrinc->tablcar,"f_thermo_dependant_")!=0) { thermo_dependant=true; *(entreePrinc->entree) >> nom; if (nom != "f_thermo_dependant_") { cout << "\n erreur en lecture de la thermodependance de f, on aurait du lire " << " le mot cle f_thermo_dependant_" << " suivi du nom d'une courbe de charge ou de la courbe elle meme "; entreePrinc->MessageBuffer("**erreur3 Hypo_hooke1D::LectureDonneesParticulieres (...**"); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; // lecture de la loi d'évolution en fonction de la température *(entreePrinc->entree) >> nom; // on regarde si la courbe existe, si oui on récupère la référence if (lesCourbes1D.Existe(nom)) { f_temperature = lesCourbes1D.Trouve(nom); } else { // sinon il faut la lire maintenant string non_courbe("_"); f_temperature = Courbe1D::New_Courbe1D(non_courbe,Id_Nom_Courbe1D (nom.c_str())); // lecture de la courbe f_temperature->LectDonnParticulieres_courbes (non_courbe,entreePrinc); }; // prepa du flot de lecture if(strstr(entreePrinc->tablcar,"fin_loi_HYPO_ELAS1D")==0) entreePrinc->NouvelleDonnee(); } // sinon on regarde si la fonction dépend d'une fonction nD else if(strstr(entreePrinc->tablcar,"f_fonction_nD:")!=0) { string nom; string mot_cle2="f_fonction_nD:"; // on lit le nom de la fonction string nom_fonct; bool lec = entreePrinc->Lecture_mot_cle_et_string(nom_class_methode,mot_cle2,nom_fonct); if (!lec ) { entreePrinc->MessageBuffer("**erreur02 en lecture** "+mot_cle2); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; // maintenant on définit la fonction if (lesFonctionsnD.Existe(nom_fonct)) {f_nD = lesFonctionsnD.Trouve(nom_fonct); } else {// sinon il faut la lire maintenant string non("_"); f_nD = Fonction_nD::New_Fonction_nD(non, Id_Nom_Fonction_nD(nom_fonct)); // lecture de la courbe f_nD->LectDonnParticulieres_Fonction_nD (non,entreePrinc); // maintenant on vérifie que la fonction est utilisable if (f_nD->NbComposante() != 1 ) { cout << "\n erreur en lecture de f , la fonction " << nom_fonct << " est une fonction vectorielle a " << f_nD->NbComposante() << " composante alors qu'elle devrait etre scalaire ! " << " elle n'est donc pas utilisable !! "; string message("\n**erreur03** \n"+nom_class_methode+"(..."); entreePrinc->MessageBuffer(message); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; }; // on regarde si la fonction nD intègre la température const Tableau & tab_enu = f_nD->Tab_enu_etendu(); if (tab_enu.Contient(TEMP)) thermo_dependant=true; // prepa du flot de lecture entreePrinc->NouvelleDonnee(); } else { // lecture de f *(entreePrinc->entree) >> f ; }; // on regarde si la fonction f dépend multiplicativement et directement du deuxième invariant if(strstr(entreePrinc->tablcar,"f_IIeps_")!=0) { *(entreePrinc->entree) >> nom; if (nom != "f_IIeps_") { cout << "\n erreur en lecture de la dependance de f au deuxieme invariant d'epsilon, " << " on aurait du lire le mot cle f_IIeps_" << " suivi du nom d'une courbe de charge ou de la courbe elle meme "; entreePrinc->MessageBuffer("**erreur4 Hypo_hooke1D::LectureDonneesParticulieres (...**"); throw (UtilLecture::ErrNouvelleDonnee(-1)); Sortie(1); }; // lecture de la loi d'évolution en fonction du deuxième invariant *(entreePrinc->entree) >> nom; // on regarde si la courbe existe, si oui on récupère la référence if (lesCourbes1D.Existe(nom)) { f_IIeps = lesCourbes1D.Trouve(nom); } else { // sinon il faut la lire maintenant string non_courbe("_"); f_IIeps = Courbe1D::New_Courbe1D(non_courbe,Id_Nom_Courbe1D (nom.c_str())); // lecture de la courbe f_IIeps->LectDonnParticulieres_courbes (non_courbe,entreePrinc); } // prepa du flot de lecture if(strstr(entreePrinc->tablcar,"fin_loi_HYPO_ELAS1D")==0) entreePrinc->NouvelleDonnee(); }; // ---- on regarde ensuite si le type de dérivée est indiqué string toto; if (strstr(entreePrinc->tablcar,"type_derivee")!=NULL) { // lecture du type *(entreePrinc->entree) >> toto >> type_derive; if ((type_derive!=0)&&(type_derive!=-1)&&(type_derive!=1)) { cout << "\n le type de derivee indique pour la loi de Hypo_hooke1D: "<< type_derive << " n'est pas acceptable (uniquement -1 ou 0 ou 1), on utilise le type par defaut (-1)" << " qui correspond a la derivee de Jauman "; type_derive = -1; }; }; // prepa du flot de lecture if(strstr(entreePrinc->tablcar,"fin_loi_HYPO_ELAS1D")==0) entreePrinc->NouvelleDonnee(); // ---- on regarde ensuite s'il y a une restriction sur la contrainte if (strstr(entreePrinc->tablcar,"restriction_traction_compression_")!=NULL) { // lecture du type *(entreePrinc->entree) >> toto >> restriction_traction_compression; if ( (restriction_traction_compression!=0) &&(restriction_traction_compression!=-1) &&(restriction_traction_compression!=1) ) { cout << "\n le type de restriction lu pour la loi de Hypo_hooke1D: "<< restriction_traction_compression << " n'est pas acceptable (uniquement -1 ou 0 ou 1), on utilise le type par defaut (0)" << " qui correspond a 'aucune restriction' "; restriction_traction_compression = 0; }; }; // prepa du flot de lecture if(strstr(entreePrinc->tablcar,"fin_loi_HYPO_ELAS1D")==0) entreePrinc->NouvelleDonnee(); // ---------- appel au niveau de la classe mère Loi_comp_abstraite::Lecture_type_deformation_et_niveau_commentaire (*entreePrinc,lesFonctionsnD); }; // affichage de la loi void Hypo_hooke1D::Affiche() const { cout << " \n loi de comportement HYPO_ELAS 1D "; if (Kc_temperature != NULL) { cout << " coef Kc thermo dependant explicitement " << " courbe Kc=f(T): " << Kc_temperature->NomCourbe() <<" ";} else if(compress_thermophysique) { cout << " compressibilite calculee avec une loi thermophysique " ;} else if (Kc_nD != NULL) {cout << "\n compressibilite calculee avec une fonction nD: "; cout << "Kc_fonction_nD:" << " "; if (Kc_nD->NomFonction() != "_") cout << Kc_nD->NomFonction(); else Kc_nD->Affiche(); } else { cout << " coef Kc= " << Kc ;} if ( Kc_IIeps != NULL) { cout << " compressibilite dependant multiplicativement du deuxieme invariant d'epsilon " << " courbe Kc_final=Kc*f(eps:eps): " << Kc_IIeps->NomCourbe() <<" ";}; if ( f_temperature != NULL) { cout << " fonction f thermo dependante explicite " << " courbe f=f(T): " << f_temperature->NomCourbe() <<" ";} else if (f_nD != NULL) {cout << "\n fonction f calculee avec une fonction nD: "; cout << "f_fonction_nD:" << " "; if (f_nD->NomFonction() != "_") cout << f_nD->NomFonction(); else f_nD->Affiche(); } else { cout << " cisaillement f= " << f ;} if ( f_IIeps != NULL) { cout << " cisaillement dependant multiplicativement du deuxieme invariant d'epsilon " << " courbe f_final=f*f(eps:eps): " << f_IIeps->NomCourbe() <<" ";}; switch (type_derive) { case -1: cout << ", et derivee de Jauman pour la contrainte" << endl;break; case 0: cout << ", et derivee de Lie deux fois covariantes pour la contrainte" << endl; break; case 1: cout << ", et derivee de Lie deux fois contravariantes pour la contrainte" << endl; break; }; // =0 -> pas de restriction // = -1 : traction uniquement autorisée, la compression est mise à 0 // = 1 : compression uniquement autorisée, la traction est mise à 0 switch (restriction_traction_compression) { case -1: cout << ", pas de compression possible " << endl;break; case 0: cout << ", pas de restriction traction-compression " << endl; break; case 1: cout << ", pas de traction possible " << endl; break; }; cout << endl; // appel de la classe mère Loi_comp_abstraite::Affiche_don_classe_abstraite(); }; // affichage et definition interactive des commandes particulières à chaques lois void Hypo_hooke1D::Info_commande_LoisDeComp(UtilLecture& entreePrinc) { ofstream & sort = *(entreePrinc.Commande_pointInfo()); // pour simplifier cout << "\n definition standart (rep o) ou exemples exhaustifs (rep n'importe quoi) ? "; string rep = "_"; // procédure de lecture avec prise en charge d'un retour chariot rep = lect_return_defaut(true,"o"); if (Kc == -ConstMath::tresgrand) { // on initialise à une valeur arbitraire Kc = 10.;} if (f == -ConstMath::trespetit) { // on initialise à une valeur arbitraire f = 0.15;} sort << "\n# ......................... loi de comportement HYPO_ELAS 1D ........................." << "\n# | coef compressibilite | rapport contrainte | type de derivee objective utilisee |" << "\n# | instantane | deformation | pour le calcul de la contrainte |" << "\n# | Kc | f |type_derivee (facultatif) |" << "\n#......................................................................................" << "\n Kc= " << setprecision(8) << Kc << " f= " << setprecision(8) << f << " type_derivee " << type_derive << "\n fin_loi_HYPO_ELAS1D " << endl; if ((rep != "o") && (rep != "O" ) && (rep != "0") ) {sort << "\n# --------- type de derivee ----------" << "\n# le type de derivee est optionnel: = -1 -> derivee de Jauman (valeur par defaut), " << "\n# = 0 -> derivee deux fois covariantes , " << "\n# = 1 -> derivee deux fois contravariantes" << "\n# dans le cas ou l'on veut une valeur differente de la valeur par defaut" << "\n# il faut mettre le mot cle" << "\n# suivi de la valeur -1 ou 0 ou 1" << "\n# "; sort << "\n# --------- Kc calcule via une loi thermophysique ---------" << "\n# pour le module de compressibilite: Kc, il est possible d'indiquer que" << "\n# le module sera fourni " << "\n# par une loi thermophysique (celle associee au meme element), pour ce faire" << "\n# on indique uniquement: " << "\n# Kc= Kc_avec_compressibilite_externe " << "\n# NB: ensuite aucune autre info concernant Kc ne doit etre fournie. " << "\n# "; sort << "\n# ------- dependance explicite a la temperature ----------" << "\n# \n# chacun des 2 parametres materiau Kc et f " << "\n# peut etre remplace par une fonction dependante de la temperature " << "\n# (sauf si le cas Kc_avec_compressibilite_externe est actif, si oui " << "\n# aucune info complementaire ne peut-etre fournie pour Kc) " << "\n# pour ce faire on utilise un mot cle puis une nom de courbe ou la courbe" << "\n# directement comme avec les autre loi de comportement " << "\n# exemple pour Kc: Kc= Kc_thermo_dependant_ courbe3 " << "\n# exemple pour f: f= f_thermo_dependant_ courbe2 " << "\n# IMPORTANT: a chaque fois qu'il y a une thermodependance, il faut passer une ligne " << "\n# apres la description de la grandeur thermodependante, mais pas de passage " << "\n# a la ligne si ce n'est pas thermo dependant" << "\n# " << "\n# -------- dependance a une fonction nD -------------" << "\n# Kc et f peuvent au choix dependre d'une fonction nD. " << "\n# Dans ce cas Kc et f ne doivent pas avoir ete declare thermodependant, " << "\n# mais on peut inclure une thermo-dependance dans la fonction nD, " << "\n# Exemple pour Kc : " << "\n# Kc= Kc_fonction_nD: fonction_1 " << "\n# " << "\n# Exemple pour f : " << "\n# f= f_fonction_nD: fonction_2 " << "\n# " << "\n# La declaration des fonctions suit la meme logique que pour les courbes 1D " << "\n# " << "\n# -------- dependance directe a (Eps:Eps) -------------" << "\n# il est également possible de definir un facteur multiplicatif pour Kc et f " << "\n# fonction explicite du deuxieme invariant de la deformation = (Eps:Eps) " << "\n# = eps_1^1 * eps_1^1 " << "\n# Ceci permet d'obtenir une evolution finale non lineaire dont l'operateur tangent" << "\n# prend en compte les variations de la deformation." << "\n# pour cela on indique le mot cle Kc_IIeps_ suivi du nom d'une courbe" << "\n# suivi d'un passage a la ligne." << "\n# Exemple de declaration d'un Kc thermo dependant et fonction " << "\n# egalement du deuxieme invariant " << "\n# Kc= Kc_thermo_dependant_ courbe2 " << "\n# Kc_IIeps_ courbe5 " << "\n# " << "\n# Exemple de declaration d'un Kc non thermo-dependant, " << "\n# mais fonction du deuxieme invariant" << "\n# Kc= 10000 Kc_IIeps_ courbe4 " << "\n# " << "\n# NB: 1) on peut cumuler une dependance a une fonction nD et une dependance " << "\n# explicite multiplicative a (Eps:Eps)" << "\n# 2) a la fin de la definition de Kc_IIeps_ , il faut passer une ligne " << "\n# " << "\n# de la meme maniere on peut definir un facteur multiplicatif pour f " << "\n# fonction du deuxieme invariant de la deformation = (Eps:Eps) . " << "\n# pour cela on indique le mot cle f_IIeps_ suivi du nom d'une " << "\n# courbe suivi d'un passage a la ligne" << "\n# Exemple de declaration d'un f thermo dependant et fonction egalement " << "\n# du deuxieme invariant " << "\n# f= f_thermo_dependant_ courbe2 " << "\n# f_IIeps_ courbe4 " << "\n# " << "\n# Exemple de declaration d'un f non thermo-dependant, mais fonction " << "\n# du deuxieme invariant " << "\n# f= 1200 f_IIeps_ courbe4 " << "\n# " << "\n# NB: 1) comme pour Kc, on peut cumuler une dependance a une fonction nD et " << "\n# une dependance explicite multiplicative a (Eps:Eps)" << "\n# 2) a la fin de la definition de f_IIeps_ , il faut passer une ligne " << "\n# " << "\n# " << "\n# NB: pour toutes les definition de courbe ou de fct nD il est aussi possible" << "\n# d'introduire directement la courbe a la place de son nom, comme pour " << "\n# toutes les autres lois " << "\n#" << "\n# --------- restriction traction-compression ----------" << "\n# a la suite du type de derivee on peut indiquer de maniere optionnelle " << "\n# une restriction sur la contrainte calculee, via le mot cle: " << "\n# restriction_traction_compression_ " << "\n# suivi d'un scalaire tel que : " << "\n# = 0 -> pas de restriction (valeur par defaut) " << "\n# = -1 : traction uniquement autorisee, la compression est mise a 0 " << "\n# = 1 : compression uniquement autorisée, la traction est mise a 0 " << "\n# " << "\n# ----------- fin des parametres --------" << "\n# la derniere ligne doit contenir uniquement le mot cle: " << "\n# fin_loi_HYPO_ELAS1D " << "\n# "; }; sort << endl; // appel de la classe mère Loi_comp_abstraite::Info_commande_don_LoisDeComp(entreePrinc); }; // test si la loi est complete int Hypo_hooke1D::TestComplet() { int ret = LoiAbstraiteGeneral::TestComplet(); if (!compress_thermophysique &&((Kc_temperature == NULL) && (Kc == -ConstMath::tresgrand)) && (Kc_nD == NULL) ) { cout << " \n le coef de compressibilite instantane n'est pas defini pour la loi " << Nom_comp(id_comp) << '\n'; ret = 0; } if ((f_temperature == NULL) && (f == -ConstMath::trespetit) && (f_nD == NULL) ) { cout << " \n la fonction f n'est pas defini pour la loi " << Nom_comp(id_comp) << '\n'; ret = 0; } return ret; }; //----- lecture écriture de restart ----- // 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 Hypo_hooke1D::Lecture_base_info_loi(ifstream& ent,const int cas,LesReferences& lesRef ,LesCourbes1D& lesCourbes1D ,LesFonctions_nD& lesFonctionsnD) { string nom; if (cas == 1) { ent >> nom; if (nom != "HYPO_ELAS1D") { cout << "\n erreur en lecture de la loi : Hypo_hooke1D, on attendait le mot cle : HYPO_ELAS1D " << "\n Hypo_hooke1D::Lecture_base_info_loi(..."; Sortie(1); }; // ensuite normalement il n'y a pas de pb de lecture puisque c'est écrit automatiquement (sauf en debug) ent >> nom >> type_derive >> nom >> restriction_traction_compression; // --- compressibilité ent >> nom >> compress_thermophysique ; int test; // sert pour le test des courbes if (!compress_thermophysique) { ent >> test; switch (test) {case 1: // cas d'une valeur numérique ent >> Kc; break; case 2: // cas d'une fonction nD {ent >> nom ; if (nom != " Kc_fonction_nD: ") { cout << "\n erreur en lecture de la fonction nD, on attendait " << " Kc_fonction_nD: et on a lue " << nom << "\n Hypo_hooke1D::Lecture_base_info_loi(..."; Sortie(1); }; Kc_nD = lesFonctionsnD.Lecture_pour_base_info(ent,cas,Kc_nD); break; } case 3: // cas d'une fonction_Kc_temperature {ent >> nom; if (nom != " fonction_Kc_temperature ") { cout << "\n erreur en lecture de la fonction Kc_temperature, on attendait " << " fonction_Kc_temperature et on a lue " << nom << "\n Hypo_hooke1D::Lecture_base_info_loi(..."; Sortie(1); }; Kc_temperature = lesCourbes1D.Lecture_pour_base_info(ent,cas,Kc_temperature); break; } default: { cout << "\n erreur en lecture de Kc, on attendait un nombre 1,2 ou 3 " << " et on a lue " << test << "\n Hypo_hooke1D::Lecture_base_info_loi(..."; Sortie(1); }; }; // dépendance multiplicative à IIeps éventuelle ent >> nom >> test; if (!test) { if (Kc_IIeps != NULL) {if (Kc_IIeps->NomCourbe() == "_") delete Kc_IIeps; Kc_IIeps = NULL;}; } else { ent >> nom; Kc_IIeps = lesCourbes1D.Lecture_pour_base_info(ent,cas,Kc_IIeps); }; }; // --- fonction f ent >> nom >> test; switch (test) {case 1: // cas d'une valeur numérique ent >> f; break; case 2: // cas d'une fonction nD {ent >> nom ; if (nom != " f_fonction_nD: ") { cout << "\n erreur en lecture de la fonction nD, on attendait " << " f_fonction_nD: et on a lue " << nom << "\n Hypo_hooke1D::Lecture_base_info_loi(..."; Sortie(1); }; f_nD = lesFonctionsnD.Lecture_pour_base_info(ent,cas,f_nD); break; } case 3: // cas d'une fonction_f_temperature {ent >> nom; if (nom != " fonction_f_temperature ") { cout << "\n erreur en lecture de la fonction f_temperature, on attendait " << " fonction_f_temperature et on a lue " << nom << "\n Hypo_hooke1D::Lecture_base_info_loi(..."; Sortie(1); }; f_temperature = lesCourbes1D.Lecture_pour_base_info(ent,cas,f_temperature); break; } default: { cout << "\n erreur en lecture de f, on attendait un nombre 1,2 ou 3 " << " et on a lue " << test << "\n Hypo_hooke1D::Lecture_base_info_loi(..."; Sortie(1); }; }; // dépendance multiplicative à IIeps éventuelle ent >> nom >> test; if (!test) { if (f_IIeps != NULL) {if (f_IIeps->NomCourbe() == "_") delete f_IIeps; f_IIeps = NULL;}; } else { ent >> nom; f_IIeps = lesCourbes1D.Lecture_pour_base_info(ent,cas,f_IIeps); }; } Loi_comp_abstraite::Lecture_don_base_info(ent,cas,lesRef,lesCourbes1D,lesFonctionsnD); }; // cas donne le niveau de sauvegarde // = 1 : on sauvegarde tout // = 2 : on sauvegarde uniquement les données variables (supposées comme telles) void Hypo_hooke1D::Ecriture_base_info_loi(ofstream& sort,const int cas) { if (cas == 1) { sort << " HYPO_ELAS1D " << " type_derivee " << type_derive << " restriction_traction_compression "<& liTQ,Loi_comp_abstraite::SaveResul * saveDon,list& decal ) const {// tout d'abord on récupère le conteneur SaveResulLoi_Hypo1D & save_resul = *((SaveResulLoi_Hypo1D*) saveDon); // on passe en revue la liste List_io::iterator itq,itqfin=liTQ.end(); list::iterator idecal=decal.begin(); for (itq=liTQ.begin();itq!=itqfin;itq++,idecal++) {TypeQuelconque& tipParticu = (*itq); // pour simplifier if (tipParticu.EnuTypeQuelconque().Nom_vide()) // veut dire que c'est un enum pur { EnumTypeQuelconque enuTQ = tipParticu.EnuTypeQuelconque().EnumTQ(); switch (enuTQ) {// 1) -----cas du module de compressibilité dépendant de la température case COMPRESSIBILITE_TANGENTE: { Tab_Grandeur_scalaire_double& tyTQ= *((Tab_Grandeur_scalaire_double*) (*itq).Grandeur_pointee()); // pour simplifier tyTQ(1+(*idecal))=save_resul.Kc/3.;(*idecal)++; break; } case MODULE_TANGENT_1D: { Tab_Grandeur_scalaire_double& tyTQ= *((Tab_Grandeur_scalaire_double*) (*itq).Grandeur_pointee()); // pour simplifier tyTQ(1+(*idecal))=save_resul.f;(*idecal)++; break; } case DEF_EPAISSEUR: { Tab_Grandeur_scalaire_double& tyTQ= *((Tab_Grandeur_scalaire_double*) (*itq).Grandeur_pointee()); // pour simplifier tyTQ(1+(*idecal))=save_resul.eps33;(*idecal)++; break; } case DEF_LARGEUR: { Tab_Grandeur_scalaire_double& tyTQ= *((Tab_Grandeur_scalaire_double*) (*itq).Grandeur_pointee()); // pour simplifier tyTQ(1+(*idecal))=save_resul.eps22;(*idecal)++; break; } case DEF_ASSO_LOI: { Tab_Grandeur_TenseurBB& tyTQ= *((Tab_Grandeur_TenseurBB*) (*itq).Grandeur_pointee()); // pour simplifier if (ParaGlob::Dimension() != 1) // il faut alors affecter en le 1D à la bonne dim {tyTQ(1+(*idecal)).Affectation_trans_dimension(save_resul.eps_cumulBB,true); } else // cas même dimension {tyTQ(1+(*idecal)) = save_resul.eps_cumulBB; }; (*idecal)++; break; } default: break; // sinon rien à faire }; }; // fin du if }; // fin de la boucle }; // récupération et création de la liste de tous les grandeurs particulières // ces grandeurs sont ajoutées à la liste passées en paramètres // absolue: indique si oui ou non on sort les tenseurs dans la base absolue ou une base particulière void Hypo_hooke1D::ListeGrandeurs_particulieres(bool absolue,List_io& liTQ) const { //on commence par définir une grandeur_scalaire_double Tableau tab_1(1); Tab_Grandeur_scalaire_double grand_courant(tab_1); // def d'un type quelconque représentatif à chaque grandeur // a priori ces grandeurs sont défini aux points d'intégration identique à la contrainte par exemple // enu_ddl_type_pt est définit dans la loi Abtraite générale //on regarde si ce type d'info existe déjà: si oui on augmente la taille du tableau, si non on crée // $$$ cas de MODULE_COMPRESSIBILITE: {List_io::iterator itq,itqfin=liTQ.end(); bool nexistePas = true; for (itq=liTQ.begin();itq!=itqfin;itq++) if ((*itq).EnuTypeQuelconque() == COMPRESSIBILITE_TANGENTE) { Tab_Grandeur_scalaire_double& tyTQ= *((Tab_Grandeur_scalaire_double*) (*itq).Grandeur_pointee()); // pour simplifier int taille = tyTQ.Taille()+1; tyTQ.Change_taille(taille); nexistePas = false; }; if (nexistePas) {TypeQuelconque typQ1(COMPRESSIBILITE_TANGENTE,enu_ddl_type_pt,grand_courant); liTQ.push_back(typQ1); }; }; // $$$ cas du module tangent {List_io::iterator itq,itqfin=liTQ.end(); bool nexistePas = true; for (itq=liTQ.begin();itq!=itqfin;itq++) if ((*itq).EnuTypeQuelconque() == MODULE_TANGENT_1D) { Tab_Grandeur_scalaire_double& tyTQ= *((Tab_Grandeur_scalaire_double*) (*itq).Grandeur_pointee()); // pour simplifier int taille = tyTQ.Taille()+1; tyTQ.Change_taille(taille); nexistePas = false; }; if (nexistePas) {TypeQuelconque typQ2(MODULE_TANGENT_1D,enu_ddl_type_pt,grand_courant); liTQ.push_back(typQ2); }; }; // $$$ cas de la déformation en épaisseur (suivant 3) {List_io::iterator itq,itqfin=liTQ.end(); bool nexistePas = true; for (itq=liTQ.begin();itq!=itqfin;itq++) if ((*itq).EnuTypeQuelconque() == DEF_EPAISSEUR) { Tab_Grandeur_scalaire_double& tyTQ= *((Tab_Grandeur_scalaire_double*) (*itq).Grandeur_pointee()); // pour simplifier int taille = tyTQ.Taille()+1; tyTQ.Change_taille(taille); nexistePas = false; }; if (nexistePas) {TypeQuelconque typQ1(DEF_EPAISSEUR,enu_ddl_type_pt,grand_courant); liTQ.push_back(typQ1); }; }; // $$$ cas de la déformation en largeur (suivant 2) {List_io::iterator itq,itqfin=liTQ.end(); bool nexistePas = true; for (itq=liTQ.begin();itq!=itqfin;itq++) if ((*itq).EnuTypeQuelconque() == DEF_LARGEUR) { Tab_Grandeur_scalaire_double& tyTQ= *((Tab_Grandeur_scalaire_double*) (*itq).Grandeur_pointee()); // pour simplifier int taille = tyTQ.Taille()+1; tyTQ.Change_taille(taille); nexistePas = false; }; if (nexistePas) {TypeQuelconque typQ1(DEF_LARGEUR,enu_ddl_type_pt,grand_courant); liTQ.push_back(typQ1); }; }; // $$$ cas de la déformation cumulée associée à la loi, c-a-d au type d'intégration {List_io::iterator itq,itqfin=liTQ.end(); bool nexistePas = true; for (itq=liTQ.begin();itq!=itqfin;itq++) if ((*itq).EnuTypeQuelconque() == DEF_ASSO_LOI) {Tab_Grandeur_TenseurBB& tyTQ= *((Tab_Grandeur_TenseurBB*) (*itq).Grandeur_pointee()); // pour simplifier int taille = tyTQ.Taille()+1; tyTQ.Change_taille(taille); nexistePas = false; }; if (nexistePas) {int dim_espace = ParaGlob::Dimension(); TenseurBB* tens = NevezTenseurBB(dim_espace); // un tenseur typique Tab_Grandeur_TenseurBB epsassoBB(*tens,1); // def d'un type quelconque représentatif TypeQuelconque typQ(DEF_ASSO_LOI,EPS11,epsassoBB); liTQ.push_back(typQ); delete tens; // car on n'en a plus besoin }; }; }; // calcul d'un module d'young équivalent à la loi double Hypo_hooke1D::Module_young_equivalent(Enum_dure temps,const Deformation & ,SaveResul * saveDon) { // ici le module d'Young correspond à la fonction f, cependant f est le module tangent et non // sécant ... peut-être sera changé par la suite // compte tenu du fait que l'on ne connait pas la métrique etc... on ramène le module en cours // on récupère le conteneur SaveResulLoi_Hypo1D & save_resul = *((SaveResulLoi_Hypo1D*) saveDon); // au niveau instantané on a: E=(3*Kc*mu)/(2*Kc+mu); // et nu = (Kc-mu)/(2*Kc+mu) double E=0.; switch (temps) {case TEMPS_0 : // rien n'a été calculé on utilise les grandeurs initiales E=f; break; case TEMPS_t : // on utilise les grandeurs stockées à t {E= save_resul.f_t; break; } case TEMPS_tdt : // on utilise les grandeurs stockées à tdt {E= save_resul.f; break; } }; return E; }; // récupération d'un module de compressibilité équivalent à la loi, ceci pour un chargement nul // il s'agit ici de la relation -pression = sigma_trace/3. = module de compressibilité * I_eps // >>> en fait ici il s'agit du dernier module tangent calculé !! double Hypo_hooke1D::Module_compressibilite_equivalent (Enum_dure temps,const Deformation & def,SaveResul * saveDon) { // ici le module correspond à Kc, cependant Kc est un module tangent et non // sécant ... peut-être sera changé par la suite // compte tenu du fait que l'on ne connait pas la métrique etc... on ramène le module en cours // on récupère le conteneur SaveResulLoi_Hypo1D & save_resul = *((SaveResulLoi_Hypo1D*) saveDon); double K=0.; switch (temps) {case TEMPS_0 : // rien n'a été calculé on utilise les grandeurs initiales K=Kc/3.; break; case TEMPS_t : // on utilise les grandeurs stockées à t {K=save_resul.Kc_t/3.; break; } case TEMPS_tdt : // on utilise les grandeurs stockées à tdt {K = save_resul.Kc/3.; break; } }; return K; }; // récupération de la variation relative d'épaisseur calculée: h/h0 // cette variation n'est utile que pour des lois en contraintes planes // - pour les lois 3D : retour d'un nombre très grand, indiquant que cette fonction est invalide // - pour les lois 2D def planes: retour de 0 // les infos nécessaires à la récupération , sont stockées dans saveResul // qui est le conteneur spécifique au point où a été calculé la loi double Hypo_hooke1D::HsurH0(SaveResul * saveDon) const { // on récupère le conteneur SaveResulLoi_Hypo1D & save_resul = *((SaveResulLoi_Hypo1D*) saveDon); // la déformation d'épaisseur obtenue dépend du type de dérivée retenue double h_sur_h0 = 0.; // init switch (type_derive) {case 0: // cas d'une dérivée de Lie deux fois covariantes // cas d'une déformation d'Almansi { // epsBB33 = 1/2 * (1. - (h0/h)^2)= 1/2 * (1. - 1./(var_epai)^2), en orthonormee h_sur_h0 = sqrt(1./(1.-2.* save_resul.eps33)); }; break; case -1: // cas d'une dérivée de jauman: // cas d'une def logarithmique cumulée que l'on approxime à une def logarithmique { // eps_33 = log(var_epai); h_sur_h0 = exp(save_resul.eps33); }; break; //**** à finir // case 1: // cas d'une dérivée de Lie deux fois contravariantes // // cas d'une déformation e_{2} // { // epsBB33 = 1/2 * (1. - (h0/h)^2)= 1/2 * (1. - 1./(var_epai)^2), en orthonormee // h_sur_h0 = sqrt(1./(1.-2.* save_resul.eps33)); // }; // break; default : cout << "\nErreur : type de derivee qui n'est pas actuellement pris en compte, type= " << type_derive; cout << "\n Hypo_hooke1D::HsurH0 (... \n"; Sortie(1); }; return h_sur_h0; }; // récupération de la variation relative de largeur calculée: b/b0 // cette variation n'est utile que pour des lois en contraintes planes double // - pour les lois 3D et 2D : retour d'un nombre très grand, indiquant que cette fonction est invalide // les infos nécessaires à la récupération , sont stockées dans saveResul // qui est le conteneur spécifique au point où a été calculé la loi double Hypo_hooke1D::BsurB0(SaveResul * saveDon) const { // on récupère le conteneur SaveResulLoi_Hypo1D & save_resul = *((SaveResulLoi_Hypo1D*) saveDon); // la déformation d'épaisseur obtenue dépend du type de dérivée retenue double b_sur_b0 = 0.; // init switch (type_derive) {case 0: // cas d'une dérivée de Lie deux fois covariantes // cas d'une déformation d'Almansi { // epsBB33 = 1/2 * (1. - (h0/h)^2)= 1/2 * (1. - 1./(var_epai)^2), en orthonormee b_sur_b0 = sqrt(1./(1.-2.* save_resul.eps22)); }; break; case -1: // cas d'une dérivée de jauman: // cas d'une def logarithmique cumulée que l'on approxime à une def logarithmique { // eps_33 = log(var_epai); b_sur_b0 = exp(save_resul.eps22); }; break; //**** à finir // case 1: // cas d'une dérivée de Lie deux fois contravariantes // // cas d'une déformation e_{2} // { // epsBB33 = 1/2 * (1. - (h0/h)^2)= 1/2 * (1. - 1./(var_epai)^2), en orthonormee // b_sur_b0 = sqrt(1./(1.-2.* save_resul.eps33)); // }; // break; default : cout << "\nErreur : type de derivee qui n'est pas actuellement pris en compte, type= " << type_derive; cout << "\n Hypo_hooke1D::HsurH0 (... \n"; Sortie(1); }; return b_sur_b0; }; // ========== codage des METHODES VIRTUELLES protegees:================ // calcul des contraintes a t+dt void Hypo_hooke1D::Calcul_SigmaHH (TenseurHH& sigHH_t,TenseurBB& ,DdlElement & ,TenseurBB & gijBB_t_,TenseurHH & gijHH_t_,BaseB& giB,BaseH& gi_H,TenseurBB& epsBB_ ,TenseurBB& delta_epsBB_, TenseurBB& gijBB_ ,TenseurHH & gijHH_,Tableau & ,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& ex) { #ifdef MISE_AU_POINT if (delta_epsBB_.Dimension() != 1) { cout << "\nErreur : la dimension devrait etre 1 !\n"; cout << " Hypo_hooke1D::Calcul_SigmaHH\n"; Sortie(1); }; #endif const Tenseur1BB & epsBB = *((Tenseur1BB*) &epsBB_); // passage en dim 1 const Tenseur1BB & delta_epsBB = *((Tenseur1BB*) &delta_epsBB_); // passage en dim 1 const Tenseur1HH & gijHH = *((Tenseur1HH*) &gijHH_); // " " " " const Tenseur1BB & gijBB = *((Tenseur1BB*) &gijBB_); // " " " " const Tenseur1HH & gijHH_t = *((Tenseur1HH*) &gijHH_t_); // " " " " const Tenseur1BB & gijBB_t = *((Tenseur1BB*) &gijBB_t_); // " " " " Tenseur1HH & sigHH = *((Tenseur1HH*) &sigHH_); // " " " " Tenseur1HH & sigHH_nn = *((Tenseur1HH*) &sigHH_t); // " " " " SaveResulLoi_Hypo1D & save_resul = *((SaveResulLoi_Hypo1D*) saveResul); // tenseur intermédiaires utilisées selon les cas (par forcément dans tous les cas !!) Tenseur1BH sigBH_n;Tenseur1HH sigHH_n; Tenseur1BB sigBB_n;Tenseur1BB sig_interBB_n; switch (type_derive) {case -1: // cas d'une dérivée de jauman: 1/2 deux fois covariant + deux fois contra {sig_interBB_n = gijBB_t * sigHH_nn * gijBB_t; sigBH_n = 0.5*( sig_interBB_n * gijHH + gijBB * sigHH_nn) ; sigHH_n = gijHH * sigBH_n ; sigBB_n = sigBH_n * gijBB; // pour la déformation cumulée associée, Tenseur1HH eps_interHH_n = gijHH_t * save_resul.eps_cumulBB_t * gijHH_t; Tenseur1BB epsBB_n = 0.5*( gijBB * eps_interHH_n * gijBB + save_resul.eps_cumulBB_t); save_resul.eps_cumulBB = delta_epsBB + epsBB_n; break;} case 0: //cas d'une dérivée de Lie deux fois covariantes {sigBB_n = gijBB_t * sigHH_nn * gijBB_t; sigBH_n = sigBB_n * gijHH ; sigHH_n = gijHH * sigBH_n ; // pour la déformation cumulée associée, on sait que cela donne directement la def d'Almansi // mais la def passée en paramètre pourrait ne pas être d'Almansi, donc on calcule // quand même la def cumulée save_resul.eps_cumulBB = save_resul.eps_cumulBB_t + delta_epsBB; break; } case 1: // cas d'une dérivée de Lie deux fois contravariante {sigHH_n = sigHH_nn; sigBH_n = gijBB * sigHH_n; sigBB_n = sigBH_n * gijBB; // pour la déformation cumulée associée, Tenseur1HH eps_interHH_n = gijHH_t * save_resul.eps_cumulBB_t * gijHH_t; Tenseur1BB epsBB_n = gijBB * eps_interHH_n * gijBB ; save_resul.eps_cumulBB = delta_epsBB + epsBB_n; break; } }; // opération de transmission de la métrique const Met_abstraite::Impli* ex_impli = NULL; const Met_abstraite::Expli_t_tdt* ex_expli_tdt = &ex; const Met_abstraite::Umat_cont* ex_expli = NULL; // cas de la thermo dépendance, on calcul les grandeurs if (f_temperature != NULL) { f = f_temperature->Valeur(*temperature);} else if (f_nD != NULL) // là il faut calculer la fonction nD { // on utilise la méthode générique de loi abstraite list list_save; // inter pour l'appel de la fonction list_save.push_back(saveResul); Tableau & tab_val = Loi_comp_abstraite::Loi_comp_Valeur_FnD_Evoluee (f_nD,1 // une seule valeur attendue en retour ,ex_impli,ex_expli_tdt,ex_expli ,NULL ,NULL ,&list_save ); /* // ici on utilise les variables connues aux pti, ou calculées à partir de // on commence par récupérer les conteneurs des grandeurs à fournir List_io & li_enu_scal = f_nD->Li_enu_etendu_scalaire(); List_io & li_quelc = f_nD->Li_equi_Quel_evolue(); bool absolue = true; // on se place systématiquement en absolu // on va utiliser la méhode Valeur_multi_interpoler_ou_calculer // pour les grandeurs strictement scalaire Tableau val_ddl_enum(Valeur_multi_interpoler_ou_calculer (absolue,TEMPS_tdt,li_enu_scal,ex_impli,ex_expli_tdt,ex_expli,NULL) ); // on utilise la méthode Valeurs_Tensorielles_interpoler_ou_calculer // pour les Coordonnees et Tenseur Valeurs_Tensorielles_interpoler_ou_calculer (absolue,TEMPS_tdt,li_quelc,ex_impli,ex_expli_tdt,ex_expli,NULL); // calcul de la valeur et retour dans tab_ret Tableau & tab_val = f_nD->Valeur_FnD_Evoluee(&val_ddl_enum,&li_enu_scal,&li_quelc,NULL,NULL); #ifdef MISE_AU_POINT if (tab_val.Taille() != 1) { cout << "\nErreur : la fonction nD relative a la fonction f " << " doit calculer un scalaire or le tableau de retour est de taille " << tab_val.Taille() << " ce n'est pas normal !\n"; cout << " Hypo_hooke1D::Calcul_SigmaHH\n"; Sortie(1); }; #endif */ // on récupère le premier élément du tableau uniquement f = tab_val(1); }; if (!compress_thermophysique) // sinon kc est mis à jour avec la méthode CalculGrandeurTravail(..) {if (Kc_temperature != NULL) {Kc = Kc_temperature->Valeur(*temperature);} else if (Kc_nD != NULL) // là il faut calculer la fonction nD { // on utilise la méthode générique de loi abstraite list list_save; // inter pour l'appel de la fonction list_save.push_back(saveResul); Tableau & tab_val = Loi_comp_abstraite::Loi_comp_Valeur_FnD_Evoluee (Kc_nD,1 // une seule valeur attendue en retour ,ex_impli,ex_expli_tdt,ex_expli ,NULL ,NULL ,&list_save ); /* // ici on utilise les variables connues aux pti, ou calculées à partir de // on commence par récupérer les conteneurs des grandeurs à fournir List_io & li_enu_scal = Kc_nD->Li_enu_etendu_scalaire(); List_io & li_quelc = Kc_nD->Li_equi_Quel_evolue(); bool absolue = true; // on se place systématiquement en absolu // on va utiliser la méhode Valeur_multi_interpoler_ou_calculer // pour les grandeurs strictement scalaire Tableau val_ddl_enum(Valeur_multi_interpoler_ou_calculer (absolue,TEMPS_tdt,li_enu_scal,ex_impli,ex_expli_tdt,ex_expli,NULL) ); // on utilise la méthode Valeurs_Tensorielles_interpoler_ou_calculer // pour les Coordonnees et Tenseur Valeurs_Tensorielles_interpoler_ou_calculer (absolue,TEMPS_tdt,li_quelc,ex_impli,ex_expli_tdt,ex_expli,NULL); // calcul de la valeur et retour dans tab_ret Tableau & tab_val = Kc_nD->Valeur_FnD_Evoluee(&val_ddl_enum,&li_enu_scal,&li_quelc,NULL,NULL); #ifdef MISE_AU_POINT if (tab_val.Taille() != 1) { cout << "\nErreur : la fonction nD relative a la fonction Kc " << " doit calculer un scalaire or le tableau de retour est de taille " << tab_val.Taille() << " ce n'est pas normal !\n"; cout << " Hypo_hooke1D::Calcul_SigmaHH\n"; Sortie(1); }; #endif */ // on récupère le premier élément du tableau uniquement Kc = tab_val(1); }; }; // cas d'une non-linéarité multiplicative en fonction de la déformation double Kc_use=Kc;// pour ne pas changer les valeurs à chaque passage !! double f_use=f;// " if ((Kc_IIeps != NULL)||(f_IIeps != NULL)) { Tenseur1BH epsBH = epsBB * gijHH;double II_eps=epsBH.II(); if (Kc_IIeps != NULL) { double coef1 = Kc_IIeps->Valeur(II_eps); Kc_use *= coef1; }; if (f_IIeps != NULL) { double coef2 = f_IIeps->Valeur(II_eps); f_use *= coef2; }; }; // sauvegarde des paramètres matériau save_resul.Kc = Kc_use; save_resul.f = f_use; // calcul de la contrainte unidirectionnelle Tenseur1BH delta_epsBH = delta_epsBB * gijHH; Tenseur1BH delta_sigmaBH(delta_epsBH*f_use); Tenseur1BH sigBH = sigBH_n + delta_sigmaBH; sigHH = gijHH * sigBH; // on applique éventuellement une restriction sur la contrainte // =0 -> pas de restriction // = -1 : traction uniquement autorisée, la compression est mise à 0 // = 1 : compression uniquement autorisée, la traction est mise à 0 switch (restriction_traction_compression) { case 0: break; // cas sans restriction case -1: // traction uniquement autorisée, la compression est mise à 0 {if (sigHH(1,1)<0.) sigHH.Coor(1,1) = 0.; break; } case 1: // compression uniquement autorisée, la traction est mise à 0 {if (sigHH(1,1)>0.) sigHH.Coor(1,1) = 0.; break; } default: cout << "\n erreur Hypo_hooke1D::Calcul_SigmaHH " << " restriction_traction_compression= "<< restriction_traction_compression << " ce qui n'est pas exploitable !! "; Sortie(1); }; // traitement des énergies // on incrémente l'énergie élastique energ.ChangeEnergieElastique(energ_t.EnergieElastique() + 0.5 * ((sigHH + sigHH_nn) && delta_epsBB)); // récup de la compressibilité (-p_point = compress * I_D, S_point = cisaille * D_barre) module_compressibilite = Kc_use/3.; // en fait il s'agit ici de la compressibilité tangente // en faisant l'analogie : K == Kc/3. = E/3/(1-2nu) et f==E // on obtient: nu = 1/2 (1-f/Kc) et G = Kc*f/(3*Kc - f) module_cisaillement = Kc_use * f_use/(3. * Kc_use - f_use); // calcul des déformations transversales // \varepsilon_i^{.i}~(t+\Delta t) = \varepsilon_i^{.i} (t) + //\frac{1}{2} \left ( \frac{{\Delta_t^{t+\Delta t} \sigma_1^{.1}}}{3~K_c - \Delta_t^{t+\Delta t} \varepsilon_1^{.1} } \right ) save_resul.eps22 = save_resul.eps22_t + 0.5 *(delta_sigmaBH(1,1)/(Kc_use) - delta_epsBH(1,1) ); save_resul.eps33 = save_resul.eps22; // ******** arrête ici ******** // il faut également traiter les déformations eps22 et 33 qui sont stockées et qui dépendent du type de // transport (contrairement à ce qui a été défini dans Hypo_hooke1D::HsurH0 qu'il faut revoir // du coup il faut également revoir la def du type de def, ou plutôt sa cohérence avec le type de transport // ici, c'est le type de transport qui va définir la def utilisée et non l'utilisateur // donc au moment de la lecture il faut mettre un message d'erreur si l'utilisateur choisit une def différente du type de transport -> donc modifier : Hypo_hooke1D::LectureDonneesParticulieres // ... à faire LibereTenseur(); }; // calcul des contraintes a t+dt et de ses variations void Hypo_hooke1D::Calcul_DsigmaHH_tdt (TenseurHH& sigHH_t,TenseurBB& ,DdlElement & tab_ddl ,BaseB& giB_t,TenseurBB & gijBB_t_,TenseurHH & gijHH_t_,BaseB& giB_tdt,Tableau & d_giB_tdt ,BaseH& giH_tdt,Tableau & d_giH_tdt, TenseurBB & epsBB_tdt,Tableau & d_epsBB,TenseurBB & delta_epsBB_,TenseurBB & gijBB_tdt, TenseurHH & gijHH_tdt,Tableau & d_gijBB_tdt, Tableau & d_gijHH_tdt,double& jacobien_0,double& jacobien, Vecteur& d_jacobien_tdt,TenseurHH& sigHH_,Tableau & d_sigHH ,EnergieMeca & energ,const EnergieMeca & energ_t,double& module_compressibilite,double& module_cisaillement ,const Met_abstraite::Impli& ex) { #ifdef MISE_AU_POINT if (delta_epsBB_.Dimension() != 1) { cout << "\nErreur : la dimension devrait etre 1 !\n"; cout << " Hypo_hooke1D::Calcul_DsigmaHH_tdt\n"; Sortie(1); }; if (tab_ddl.NbDdl() != d_gijBB_tdt.Taille()) { cout << "\nErreur : le nb de ddl est != de la taille de d_gijBB_tdt !\n"; cout << " Hypo_hooke1D::Calcul_DsigmaHH_tdt\n"; Sortie(1); }; #endif const Tenseur1BB & epsBB = *((Tenseur1BB*) &epsBB_tdt); // passage en dim 1 const Tenseur1BB & delta_epsBB = *((Tenseur1BB*) &delta_epsBB_); // passage en dim 1 const Tenseur1HH & gijHH = *((Tenseur1HH*) &gijHH_tdt); // " " " " const Tenseur1BB & gijBB = *((Tenseur1BB*) &gijBB_tdt); // " " " " Tenseur1HH & sigHH = *((Tenseur1HH*) &sigHH_); // " " " " Tenseur1HH & sigHH_nn = *((Tenseur1HH*) &sigHH_t); // " " " " const Tenseur1HH & gijHH_t = *((Tenseur3HH*) &gijHH_t_); // " " " " const Tenseur1BB & gijBB_t = *((Tenseur3BB*) &gijBB_t_); // " " " " SaveResulLoi_Hypo1D & save_resul = *((SaveResulLoi_Hypo1D*) saveResul); // tenseur intermédiaires utilisées selon les cas (par forcément dans tous les cas !!) Tenseur1BH sigBH_n;Tenseur1HH sigHH_n; Tenseur1BB sigBB_n;Tenseur1BB sig_interBB_n; switch (type_derive) {case -1: // cas d'une dérivée de jauman: 1/2 deux fois covariant + deux fois contra {sig_interBB_n = gijBB_t * sigHH_nn * gijBB_t; sigBH_n = 0.5*( sig_interBB_n * gijHH + gijBB * sigHH_nn) ; sigHH_n = gijHH * sigBH_n ; sigBB_n = sigBH_n * gijBB; // pour la déformation cumulée associée, Tenseur1HH eps_interHH_n = gijHH_t * save_resul.eps_cumulBB_t * gijHH_t; Tenseur1BB epsBB_n = 0.5*( gijBB * eps_interHH_n * gijBB + save_resul.eps_cumulBB_t); save_resul.eps_cumulBB = delta_epsBB + epsBB_n; break; } case 0: // cas d'une dérivée de Lie deux fois covariantes {sigBB_n = gijBB_t * sigHH_nn * gijBB_t; sigBH_n = sigBB_n * gijHH ; sigHH_n = gijHH * sigBH_n ; // pour la déformation cumulée associée, on sait que cela donne directement la def d'Almansi // mais la def passée en paramètre pourrait ne pas être d'Almansi, donc on calcule // quand même la def cumulée save_resul.eps_cumulBB = save_resul.eps_cumulBB_t + delta_epsBB; break; } case 1: // cas d'une dérivée de Lie deux fois contravariantes {sigHH_n = sigHH_nn; sigBH_n = gijBB * sigHH_n; sigBB_n = sigBH_n * gijBB; // pour la déformation cumulée associée, Tenseur1HH eps_interHH_n = gijHH_t * save_resul.eps_cumulBB_t * gijHH_t; Tenseur1BB epsBB_n = gijBB * eps_interHH_n * gijBB ; save_resul.eps_cumulBB = delta_epsBB + epsBB_n; break; } }; // opération de transmission de la métrique const Met_abstraite::Impli* ex_impli = &ex; const Met_abstraite::Expli_t_tdt* ex_expli_tdt = NULL; const Met_abstraite::Umat_cont* ex_expli = NULL; // cas de la thermo dépendance, on calcul les grandeurs if (f_temperature != NULL) { f = f_temperature->Valeur(*temperature);} else if (f_nD != NULL) // là il faut calculer la fonction nD { // on utilise la méthode générique de loi abstraite list list_save; // inter pour l'appel de la fonction list_save.push_back(saveResul); Tableau & tab_val = Loi_comp_abstraite::Loi_comp_Valeur_FnD_Evoluee (f_nD,1 // une seule valeur attendue en retour ,ex_impli,ex_expli_tdt,ex_expli ,NULL ,NULL ,&list_save ); /* // ici on utilise les variables connues aux pti, ou calculées à partir de // on commence par récupérer les conteneurs des grandeurs à fournir List_io & li_enu_scal = f_nD->Li_enu_etendu_scalaire(); List_io & li_quelc = f_nD->Li_equi_Quel_evolue(); bool absolue = true; // on se place systématiquement en absolu // on va utiliser la méhode Valeur_multi_interpoler_ou_calculer // pour les grandeurs strictement scalaire Tableau val_ddl_enum(Valeur_multi_interpoler_ou_calculer (absolue,TEMPS_tdt,li_enu_scal,ex_impli,ex_expli_tdt,ex_expli,NULL) ); // on utilise la méthode Valeurs_Tensorielles_interpoler_ou_calculer // pour les Coordonnees et Tenseur Valeurs_Tensorielles_interpoler_ou_calculer (absolue,TEMPS_tdt,li_quelc,ex_impli,ex_expli_tdt,ex_expli,NULL); // calcul de la valeur et retour dans tab_ret Tableau & tab_val = f_nD->Valeur_FnD_Evoluee(&val_ddl_enum,&li_enu_scal,&li_quelc,NULL,NULL); #ifdef MISE_AU_POINT if (tab_val.Taille() != 1) { cout << "\nErreur : la fonction nD relative a la fonction f " << " doit calculer un scalaire or le tableau de retour est de taille " << tab_val.Taille() << " ce n'est pas normal !\n"; cout << " Hypo_hooke1D::Calcul_DsigmaHH_tdt \n"; Sortie(1); }; #endif */ // on récupère le premier élément du tableau uniquement f = tab_val(1); }; if (!compress_thermophysique) // sinon kc est mis à jour avec la méthode CalculGrandeurTravail(..) {if (Kc_temperature != NULL) {Kc = Kc_temperature->Valeur(*temperature);} else if (Kc_nD != NULL) // là il faut calculer la fonction nD { // on utilise la méthode générique de loi abstraite list list_save; // inter pour l'appel de la fonction list_save.push_back(saveResul); Tableau & tab_val = Loi_comp_abstraite::Loi_comp_Valeur_FnD_Evoluee (Kc_nD,1 // une seule valeur attendue en retour ,ex_impli,ex_expli_tdt,ex_expli ,NULL ,NULL ,&list_save ); /* // ici on utilise les variables connues aux pti, ou calculées à partir de // on commence par récupérer les conteneurs des grandeurs à fournir List_io & li_enu_scal = Kc_nD->Li_enu_etendu_scalaire(); List_io & li_quelc = Kc_nD->Li_equi_Quel_evolue(); bool absolue = true; // on se place systématiquement en absolu // on va utiliser la méhode Valeur_multi_interpoler_ou_calculer // pour les grandeurs strictement scalaire Tableau val_ddl_enum(Valeur_multi_interpoler_ou_calculer (absolue,TEMPS_tdt,li_enu_scal,ex_impli,ex_expli_tdt,ex_expli,NULL) ); // on utilise la méthode Valeurs_Tensorielles_interpoler_ou_calculer // pour les Coordonnees et Tenseur Valeurs_Tensorielles_interpoler_ou_calculer (absolue,TEMPS_tdt,li_quelc,ex_impli,ex_expli_tdt,ex_expli,NULL); // calcul de la valeur et retour dans tab_ret Tableau & tab_val = Kc_nD->Valeur_FnD_Evoluee(&val_ddl_enum,&li_enu_scal,&li_quelc,NULL,NULL); #ifdef MISE_AU_POINT if (tab_val.Taille() != 1) { cout << "\nErreur : la fonction nD relative a la fonction Kc " << " doit calculer un scalaire or le tableau de retour est de taille " << tab_val.Taille() << " ce n'est pas normal !\n"; cout << " Hypo_hooke1D::Calcul_SigmaHH\n"; Sortie(1); }; #endif */ // on récupère le premier élément du tableau uniquement Kc = tab_val(1); }; }; // cas d'une non-linéarité en fonction de la déformation Tenseur1BH epsBH; Courbe1D::ValDer valder_Kc,valder_f; double Kc_base=Kc;double Kc_use=Kc;// pour ne pas changer les valeurs à chaque passage !! double f_base=f;double f_use=f;// idem if ((Kc_IIeps != NULL)||(f_IIeps != NULL)) { epsBH = epsBB * gijHH;double II_eps=epsBH.II(); if (Kc_IIeps != NULL) { valder_Kc = Kc_IIeps->Valeur_Et_derivee(II_eps); Kc_use *= valder_Kc.valeur; }; if (f_IIeps != NULL) { valder_f = f_IIeps->Valeur_Et_derivee(II_eps); f_use *= valder_f.valeur; }; }; // sauvegarde des paramètres matériau save_resul.Kc = Kc_use; save_resul.f = f_use; // calcul de la contrainte unidirectionnelle Tenseur1BH delta_epsBH = delta_epsBB * gijHH; Tenseur1BH delta_sigmaBH(delta_epsBH*f_use); Tenseur1BH sigBH = sigBH_n + delta_sigmaBH; sigHH = gijHH * sigBH; // on applique éventuellement une restriction sur la contrainte // =0 -> pas de restriction // = -1 : traction uniquement autorisée, la compression est mise à 0 // = 1 : compression uniquement autorisée, la traction est mise à 0 int restriction_effectuee=0; switch (restriction_traction_compression) { case 0: break; // cas sans restriction case -1: // traction uniquement autorisée, la compression est mise à 0 {if (sigHH(1,1)<0.) sigHH.Coor(1,1) = 0.; restriction_effectuee = -1; break; } case 1: // compression uniquement autorisée, la traction est mise à 0 {if (sigHH(1,1)>0.) sigHH.Coor(1,1) = 0.; restriction_effectuee = 1; break; } default: cout << "\n erreur Hypo_hooke1D::Calcul_DsigmaHH_tdt " << " restriction_traction_compression= "<< restriction_traction_compression << " ce qui n'est pas exploitable !! "; Sortie(1); }; // traitement des énergies // on incrémente l'énergie élastique energ.ChangeEnergieElastique(energ_t.EnergieElastique() + 0.5 * ((sigHH + sigHH_nn) && delta_epsBB)); // récup de la compressibilité (-p_point = compress * I_D, S_point = cisaille * D_barre) module_compressibilite = Kc_use/3.; // en fait il s'agit ici de la compressibilité tangente // en faisant l'analogie : K == Kc/3. = E/3/(1-2nu) et f==E // on obtient: nu = 1/2 (1-f/Kc) et G = Kc*f/(3*Kc - f) module_cisaillement = Kc_use * f_use/(3. * Kc_use - f_use); // calcul des déformations transversales // \varepsilon_i^{.i}~(t+\Delta t) = \varepsilon_i^{.i} (t) + //\frac{1}{2} \left ( \frac{{\Delta_t^{t+\Delta t} \sigma_1^{.1}}}{3~K_c - \Delta_t^{t+\Delta t} \varepsilon_1^{.1} } \right ) save_resul.eps22 = save_resul.eps22_t + 0.5 *(delta_sigmaBH(1,1)/(Kc_use) - delta_epsBH(1,1) ); save_resul.eps33 = save_resul.eps22; ////--- debug //cout << "\n Hypo_hooke1D::Calcul_DsigmaHH_tdt"; //cout << "\n save_resul.eps22= "<Valeur(*temperature); if (!compress_thermophysique) // sinon kc est mis à jour avec la méthode CalculGrandeurTravail(..) {if (Kc_temperature != NULL) Kc = Kc_temperature->Valeur(*temperature);}; // cas d'une non-linéarité en fonction de la déformation Tenseur3BH DepsBH; Tenseur3BB Deps_barre_BB;Tenseur3HH Deps_barre_HH; // init Tenseur3BH Deps_barre_BH; double IDeps=0.; // " Tenseur3BH epsBH; Courbe1D::ValDer valder_Kc,valder_f; double Kc_base=Kc;double Kc_use=Kc;// pour ne pas changer les valeurs à chaque passage !! double f_base=f;double f_use=f;// idem static const double untier=1./3.; if (en_base_orthonormee) {DepsBH = DepsBB.MonteDernierIndice(); IDeps=DepsBH.Trace(); Deps_barre_BH = DepsBH - (untier * IDeps) * IdBH3; Deps_barre_BB = Deps_barre_BH * IdBB3; Deps_barre_HH = IdHH3 * Deps_barre_BH ; epsBH = epsBB.MonteDernierIndice(); // deformation en mixte } else {DepsBH = DepsBB * gijHH; IDeps=DepsBH.Trace(); Deps_barre_BH = DepsBH - (untier * IDeps) * IdBH3; Deps_barre_BB = Deps_barre_BH * gijBB; Deps_barre_HH = gijHH * Deps_barre_BH ; epsBH = epsBB * gijHH; // deformation en mixte }; if ((Kc_IIeps != NULL)||(f_IIeps != NULL)) { double Ieps = epsBH.Trace(); Tenseur3BH eps_barre_BH = epsBH - (untier * Ieps) * IdBH3; double II_eps=epsBH.II(); if (Kc_IIeps != NULL) { valder_Kc = Kc_IIeps->Valeur_Et_derivee(II_eps); Kc_use *= valder_Kc.valeur; }; if (f_IIeps != NULL) { valder_f = f_IIeps->Valeur_Et_derivee(II_eps); f_use *= valder_f.valeur; }; }; // recup de l'incrément de temps double deltat=ParaGlob::Variables_de_temps().IncreTempsCourant(); // le calcul de la contrainte correspond à l'intégration d'une équation différencielle // delta_eps_barre=1/f delta_sigma/dt pour la partie déviatoire // et pour la partie sphérique // I_D = 1/(Kc) I_Sig_point double Isig_n = sigBH_n.Trace(); Tenseur3BH SBH_n = sigBH_n - (untier * Isig_n) * IdBH3; double unSurDeltat=0; if (Abs(deltat) >= ConstMath::trespetit) {unSurDeltat = 1./deltat;} else // si l'incrément de temps est tres petit on remplace 1/deltat par un nombre tres grand { // non un pas de temps doit être positif !! or certaine fois il peut y avoir des pb if (unSurDeltat < 0) { cout << "\n le pas de temps est négatif !! "; }; unSurDeltat = ConstMath::tresgrand; }; // cas de la partie sphérique scalaire double Isigma = Isig_n + deltat * Kc_use * IDeps; // cas de la partie déviatorique // if (en_base_orthonormee) // { switch (cas_calcul) // { case 0: // calcul normal (tous les termes) // { // la partie sphérique est déjà calculé, cas de la partie déviatorique // Tenseur3HH SHH = IdHH3 * (SBH_n + (deltat * f_use) * Deps_barre_BH); // sigHH = SHH + (untier * Isigma) * IdHH3; // break; // } // case 1: // calcul de la partie déviatorique seule // { sigHH = IdHH3 * (SBH_n + (deltat * f_use) * Deps_barre_BH); // break; // } // case 2: // calcul de la partie sphérique seule // { sigHH = (untier * Isigma) * IdHH3; // break; // } // default: // { cout << "\n erreur l'indicateur cas_calcul= " << cas_calcul << " n'a pas une valeur correcte !! " // << "\n Hypo_hooke1D::Calcul_dsigma_deps (.... "; // Sortie(1); // } // }; // } // else // { switch (cas_calcul) // { case 0: // calcul normal (tous les termes) // { // la partie sphérique est déjà calculé, cas de la partie déviatorique // Tenseur3HH SHH = gijHH * (SBH_n + deltat * f_use * Deps_barre_BH); // sigHH = SHH + (untier * Isigma) * gijHH ; // break; // } // case 1: // calcul de la partie déviatorique seule // { sigHH = gijHH * (SBH_n + deltat * f_use * Deps_barre_BH); // break; // } // case 2: // calcul de la partie sphérique seule // { sigHH = (untier * Isigma) * gijHH ; // *IdBH3; // break; // } // default: // { cout << "\n erreur l'indicateur cas_calcul= " << cas_calcul << " n'a pas une valeur correcte !! " // << "\n Hypo_hooke1D::Calcul_dsigma_deps (.... "; // Sortie(1); // } // }; // }; // // // ////--- debug //cout << "\n Hypo_hooke1D::Calcul_dsigma_deps"; //cout << "\n IDeps= "<Monte4Indices(); // switch (cas_calcul) // { case 0: // calcul normal (tous les termes) // { d_sigma_depsHHHH = f_use * d_deltat_Deps_barre_HHHH; // if (f_IIeps != NULL) // d_sigma_depsHHHH += ((Tenseur3BBBB*) &((f_base * valder_f.derivee * 2. * deltat) // * Tenseur3BBBB::Prod_tensoriel((Deps_barre_BH * IdBB3),epsBB)))->Monte4Indices(); // d_sigma_depsHHHH += (untier * Kc) * IdHHHH3; // if (Kc_IIeps != NULL) // d_sigma_depsHHHH += ((Tenseur3BBBB*) &((untier * Kc_base * valder_Kc.derivee * 2. ) // * Tenseur3BBBB::Prod_tensoriel(IdBB3, epsBB)))->Monte4Indices(); // break; // } // case 1: // calcul de la partie déviatorique seule // { d_sigma_depsHHHH = f_use * d_deltat_Deps_barre_HHHH; // if (f_IIeps != NULL) // d_sigma_depsHHHH += ((Tenseur3BBBB*) &((f_base * valder_f.derivee * 2. * deltat) // * Tenseur3BBBB::Prod_tensoriel((Deps_barre_BH * IdBB3),epsBB)))->Monte4Indices(); // break; // } // case 2: // calcul de la partie sphérique seule // { d_sigma_depsHHHH = (untier * Kc) * IdHHHH3; // if (Kc_IIeps != NULL) // d_sigma_depsHHHH += ((Tenseur3BBBB*) &((untier * Kc_base * valder_Kc.derivee * 2. ) // * Tenseur3BBBB::Prod_tensoriel(IdBB3, epsBB)))->Monte4Indices(); // break; // } // default: // { cout << "\n erreur l'indicateur cas_calcul= " << cas_calcul << " n'a pas une valeur correcte !! " // << "\n Hypo_hooke1D::Calcul_DsigmaHH_tdt (.... "; // Sortie(1); // } // }; } else // sinon cas où les bases sont curvilignes { // calcul de variables intermédiaires I_x_I_HHHH=Tenseur3HHHH::Prod_tensoriel(gijHH,gijHH); I_xbarre_I_HHHH=Tenseur3HHHH::Prod_tensoriel_barre(gijHH,gijHH); Tenseur3HH epsHH(gijHH * epsBH); I_x_eps_HHHH=Tenseur3HHHH::Prod_tensoriel(gijHH,epsHH); Tenseur3HH Deps_HH(gijHH * DepsBH); I_x_D_HHHH=Tenseur3HHHH::Prod_tensoriel(gijHH,Deps_HH); I_xbarre_D_HHHH=Tenseur3HHHH::Prod_tensoriel_barre(gijHH,Deps_HH); // variation de sigma_n et de la trace en fonction du transport switch (type_derive) { case -1: // // cas d'une dérivée de jauman: 1/2 deux fois covariant + deux fois contra {// pour info sigBH_n = 0.5*(gijBB * sigHH_n + (gijBB_t * sigHH_n * gijBB_t) * gijHH) d_sig_t_HHHH = Tenseur3HHHH::Prod_tensoriel_barre(gijHH,sigHH_n); d_spherique_sig_t_HHHH = (1. /6.) * (Tenseur3HHHH::Prod_tensoriel(gijHH,sigHH_n)+ Isig_n * I_xbarre_I_HHHH); break; } case 0: // cas d'une dérivée de Lie deux fois covariantes {// pour info sigBH_n = (gijBB_t * sigHH_n * gijBB_t) * gijHH d_sig_t_HHHH = 2.*Tenseur3HHHH::Prod_tensoriel_barre(gijHH,sigHH_n); d_spherique_sig_t_HHHH = untier * (Tenseur3HHHH::Prod_tensoriel(gijHH,sigHH_n)+ Isig_n * I_xbarre_I_HHHH); break; } case 1: // cas d'une dérivée de Lie deux fois contravariantes {// pour info sigBH_n = gijBB * sigHH_n, ici il n'y a aucune conséquence // sur sigma(t), par contre il y en a une sur la trace d_sig_t_HHHH.Inita(0.); d_spherique_sig_t_HHHH = untier * (Tenseur3HHHH::Prod_tensoriel(gijHH,sigHH_n)+ Isig_n * I_xbarre_I_HHHH); break; } }; // on choisit entre les différents cas double e01,e02;e01=e02=0.; // les coefficients pour la partie transportée double b01,b02,b03,b04,b05;b01=b02=b03=b04=b05=0.; // les coefficients pour le reste de l'équation // l'expression générique est: // d_sigma_depsHHHH = b01 * I_x_I_HHHH + b02 * I_x_D_HHHH + b03 * I_xbarre_I_HHHH // + b04 * I_xbarre_D_HHHH // + e01 * d_sig_t_HHHH + e02 * d_spherique_sig_t_HHHH; // en fonction des coefficients nuls on simplifie // switch (cas_calcul) // { case 0: // calcul normal (tous les termes) // cas complet // { b01= (Kc_use - f_use) * untier; // b02= -2.* (Kc_use - f_use) * untier * deltat; // b03= -2. * (Kc_use - f_use) * untier * deltat * IDeps + f_use; // b04= -4. * f_use * deltat; // e01= f_use; // e02= untier; // d_sigma_depsHHHH = b01 * I_x_I_HHHH + b02 * I_x_D_HHHH // + b03 * I_xbarre_I_HHHH + b04 * I_xbarre_D_HHHH // + e01 * d_sig_t_HHHH + e02 * d_spherique_sig_t_HHHH; // break; // } // case 1: // calcul de la partie déviatorique seule // { b01= - f_use * untier; // b02= 2.* f_use * untier * deltat; // b03= 2. * f_use * untier * deltat * IDeps + f_use; // b04= -4. * f_use * deltat; // e01= f_use; // e02= untier; // d_sigma_depsHHHH = b01 * I_x_I_HHHH + b02 * I_x_D_HHHH // + b03 * I_xbarre_I_HHHH + b04 * I_xbarre_D_HHHH // + e01 * d_sig_t_HHHH + e02 * d_spherique_sig_t_HHHH; // break; // } // case 2: // calcul de la partie sphérique seule // { b01= Kc_use * untier; // b02= -2.*Kc_use * untier * deltat; // b03= -2.*Kc_use * untier * deltat * IDeps; // b04= 0.; // e01= 0.; // e02= untier; // d_sigma_depsHHHH = b01 * I_x_I_HHHH + b02 * I_x_D_HHHH // + b03 * I_xbarre_I_HHHH // + e02 * d_spherique_sig_t_HHHH; // break; // } // default: // { cout << "\n erreur l'indicateur cas_calcul= " << cas_calcul << " n'a pas une valeur correcte !! " // << "\n Hypo_hooke1D::Calcul_dsigma_deps (.... "; // Sortie(1); // } // }; }; // traitement des énergies // on incrémente l'énergie élastique energ.ChangeEnergieElastique(energ_t.EnergieElastique() + 0.5 * deltat * ((sigHH + sigHH_nn) && DepsBB)); // récup de la compressibilité (-p_point = compress * I_D, S_point = cisaille * D_barre) module_compressibilite = Kc_use/3.; module_cisaillement = 0.5 * f_use; LibereTenseur(); LibereTenseurQ(); };