// 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: . /************************************************************************ * DATE: 23/01/97 * * $ * * AUTEUR: G RIO (mailto:gerardrio56@free.fr) * * $ * * PROJET: Herezh++ * * $ * ************************************************************************ * BUT: Classe de base de Defintion des differents algoritmes * * de resolution. * * $ * * '''''''''''''''''''''''''''''''''''''''''''''''''''''''''''' * * * VERIFICATION: * * * * ! date ! auteur ! but ! * * ------------------------------------------------------------ * * ! ! ! ! * * $ * * '''''''''''''''''''''''''''''''''''''''''''''''''''''''''''' * * MODIFICATIONS: * * ! date ! auteur ! but ! * * ------------------------------------------------------------ * * $ * ************************************************************************/ #ifndef ALGO_H #define ALGO_H #include #include "UtilLecture.h" #include "EnumTypeCalcul.h" #include "MotCle.h" // #include "bool.h" #include "string" #include "Tableau_T.h" #include "LesMaillages.h" #include "LesReferences.h" #include "LesCourbes1D.h" #include "LesFonctions_nD.h" #include "DiversStockage.h" //#include "Charge.h" #include "LesCondLim.h" #include "ParaGlob.h" #include "ParaAlgoControle.h" #include "LesContacts.h" #include "List_io.h" #include "Visualisation.h" #include "Visualisation_maple.h" #include "Visualisation_geomview.h" #include "Visualisation_Gid.h" #include "Visualisation_Gmsh.h" #include "Nb_assemb.h" #include "Temps_CPU_HZpp.h" #include "TypeQuelconqueParticulier.h" #include "VariablesExporter.h" #ifdef UTILISATION_MPI #include "Distribution_CPU.h" #include "ResRaid_MPI.h" #endif // peut-être à modifier sur linux #ifdef BIBLIOTHEQUE_STL_CODE_WARRIOR #include // cas code warrior #else #ifdef SYSTEM_MAC_OS_X_unix #include // cas unix sur osX #else #include // cas linux #endif #endif class Resultats; class Charge; class AlgoriCombine; /** @defgroup Les_algorithmes_de_resolutions_globales * * BUT: groupe des algorithmes de résolutions globales * * * \author Gérard Rio * \version 1.0 * \date 10/02/2004 * \brief groupe des algorithmes de résolutions globales * */ /// @addtogroup Les_algorithmes_de_resolutions_globales /// @{ /// /// BUT: Classe de base de Defintion des differents algoritmes /// de resolution. class Algori { friend class Charge;friend class AlgoriCombine; public : // deux container de base qui servent pour les données externes class DeuxString { // surcharge de l'operator de lecture /écriture friend istream & operator >> (istream &, Algori::DeuxString &); friend ostream & operator << (ostream &, const Algori::DeuxString &); public : string nomFichier; string nomMaillage; }; class ListDeuxString { // surcharge de l'operator de lecture /écriture friend istream & operator >> (istream & , Algori::ListDeuxString & ); friend ostream & operator << (ostream &, const Algori::ListDeuxString &); public : List_io infos; }; // surcharge de l'operator d'ecriture pour avoir des tableaux mais // non opérationnelle friend ostream & operator << (ostream & sort, const Algori & ) { // tout d'abord un indicateur cout << "\n **** ecriture directe d'un algorithme non implante, il faut utiliser Affiche()"; Sortie(1); return sort; }; // CONSTRUCTEURS : Algori (); // par defaut // constructeur en fonction du type de calcul // du sous type (pour les erreurs, remaillage etc...) // il y a ici lecture des parametres attaches au type Algori (EnumTypeCalcul type,const bool avec_typeDeCalcul ,const list & soustype ,const list & avec_soustypeDeCalcul ,UtilLecture& entreePrinc ); // constructeur de copie Algori (const Algori& algo); // constructeur de copie à partie d'une instance indifférenciée virtual Algori * New_idem(const Algori* algo) const =0 ; // DESTRUCTEUR : virtual ~Algori () ; // METHODES PUBLIQUES : // ---------- static pour la création d'un algorithme particulier --------- // ramène un pointeur sur l'algorithme spécifique correspondant aux paramètres // IMPORTANT : il y a création de l'algorithme correspondant (utilisation d'un new) static Algori* New_Agori(EnumTypeCalcul type,const bool avec_typeDeCalcul ,const list & soustype ,const list & avec_soustypeDeCalcul ,UtilLecture& entreePrinc ); // ramène un tableau de pointeurs sur tous les algorithmes spécifique existants // IMPORTANT : il y a création des algorithmes (utilisation d'un new) static Tableau New_tous_les_Algo (const bool avec_typeDeCalcul ,const list & soustype ,const list & avec_soustypeDeCalcul ,UtilLecture& entreePrinc ); //lecture des parametres de controle // peut être modifié dans le cas d'un algo particulier virtual void Lecture(UtilLecture & entreePrinc,ParaGlob & paraGlob,LesMaillages& lesMail); // cohérence du type de convergence adopté en fonction de l'algorithme: par exemple on vérifie que l'on n'utilise // pas une convergence sur l'énergie cinétique avec un algo statique !! void Coherence_Algo_typeConvergence(); // définition du nombre de cas d'assemblage, ce qui dépend du cas de calcul traité int NbCasAssemblage() const {return nb_CasAssemblage;}; // execution de l'algorithme en fonction du type de calcul virtual void Execution(ParaGlob *,LesMaillages *,LesReferences*, LesCourbes1D*,LesFonctions_nD*,VariablesExporter*, LesLoisDeComp*, DiversStockage*, Charge*,LesCondLim*,LesContacts*,Resultats* ) = 0; //------- décomposition en 3 du calcul d'équilibre ------------- // a priori : InitAlgorithme et FinCalcul ne s'appellent qu'une fois, // par contre : CalEquilibre peut s'appeler plusieurs fois, le résultat sera différent si entre deux calcul // certaines variables ont-été changés // initialisation virtual void InitAlgorithme(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* ,LesFonctions_nD* ,VariablesExporter*,LesLoisDeComp* ,DiversStockage*,Charge*,LesCondLim*,LesContacts*,Resultats* ) = 0; // mise à jour virtual void MiseAJourAlgo(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* ,LesFonctions_nD* ,VariablesExporter* ,LesLoisDeComp* ,DiversStockage*,Charge*,LesCondLim*,LesContacts*,Resultats* ) = 0; // calcul de l'équilibre // si tb_combiner est non null -> un tableau de 2 fonctions // - la première fct dit si on doit valider ou non le calcul à convergence ok, // - la seconde dit si on doit sortir de la boucle ou non à convergence ok // // si la validation est effectuée, la sauvegarde pour le post-traitement est également effectuée // en fonction de la demande de sauvegard, // sinon pas de sauvegarde pour le post-traitement à moins que l'on a demandé un mode debug // qui lui fonctionne indépendamment virtual void CalEquilibre(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* ,LesFonctions_nD* ,VariablesExporter* ,LesLoisDeComp* ,DiversStockage*,Charge*,LesCondLim*,LesContacts*,Resultats* ,Tableau < Fonction_nD* > * tb_combiner) = 0; // dernière passe virtual void FinCalcul(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* ,LesFonctions_nD* ,VariablesExporter* ,LesLoisDeComp* ,DiversStockage*,Charge*,LesCondLim*,LesContacts*,Resultats* ) = 0; // sortie du schemaXML: en fonction de enu virtual void SchemaXML_Algori(ofstream& sort,const Enum_IO_XML enu) const = 0; // affichage des parametres void Affiche() const ; // affichage de tous les parametres void Affiche1() const ; // affichage du type de calcul void Affiche2() const ; // affichage des parametres de controle // définition interactive exhaustive des paramètres de contrôle indépendamment de l'algorithme // écriture du fichier de commande void Info_commande_ParaAlgoControle(UtilLecture& lec) {pa.Info_commande_ParaAlgoControle(lec);} // sortie d'info sur un increment de calcul pour les ddl // nb_casAssemb : donne le numéro du cas d'assemblage void InfoIncrementDdl(LesMaillages * lesMail, int inSol,double maxDeltaDdl,const Nb_assemb& nb_casAssemb); // sortie d'info sur un increment de calcul pour les réaction // dans le cas d'un compteur, par exemple dans le cas implicite // nb_casAssemb : donne le numéro du cas d'assemblage void InfoIncrementReac(LesMaillages * lesMail,int compteur, int inreaction,double maxreaction,const Nb_assemb& nb_casAssemb); // idem mais dans le cas où il n'y a pas de compteur, par exemple // dans le cas explicite // nb_casAssemb : donne le numéro du cas d'assemblage void InfoIncrementReac(LesMaillages * lesMail, int inreaction,double maxreaction,const Nb_assemb& nb_casAssemb); // retourne le type de calcul EnumTypeCalcul TypeDeCalcul() const { return typeCalcul;}; // dans le cas où un grand nombre d'information nécessaire pour l'algorithme // est stocké de manière externe sur un fichier // exemple : calcul d'erreur a partir // d'un champ d'information stocké de manière externe sur un fichier // l'acquisition du (ou des) lieu(s) (flot) où sont stocké ces infos est effectuée // ce lieu fait parti des paramètres de l'algorithme void DefFlotExterne(UtilLecture* entreePrinc); // indique en retour s'il y a des données externes ou non bool ExisteDonneesExternes() const { if (noms_fichier.Taille() != 0) return true; else return false;}; // retourne le nombre total de flot externe où il faut lire // des données externes int NbTypeFlotExt() const {return noms_fichier.Taille();}; // retourne la liste des noms du fichier et des noms de maillage correspondant // correspondant aux données externes nbf const ListDeuxString& Noms_fichier(int nbf) const // {return ((char*)(nom_fichier(nbf)).c_str());}; {return noms_fichier(nbf);}; // retourne le type de données externe pour le flot nbf const int TypeDonneesExternes(int nbf) const { return typeFlotExterne(nbf);}; // retourne le numéro demandé de restart const int Num_restart() const {return pa.Restart();}; // retourne le numéro d'incrément const int Num_increment() const {return icharge;} ; // modif d'icharge, utilisé par les classes externes amies éventuellement void Affectation_icharge(int icharge_) {icharge = icharge_;}; // retourne vraie si la sauvegarde est activée const bool Active_sauvegarde() const {return pa.Sauvegarde();}; // sauvegarde générale sur base info // cas donne le niveau de sauvegarde // = 0 : initialisation de la sauvegarde -> c'est-à-dire de la sortie base info // = 1 : on sauvegarde tout // = 2 : on sauvegarde uniquement les données variables (supposées comme telles) // incre : numero d'incrément auquel on sauvegarde // éventuellement est définit de manière spécifique pour chaque algorithme // dans les classes filles // type_incre :signal permet de localiser le dernier incrément virtual void Ecriture_base_info (int cas,LesMaillages *lesMaillages, LesReferences* lesReferences,LesCourbes1D* lesCourbes1D, LesFonctions_nD* lesFonctionsnD, LesLoisDeComp* lesLoisDeComp, DiversStockage* diversStockage,Charge* charge, LesCondLim* lesCondlim,LesContacts* lesContacts,Resultats* resultats, OrdreVisu::EnumTypeIncre type_incre,int incre = 0); // cas d'un démarrage à partir de base_info // éventuellement est définit de manière spécifique pour chaque algorithme // dans les classes filles // si inc_voulu est négatif cela signifie que l'on est déjà positionné sur un // incrément voulu et qu'il faut simplement faire la lecture des infos // cas = 3 : on met à jour uniquement les données variables, il y a modification des grandeurs, // mais pas redimentionnement lorsque cela est viable (sinon on redimentionne) virtual void Lecture_base_info (int cas,LesMaillages *lesMaillages, LesReferences* lesReferences,LesCourbes1D* lesCourbes1D, LesFonctions_nD* lesFonctionsnD, LesLoisDeComp* lesLoisDeComp, DiversStockage* diversStockage,Charge* charge, LesCondLim* lesCondlim,LesContacts* lesContacts,Resultats* resultats, int inc_voulu=0); // visualisation intéractive via le standard vrml // la fonction est virtuelle ce qui lui permet d'être éventuellement facilement surchargée // dans les algorithmes dérivées virtual void Visu_vrml(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* , LesFonctions_nD* ,LesLoisDeComp* ,DiversStockage*, Charge*,LesCondLim*,LesContacts*,Resultats* ); // visualisation intéractive par fichier de points, utilisant le format maple ou gnuplot // la fonction est virtuelle ce qui lui permet d'être éventuellement facilement surchargée // dans les algorithmes dérivées virtual void Visu_maple(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* , LesFonctions_nD* ,LesLoisDeComp* ,DiversStockage*, Charge*,LesCondLim*,LesContacts*,Resultats* ); // visualisation intéractive via geomview // la fonction est virtuelle ce qui lui permet d'être éventuellement facilement surchargée // dans les algorithmes dérivées virtual void Visu_geomview(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* , LesFonctions_nD* ,LesLoisDeComp* ,DiversStockage*, Charge*,LesCondLim*,LesContacts*,Resultats* ); // visualisation intéractive via Gid // la fonction est virtuelle ce qui lui permet d'être éventuellement facilement surchargée // dans les algorithmes dérivées virtual void Visu_Gid(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* , LesFonctions_nD* ,LesLoisDeComp* ,DiversStockage*, Charge*,LesCondLim*,LesContacts*,Resultats* ); // visualisation intéractive via Gmsh // la fonction est virtuelle ce qui lui permet d'être éventuellement facilement surchargée // dans les algorithmes dérivées virtual void Visu_Gmsh(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* , LesFonctions_nD* ,LesLoisDeComp* ,DiversStockage*, Charge*,LesCondLim*,LesContacts*,Resultats* ); // lecture d'informations de visualisation dans un fichier de commande externe void LectureCommandeVisu(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* , LesFonctions_nD* ,LesLoisDeComp* ,DiversStockage*, Charge*,LesCondLim*,LesContacts*,Resultats* ); // visualisation au fil du calcul à partir d'un fichier de commande externe // type_incre signal ce qu'il y a a faire, il est modifié dans le programme !! donc le laisser en variable // en fait après le premier passage void VisuAuFilDuCalcul(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* ,LesFonctions_nD* ,LesLoisDeComp* ,DiversStockage* ,Charge*,LesCondLim*,LesContacts*,Resultats* ,OrdreVisu::EnumTypeIncre& type_incre ,int num_incre); // Ecriture des informations de visualisation dans un fichier de commande externe void EcritureCommandeVisu(); // sortie sur fichier des temps cpu virtual void Sortie_temps_cpu(const LesCondLim& lesCondLim , const Charge& charge,const LesContacts & contact); /* // visualisation automatique sans partie interactive void VisualisationAutomatique(ParaGlob * ,LesMaillages *,LesReferences*,LesCourbes1D* ,LesLoisDeComp* ,DiversStockage*, Charge*,LesCondLim*,LesContacts*,Resultats* ); */ // une méthode qui a pour objectif de terminer tous les comptages, utile // dans le cas d'un arrêt impromptu virtual void Arret_du_comptage_CPU(); //----------------------------------------- protégé ------------------------------ protected : // METHODES PROTEGEES : // préparation des conteneurs interne, utilisés après la lecture // où dans le cas où il n'y a pas de lecture directe // (pour un sous algo via un algo combiné par exemple) void Preparation_conteneurs_interne(LesMaillages& lesMail); // sauvegarde sur base info // cas donne le niveau de sauvegarde // = 0 : initialisation de la sauvegarde -> c'est-à-dire de la sortie base info // = 1 : on sauvegarde tout // = 2 : on sauvegarde uniquement les données variables (supposées comme telles) // incre : numero d'incrément auquel on sauvegarde // éventuellement est définit de manière spécifique pour chaque algorithme // dans les classes filles // si ptalgo est non nulle, on a une sortie des temps cpu, spécifique // à la classe AlgoriCombine passé en paramètre virtual void Ecriture_base_info (ofstream& sort,const int cas); // cas de la lecture spécifique à l'algorithme dans base_info virtual void Lecture_base_info(ifstream& ent,const int cas); // cas des paramètres spécifiques // écriture des paramètres dans la base info // = 1 : on écrit tout // = 2 : on écrot uniquement les données variables (supposées comme telles) virtual void Ecrit_Base_info_Parametre(UtilLecture& entreePrinc,const int& cas) = 0; // lecture des paramètres dans la base info // = 1 : on récupère tout // = 2 : on récupère uniquement les données variables (supposées comme telles) // choix = true : fonctionnememt normal // choix = false : la méthode ne doit pas lire mais initialiser les données à leurs valeurs par défaut // car la lecture est impossible virtual void Lecture_Base_info_Parametre(UtilLecture& entreePrinc,const int& cas,bool choix) = 0; // création d'un fichier de commande: cas des paramètres spécifiques virtual void Info_commande_parametres(UtilLecture& entreePrinc) = 0; // mise à jour par une classe friend (attention) des pointeurs de sous types void Change_affectation_pointeur_sous_type (const list & soustype ,const list & avec_soustypeDeCal) {soustypeDeCalcul = &soustype; avec_soustypeDeCalcul = &avec_soustypeDeCal;}; // initialisation de tous les paramètres de contrôle aux valeurs transmises void Init_ParaAlgoControle (ParaAlgoControle& paa) {pa = paa; if (type_calcul_visqu_critique) // cas particulier d'un amortissement critique pa.Change_amort_visco_artificielle(4); // on indique qu'il s'agit d'une gestion par un amortissement visqueux critique }; // récupération des paramètres de contrôles de l'algo ParaAlgoControle& ParaAlgoControle_de_lalgo() {return pa;}; // récupération de la position lue par l'algo dans le .BI, changé lors d'un restart par exemple // il s'agit de la sauvegarde par l'algo de la position qu'il a lue en restart streampos Debut_increment() const {return debut_increment;}; // ==============VARIABLES PROTEGEES : #ifdef UTILISATION_MPI // cas d'un calcul parallèle Distribution_CPU distribution_CPU_algo; // gestion de la distribution de charge sur CPU #endif // 1)----------- protégées: spécifiques à un algo ---------- int mode_debug; // permet de spécifier des options particulières pour débugger // controle de la sortie des informations: utilisé par les classes dérivées int permet_affichage; // pour permettre un affichage spécifique dans les méthodes, pour les erreurs et warning EnumTypeCalcul typeCalcul; bool avec_typeDeCalcul; // indique si le calcul doit être effectué ou pas int nb_CasAssemblage; // indique le nombre de cas d'assemblage Tableau paraTypeCalcul; // paramètres particuliers à chaque algorithme // identificateur secondaire optionnel //du type de Calcul renseigne sur le fait de calculer l'erreur, // une relocation, un raffinement, éventuel. list const * soustypeDeCalcul; list const *avec_soustypeDeCalcul; // indique les sous types actifs ou non ParaAlgoControle pa; // parametres de controle de l'algorithme // 2) ----------- protégées: commums à tous les algos -------- double temps_derniere_sauvegarde; // variable stockant le dernier temps de la sauvegarde sur .BI // liste de tous les différents incréments et temps des sauvegardes effectuées dans le .BI List_io list_incre_temps_sauvegarder; List_io list_incre_temps_calculer; // tous les incréments calculés double tempsDerniereSortieFilCalcul; // variable stockant le dernier temps de sortie au fil du calcul // définition des flots externes éventuels Tableau noms_fichier; Tableau typeFlotExterne; UtilLecture* entreePrinc; // sauvegarde pour éviter les passages de paramètres streampos debut_increment; // sanvegarde de la position du début d'incrément // lors d'un restart Visualisation visualise; // instance sur les utilitaires de visualisation en vrml Visualisation_maple visualise_maple; // instance sur les utilitaires de visualisation en maple Visualisation_geomview visualise_geomview; // instance sur les utilitaires de visualisation en // geomview Visualisation_Gid visualise_Gid; // instance sur les utilitaires de visualisation en Gid Visualisation_Gmsh visualise_Gmsh; // instance sur les utilitaires de visualisation en Gmsh private: Vecteur * toto_masse; // pointeur intermédiaire utilisée dans la méthode Cal_matrice_masse protected: // compteur d'increments de charge: s'avère dans les faits, utilisés par tous les algorithmes // est renseigné par les classes dérivées int icharge; // --- les énergies et forces généralisées ( utilisées par certains algorithmes ) double E_cin_0,E_cin_t,E_cin_tdt,E_int_t,E_int_tdt,E_ext_t,E_ext_tdt,bilan_E; // différentes énergies et bilan double q_mov_t,q_mov_tdt; // la quantité de mouvement double P_acc_tdt,P_int_tdt,P_ext_tdt,bilan_P_tdt; // différentes puissances et bilan double P_acc_t,P_int_t,P_ext_t,bilan_P_t; // différentes puissances et bilan double E_visco_numerique_t,E_visco_numerique_tdt; // énergie due à la viscosité numérique // === vecteurs Vecteur F_int_t,F_ext_t,F_int_tdt,F_ext_tdt; // forces généralisées int et ext à t et t+dt Vecteur F_totale_tdt; // bilan des forces internes et externes, = F_int_tdt + F_ext_tdt // évite d'avoir à le construire à chaque utilisation: utilisation en post-traitement uniquement Vecteur residu_final; // vecteur intermédiaire de sauvegarde pour la visualisation uniquement Vecteur X_t,X_tdt,delta_X,var_delta_X; // les positions, variations entre t et tdt, et variation // d'une itération à l'autre en implicite, et d'un incrément à l'autre en explicite double delta_t_precedent_a_convergence; // cf. son nom ! Vecteur vitesse_t,vitesse_tdt; // vitesses Vecteur acceleration_t,acceleration_tdt ; // accélérations double E_bulk; // énergie due à l'utilisation du bulk viscosity double P_bulk; // puissance due à l'utilisation du bulk viscosity // .... les mêmes informations individuelle (ind) par maillage List_io < double> E_cin_ind_t,E_cin_ind_tdt,E_int_ind_t,E_int_ind_tdt,E_ext_ind_t,E_ext_ind_tdt,bilan_ind_E; // différentes énergies et bilan List_io < double> q_mov_ind_t,q_mov_ind_tdt; // la quantité de mouvement List_io < double> P_acc_ind_tdt,P_int_ind_tdt,P_ext_ind_tdt,bilan_P_ind_tdt; // différentes puissances et bilan List_io < double> P_acc_ind_t,P_int_ind_t,P_ext_ind_t,bilan_P_ind_t; // différentes puissances et bilan List_io < double> E_visco_numerique_ind_t,E_visco_numerique_ind_tdt; // énergie due à la viscosité numérique List_io < Vecteur> F_int_ind_t,F_ext_ind_t,F_int_ind_tdt,F_ext_ind_tdt; // forces généralisées int et ext à t et t+dt List_io < double> E_bulk_ind; // énergie due à l'utilisation du bulk viscosity List_io < double> P_bulk_ind; // puissance due à l'utilisation du bulk viscosity // liste des variables globales scalaires sous forme d'un arbre pour faciliter la recherche // contient par exemple des pointeurs sur les énergies map < string, const double * , std::less > listeVarGlob; // liste des grandeurs globales vectorielles qui représentent en fait des infos aux noeuds: exe: F_int_t et F_ext_t List_io < TypeQuelconque > listeVecGlob; int deja_lue_entete_parametre; // une variable interne qui permet de dire si l'algo spécifique a déjà lue // l'entete des paramétres ou pas : =0 : pas de procédure spécifique de lecture de paramètre // =1, procédure spécifique, mais pas de lecture effective de l'entête par la procédure spécifique // =2, procédure spécifique et lecture effective de l'entête par la procédure spécifique // variables privées internes servant au pilotage de la convergence int PhaseDeConvergence() const {return phase_de_convergence;}; // pour l'acces en lecture uniquement void Change_PhaseDeConvergence(int phase) {phase_de_convergence=phase;}; #ifdef UTILISATION_MPI // cas d'un calcul parallèle void Passage_indicConvergenceAuxProcCalcul() ; // passage des indicateurs aux process de calcul // globalisation des grandeurs globales: proc de calcul -> maitre puis transmission à ParaGlob // et transmission aux proc de calcul: // on ne demande pas à ParaGlob de faire la transmission, car il ne sait pas ce qu'il transmet // et les infos ne sont pas contigües, le transfert ne sera pas performant void Globalisation_et_transfert_auxProcCalcul_grandeurs_globales(); #endif private : int phase_de_convergence; // =0 indique si l'on est en phase de convergence // =1 on a convergé soit très bien soit correctement // =2 on a convergé, mais avec une mauvaise convergence // =-1 on a diverge dans l'algo du à un nb d'iter trop grand, =1 on a convergé // =-2, =-3, =-4 indique qu'un arrêt a été demandé par la méthode Convergence // = -5 : indique que l'on prend en compte un jacobien négatif sous forme d'une divergence // = -6 : idem pour une variation de jacobien trop grande // = -7 : idem pour une non convergence sur la loi de comportement // = -8 : idem mais la non convergence "couve" depuis plusieurs incréments !! // = -9 : divergence due au fait que la résolution du système linéaires est impossible // = -10 : problème due au chargement // = -11 : on a trouvé un nan ou un nombre infini sur la norme des efforts int et/ou ext int nombre_de_bonnes_convergences; // indique le nombre actuel de bonnes convergences en implicite int nombre_de_mauvaises_convergences; // indique le nombre actuel de mauvaises convergences consécutives en implicite int a_converge; // indique que l'on a convergé int a_converge_iterMoins1; // garde en mémoire la convergence de l'itération précédente, n'est utile // que si on veut une vérification double (deux fois) de la convergence Tableau max_var_residu; // tableau qui garde en souvenir la suite des maxi pour le pilotage int nb_cycle_test_max_var_residu; // utilisé dans Convergence() //--- gestion des stockages // -- une classe de stockage de paramètres pour le changement de largeur de bande class ParaStockage {public: Nb_assemb nbass; int demi,total; ParaStockage():nbass(),demi(0),total(0) {}; }; // dans le cas ou on ne fait que changer la largeur de bande en fonction du contact // avec une numérotation qui reste fixe (et donc des pointeurs d'assemblage fixe) // on se sert de tab_para_stockage pour ne pas recalculer la largeur de bande sans contact // par compte dès que les pointeurs d'assemlage change, il faut remettre à jour // cf. -> Mise_a_jour_Choix_matriciel_contact Tableau < ParaStockage > tab_para_stockage; // variables internes de travail Vecteur v_int; // utilisé par CalEnergieAffichage // cas de la remontée aux contraintes ou déformation, et calcul d'erreur Nb_assemb* cas_ass_sig; // assemblage pour de tous les variables contraintes définit aux noeuds Assemblage* Ass_sig; // " Nb_assemb* cas_ass_sig11; // assemblage pour des variables contraintes sigma11 définit aux noeuds Assemblage* Ass_sig11; // " Nb_assemb* cas_ass_eps; // assemblage pour de tous les variables déformations définit aux noeuds Assemblage* Ass_eps; // " Nb_assemb* cas_ass_eps11; // assemblage pour des variables déformations EPS11 définit aux noeuds Assemblage* Ass_eps11; // " Nb_assemb* cas_ass_err; // assemblage pour des variables d'erreurs définit aux noeuds Assemblage* Ass_err; // " Mat_abstraite* mat_glob_sigeps; // matrice utilisée pour la forme variationnelle du pb Tableau * vglob_sigeps; // tableau de seconds membres globaux Vecteur* vglobal_sigeps; // vecteur global de tous les ddl de sigma ou epsilon int nbddl_sigeps,nbddl_sigeps11; // nombre total de ddl de type sigma ou eps, et nb total uniquement du premier int nbddl_err; // nombre de ddl d'erreur bool avec_remont_sigma,avec_remont_eps,avec_erreur; // indicateurs intermédiaires utilisés par CalculRemont et InitRemont // -- volume matière double volume_total_matiere; // volume totale de matière // -- volume déplacé: cas d'élément 2D dans un calcul 3D // un tableau pour avoir un volume par maillage Tableau vol_total2D_avec_plan_ref; // volume total entre la surface et : yz, xz,xy // -- cas des énergies EnergieMeca energTotal; // énergie totale du système double energHourglass; // énergie totale liée au contrôle d'hourglass EnergieMeca energFrottement; // énergie liée au frottement de contact double energPenalisation; // énergie liée à la mise en place d'une méthode de pénalisation double energStabiliMembBiel; // énergie liée à la stabilisation des membranes et biels éventuelles // ---------- debut amortissement cinétique: => spécifique à un algo --------------- protected: bool amortissement_cinetique; // indique si oui ou non on utilise de l'amortissement cinétique Vecteur X_EcinMax; // le X à l'énergie cinétique maxi (ou moy glissante) private: // deux techniques pour appliquer le test, qui conduit à mettre les vitesses à 0 // 1) on regarde si l'|énergie cinétique| diminue n1 fois (pas forcément consécutives) // ou est inférieure n1 fois fraction * dernier_max_E_cin, si fraction (= test_fraction_energie) est non nulle // NB: le test ne se déclanche que lorsque l'on a scruté au moins n1 valeurs de l'énergie cinétique // 2) on regarde si la moyenne glissante de n2 |énergie cinétique| consécutive diminue n1 fois (pas forcément consécutives) // ou est inférieure n1 fois à fraction * dernier_max_E_cin, si fraction (= test_fraction_energie) est non nulle, // NB: le test ne se déclanche que lorsque l'on a au moins n2 valeurs consécutives disponibles de l'énergie cinétique int max_nb_decroit_pourRelaxDyn; // paramètre n1 int taille_moyenne_glissante; // paramètre n2 double moyenne_glissante, moy_gliss_t; // valeur de la moyenne glissante actuelle et moyennes à t int inter_nb_entre_relax; // indique un nombre mini d'iter avec la précédente relaxation, avant d'autoriser // une nouvelle relaxation Fonction_nD* fct_nD_inter_nb_entre_relax; // si non NULL: donne une valeur qui remplace inter_nb_entre_relax string nom_fct_nD_inter_nb_entre_relax; // nom éventuel de la fonction associée int nb_depuis_derniere_relax; // variable inter qui permet d'appliquer inter_nb_entre_relax Vecteur sauve_E_cin; // sauvegarde des valeurs pour le calcul de la moyenne glissante int indice_EC1; // indice dans sauve_E_cin, de la première valeur d'énergie cinétique, stockée int nb_Ecin_stocker; // indique le nombre courant de valeur d'énergies cinétiques stockées double test_fraction_energie; // si diff de 0, indique la valeur de "fraction" int compteur_decroit_pourRelaxDyn; // compteur associé à max_nb_decroit_pourRelaxDyn // =-1 lorsque la relaxation est suspendue double coef_arret_pourRelaxDyn; // on arrête la relaxation lorsque E_cin_tdt < E_cin_t * coef double dernier_pic,pic_E_cint_t; // dernier pic, et pic à l'instant t d'énergie cinétique enregistrée double max_pic_E_cin; // le max des pic d'énergie cinétique double coef_redemarrage_pourRelaxDyn; // on redémarre la relaxation lorsque // E_cin_tdt > max_pic_E_cin * coef_redemarrage_pourRelaxDyn double max_deltaX_pourRelaxDyn; // le calcul s'arrête lorsque deltaX < max_deltaX_pourRelaxDyn m fois int nb_max_dX_OK_pourRelaxDyn; // = le m de la ligne précédente int nb_dX_OK_pourRelaxDyn; // le nombre de dx OK courant int nb_deb_testfin_pourRelaxDyn; // le nombre mini de passage dans la relaxation, a partir duquel on test la fin int nb_deb_test_amort_cinetique; // le nombre mini d'iteration, a partir duquel on démarre l'algorithme int compteur_test_pourRelaxDyn; // compteur associé avec nb_deb_testfin_pourRelaxDyn int compteur_pic_energie; // compteur absolu numérotant les pic d'énergie, pas sauvegardé, donc mis à 0 à chaque // démarrage de calcul //+ cas de l'amortissement local int amortissement_cinetique_au_noeud; // si diff de 0, indique une technique d'amortissement individuelle à chaque noeud int nb_deb_test_amort_cinetique_noe; // le nombre mini d'iteration, a partir duquel on démarre l'algorithme Vecteur E_cinNoe_tdt,E_cinNoe_t; // énergie cinétique à chaque noeud, à t et tdt int max_nb_decroit_pourRelaxDyn_noe; // idem paramètre n1 sur le globale Tableau compteur_decroit_pourRelaxDyn_noe; // compteur associé à max_nb_decroit_pourRelaxDyn // =-1 lorsque la relaxation est suspendue double coef_arret_pourRelaxDyn_noe; // on arrête la relaxation lorsque eN_cin_tdt < eN_cin_t * coef Vecteur dernier_pic_noe; // pour chaque noeuds: dernier pic, et pic à l'instant t d'énergie cinétique enregistrée Vecteur pic_E_cint_t_noe; // Vecteur max_pic_E_cin_noe; // pour chaque noeud: le max des pic d'énergie cinétique double coef_redemarrage_pourRelaxDyn_noe; // on redémarre la relaxation lorsque // eN_cin_tdt > max_pic_E_cin_noe * coef_redemarrage_pourRelaxDyn_noe int compteur_pic_energie_noe; // compteur absolu numérotant les pic d'énergie, pas sauvegardé, donc mis à 0 à chaque // démarrage de calcul Vecteur Puiss_loc_t,Puiss_loc_tdt; // puissances locales à chaques noeuds (comprend les forces externes et internes) // ---------- fin amortissement cinétique: => spécifique à un algo --------------- // ---------- début amortissement visqueux critique: => spécifique à un algo --------------- int opt_cal_C_critique; // choix pour le type de calcul du quotient de Rayleight double f_; // le coefficient multiplicateur pour la borne maxi de la fréquence mini double ampli_visco; // un coefficient d'amplification arbitraire de la viscosité int type_calcul_visqu_critique; // choix entre différents calcul de la partie visqueuse Mat_abstraite* C_inter; // matrice de travail utilisé dans la méthode CalculEnContinuMatriceViscositeCritique // ---------- fin amortissement visqueux critique: => spécifique à un algo --------------- // -- pour les algo dynamiques, convergences vers une solution statique int arret_a_equilibre_statique ; // indique, si diff de 0, que l'on veut contrôler une convergence // Si = 1: le residu statique est utilise comme critere // si = 2: nécessite que les deux critères précédents soient satisfait // -- modulation de la précision de convergence via une fonction nD string nom_modulation_precision; // nom éventuel de la fonction de modulation Fonction_nD* modulation_precision; // si non nulle => modulation //-- utilisation éventuelle d'une fonction nD pour calculer le type de norme de convergence // fait un peu double emploi avec la modulation de précision, mais est plus souple d'emplois // du fait que l'on peut définir ainsi directement le type de norme qu'on veut utiliser Fonction_nD* fct_norme; protected: Vecteur F_int_tdt_prec, F_ext_tdt_prec; // forces généralisées de l'itération précédente (éventuellement) // vecteur inermédiaire, utilisé par les classes dérivées pour stocker la puissance totale statique: sans accélération et sans visqueux numérique Vecteur * vglob_stat; // dans le cas où on veut une convergence vers la solution statique, sans viscosité dynamique inline int Arret_A_Equilibre_Statique() const {return arret_a_equilibre_statique;}; Vecteur vec_trav; // vecteur de travail, de dim le nb de ddl par exe, utilisé par Pilotage_maxi_X_V // récupération de la précision en cours double Precision_equilibre() const {double pa_precision = pa.Precision(); if (modulation_precision != NULL) // cas où en plus on fait une modulation de précision pa_precision *= (modulation_precision->Valeur_pour_variables_globales())(1); return pa_precision; }; //------- temps cpu ----- on garde des temps spécifiques pour chaque algo // ainsi en sortie on pourra différencier les temps totaux et les temps partiels Temps_CPU_HZpp tempsInitialisation; // lesTempsCpu(1) Temps_CPU_HZpp tempsMiseAjourAlgo; // lesTempsCpu(2) Temps_CPU_HZpp tempsCalEquilibre; // lesTempsCpu(3) Temps_CPU_HZpp tempsRaidSmEner; // lesTempsCpu(4) Temps_CPU_HZpp tempsSecondMembreEnerg; // lesTempsCpu(5) Temps_CPU_HZpp tempsResolSystemLineaire; // lesTempsCpu(6) Temps_CPU_HZpp tempsSauvegarde; // lesTempsCpu(7) Temps_CPU_HZpp tempsSortieFilCalcul; // lesTempsCpu(8) Temps_CPU_HZpp tempsRaidSmEnerContact; // lesTempsCpu(9) Temps_CPU_HZpp tempsSecondMembreEnergContact; // lesTempsCpu(10) Temps_CPU_HZpp temps_CL; // lesTempsCpu(11) Temps_CPU_HZpp temps_CLL; // lesTempsCpu(12) Temps_CPU_HZpp temps_lois_comportement; // lesTempsCpu(13) Temps_CPU_HZpp temps_metrique_K_SM; // lesTempsCpu(14) Temps_CPU_HZpp temps_chargement; // lesTempsCpu(15) Temps_CPU_HZpp temps_rech_contact; // lesTempsCpu(16) #ifdef UTILISATION_MPI // cas d'un calcul parallèle: // passage des infos entre process Temps_CPU_HZpp temps_transfert_court_algo ; // lesTempsCpu(17) Temps_CPU_HZpp temps_transfert_long_algo ; // lesTempsCpu(18) Temps_CPU_HZpp temps_attente_algo ; // lesTempsCpu(19) Temps_CPU_HZpp temps_transfert_court_matSm ; // lesTempsCpu(20) Temps_CPU_HZpp temps_transfert_long_matSm ; // lesTempsCpu(21) Temps_CPU_HZpp temps_attente_matSm ; // lesTempsCpu(22) Temps_CPU_HZpp temps_transfert_court_charge ; // lesTempsCpu(23) Temps_CPU_HZpp temps_transfert_long_charge ; // lesTempsCpu(24) Temps_CPU_HZpp temps_attente_charge ; // lesTempsCpu(25) Temps_CPU_HZpp temps_transfert_court_contact ; // lesTempsCpu(26) Temps_CPU_HZpp temps_transfert_long_contact ; // lesTempsCpu(27) Temps_CPU_HZpp temps_attente_contact ; // lesTempsCpu(28) #endif Tableau lesTempsCpu; // un tableau intermédiaire qui récupère et globalise les temps pour les sorties // via listeVarGlob, mais c'est les variables Temps_CPU_HZpp qui stockent // réellement les temps // la petite méthode qui sert au transfert et qui doit-être appelé avant les sorties void Temps_CPU_HZpp_to_lesTempsCpu(const LesCondLim& lesCondLim, const Charge& charge,const LesContacts & contact); //========== cas particulier des algo combiner ============================================ // un pointeur dans le cas où les temps à sortir ne sont pas ceux stockés // --> utilisé uniquement par l'algo combiner, car dans ce cas il faut globaliser les temps // de tous les sous-algo // par défaut ce pointeur est NULL AlgoriCombine* ptalgocombi; // et la méthode pour changer void Change_ptalgocombi (AlgoriCombine* ptal) {ptalgocombi = ptal;}; // idem la méthode de transfert si-dessus mais concerne uniquement les temps internes à l'algo // ajoute au tableau passé en paramètre, les temps de l'algo Tableau & Ajout_Temps_CPU_HZpp_to_lesTempsCpu(Tableau & lesTsCpu); //========== fin cas particulier des algo combiner ============================================ // --- METHODES PRIVEES // METHODES PROTEGEES : protected: // -- pour les algo dynamiques, convergences vers une solution statique const int& ArretEquilibreStatique()const {return arret_a_equilibre_statique;} ; // indique, si diff de 0, que l'on veut contrôler une convergence vers la solution // initialisation du calcul de la remontée des contraintes aux noeuds // initialisation valable également dans le cas où aucun calcul n'a précédé ce calcul void InitRemontSigma(LesMaillages * lesMail,LesReferences* , DiversStockage* ,Charge* , LesCondLim* ,LesContacts* ,Resultats* ); // initialisation du calcul de la remontée des déformations aux noeuds // initialisation valable également dans le cas où aucun calcul n'a précédé ce calcul void InitRemontEps(LesMaillages * lesMail,LesReferences* , DiversStockage* ,Charge* , LesCondLim* ,LesContacts* ,Resultats* ); // initialisation du calcul d'erreur aux éléments // initialisation valable également dans le cas où aucun calcul n'a précédé ce calcul void InitErreur(LesMaillages * lesMail,LesReferences* lesRef, DiversStockage* toto,Charge* charge, LesCondLim* lesCondLim,LesContacts* titi,Resultats* tutu); // calcul de la remontée des contraintes aux noeuds void RemontSigma(LesMaillages * lesMail); // calcul de la remontée des déformations aux noeuds void RemontEps(LesMaillages * lesMail); // calcul de l'erreur et de sa remonté aux noeuds, après un calcul de mécanique. // nécessite d'avoir fait auparavant la remonté des contraintes (initialisation et remontée !!) void RemontErreur(LesMaillages * lesMail); // initialisation d'un calcul de remonté à des variables secondaires ou d'erreur // globalement, cette initialisation dépend des paramètres utilisateurs // ramène true s'il y a initialisation bool InitRemont(LesMaillages * lesMail,LesReferences* lesRef, DiversStockage* ,Charge* charge, LesCondLim* lesCondLim,LesContacts* ,Resultats* ); // calcul d'une remontée à des variables secondaires ou d'erreur, ceci pour un post traitement de visualisation // ce calcul dépend des paramètres utilisateurs, il n'est fait qu'à chaque sortie .BI // ramène true s'il y a un calcul effectif bool CalculRemont(LesMaillages * lesMail,OrdreVisu::EnumTypeIncre type_incre,int incre); // lecture d'un type externe void LectureUnTypeExterne(UtilLecture* entreePrinc,const int type); // lecture des paramètres du calcul généraux pour tous les algos (dans le cas où il y en a) // qui peuvent éventuellement ne pas servir !! virtual void lecture_Parametres(UtilLecture& entreePrinc); // création d'un fichier de commande: cas des paramètres spécifiques // par défaut on indique qu'il n'y a pas de paramètres, // cette méthode est appelée par les classes dérivées void Info_com_parametres(UtilLecture& entreePrinc); // algorithme classique de line search // le line seach agit sur une partie de la puissance et non sur la norme du résidu // le line search à pour objectif de minimiser une partie de l'acroissement de la puissance // en jeux, ainsi on peut avoir une minimisation réussi tout en ayant une norme d'équilibre // qui continue à grandir // cependant en général les deux sont lié et l'application du line_seach garanti que dans // la direction de descente choisi on minimise norme (mais ce minimum peut cependant être // supérieur à l'ancienne norme!) // la méthode renvoi si oui ou non le line search a été effecué bool Line_search1(Vecteur& sauve_deltadept,double& puis_precedente, Vecteur& Vres,LesMaillages * lesMail,Vecteur* sol ,int& compteur,Vecteur& sauve_dept_a_tdt,Charge* charge,Vecteur& vglobex ,Assemblage& Ass,Vecteur& v_travail,LesCondLim* lesCondLim,Vecteur& vglobal ,LesReferences* lesRef,Vecteur & vglobin,const Nb_assemb& nb_casAssemb ,int cas_combi_ddl,LesCourbes1D* lesCourbes1D,LesFonctions_nD* lesFonctionsnD); // algorithme de line search qui utilise une approximation cubique du résidu ( (R).) // en fonction du paramètre lambda de line search ( X+lambda. delta_ddl) // la méthode renvoi si oui ou non le line search a été effecué bool Line_search2(Vecteur& sauve_deltadept,double& puis_precedente, Vecteur& Vres,LesMaillages * lesMail,Vecteur* sol ,int& compteur,Vecteur& sauve_dept_a_tdt,Charge* charge,Vecteur& vglobex ,Assemblage& Ass,Vecteur& v_travail,LesCondLim* lesCondLim,Vecteur& vglobal ,LesReferences* lesRef,Vecteur & vglobin,const Nb_assemb& nb_casAssemb ,int cas_combi_ddl,LesCourbes1D* lesCourbes1D,LesFonctions_nD* lesFonctionsnD); // calcul de la raideur et du second membre pour tous les maillages ainsi que des énergies // si pb retour false: indique que le calcul c'est mal passé, il y a une erreur // généré, les résultats: raideur et second membres calculées sont inexploitable bool RaidSmEner(LesMaillages * lesMail,Assemblage& Ass,Vecteur & vglobin ,Mat_abstraite& matglob ); // calcul du second membre pour tous les maillages ainsi que des énergies // si pb retour false: indique que le calcul c'est mal passé, il y a une erreur // généré, les résultats: raideur et second membres calculées sont inexploitable bool SecondMembreEnerg(LesMaillages * lesMail,Assemblage& Ass,Vecteur & vglobin); // Mise à jour de la raideur et du second membre pour le contact // calcul des énergies spécifiques au contact // si pb retour false: indique que le calcul c'est mal passé, il y a une erreur // généré, les résultats: raideur et second membres calculées sont inexploitables bool RaidSmEnerContact(LesContacts * lesCont,Assemblage& Ass,Vecteur & vglobin ,Mat_abstraite& matglob ); // cas du contact: calcul et mise à jour du second membre ainsi que des énergies spécifiques bool SecondMembreEnergContact(LesContacts * lesContacts,Assemblage& Ass,Vecteur & vglobin ,bool aff_iteration); // choix de la (les) matrices de raideur du système linéaire // 1) Si le pointeur est non null, on vérifie les tailles, si c'est différent // la matrice est reconstruite // 2) Si le pointeur est null, il y a création d'une matrice // lesCondLim est utilisé pour modifier la largeur de bande par exemple // en fonction des conditions limites linéaire // lescontacts: si différent de NULL indique qu'il y a du contact // NB: la dimension du tableau dépend du paramétrage lu dans ParaAlgoGlobal Tableau Choix_matriciel (int nbddl,Tableau & tab_matglob,LesMaillages * lesMail ,LesReferences* lesRef,const Nb_assemb& nb_casAssemb,LesCondLim* lesCondLim , LesContacts*lescontacts = NULL); // Mise à jour éventuelle de la taille de la matrice de raideur ou masse du système linéaire // en fonction du contact ou en fonction de "nouvelle_largeur_imposee" dans ce cas, lescontacts // n'est pas utilisé ! donc la méthode peut servir dans ce cas sans le contact // // la matrice *matglob: doit déjà existée // en retour le pointeur désigne toujours la même matrice qu'en entrée // avec les même caractéristique de fonctionnement que l'ancienne (mêm résolution, stockage ...) // NB: la dimension du tableau dépend du paramétrage lu dans ParaAlgoGlobal // 1) si niveau_substitution =0 : -> mise à jour de toutes les matrices // 2) si niveau_substitution = i : -> uniquement mise à jour de la matrices i // par contre la dimension de tab_matglob est toujours mis à jour si c'est nécessaire // si: nouvelle_largeur_imposee : n'est pas nulle, // alors c'est directement cette valeur en noeud qui est imposée sur la matrice // nouvelle_largeur_imposee.un = la largeur totale // nouvelle_largeur_imposee.deux = la demie largeur // nouvelles_largeur_en_ddl.trois = la demie largeur maximale pour la partie éléments finis // uniquement (sans les CLL) void Mise_a_jour_Choix_matriciel_contact (Tableau & tab_matglob,const Nb_assemb& nb_casAssemb , LesContacts*lescontacts ,int niveau_substitution ,TroisEntiers* nouvelle_largeur_imposee = NULL); // choix de la matrice de la matrice de masse globale // ou mise à jour de la matrice masse si elle existe Mat_abstraite* Choix_matrice_masse (int nbddl,Mat_abstraite* matglob,LesMaillages * lesMail,LesReferences* lesRef ,const Nb_assemb& nb_casAssemb, LesContacts* lescontacts,LesCondLim* lesCondLim); // mise à jour éventuelle du type et/ou de la taille de la matrice de masse // retourne un pointeur sur la nouvelle matrice, s'il y a eu une modification // l'ancienne est supprimée // sinon retour d'un pointeur NULL Mat_abstraite* Mise_a_jour_type_et_taille_matrice_masse_en_explicite (int nbddl,Mat_abstraite* matglob,LesMaillages * lesMail,LesReferences* lesRef ,const Nb_assemb& nb_casAssemb, LesContacts* lescontacts); // calcul la matrice de masse, y compris les masses ponctuelles si elles existent // N_ddl : donne le ddl où doit être mis la masse void Cal_matrice_masse(LesMaillages * lesMail,Assemblage& Ass,Mat_abstraite& matglob ,const DiversStockage* diversStockage,LesReferences* lesRef ,const Enum_ddl & N_ddl,LesFonctions_nD* lesFonctionsnD); // uniquement ajout dans la matrice de masse des masses ponctuelles si elles existent // N_ddl : donne le ddl où doit être mis la masse // par exemple, sert pour la relaxation dynamique void Ajout_masses_ponctuelles(LesMaillages * lesMail,Assemblage& Ass,Mat_abstraite& matglob ,const DiversStockage* diversStockage,LesReferences* lesRef ,const Enum_ddl & N_ddl,LesFonctions_nD* lesFonctionsnD); // pilotage de la convergence, ceci en fonction du type de pilotage choisit // ramène true si le pas de temps a été modifié, false sinon // arret = true: indique dans le cas d'une non convergence, neanmoins rien n'a changé, il faut donc arrêter bool Pilotage_du_temps(Charge* charge,bool& arret); // pilotage de la fin des itération en implicite : utilisation conjointe avec: Pilotage_du_temps() // retourne true dans le cas normale bool Pilotage_fin_iteration_implicite(int compteur); // examine s'il y a convergence on pas en fonction des parametres // de controle et du resultat du systeme lineaire // affiche: indique si l'on veut l'affichage ou non des infos // itera : correspond au compteur d'itération dans un cas d'équilibre itératif // en sortie: ramene le maxi du residu et le pointeur d'assemblage correspondant // arrêt: indique qu'il vaut mieux arrêter le calcul bool Convergence(bool affiche, double last_var_ddl_max,Vecteur& residu,double maxPuissExt,double maxPuissInt ,double maxReaction,int itera,bool& arret); // pilotage spécifique de la fin de la relaxation dynamique, en tenant compte des différents critères demandés // relax_vit_acce : // = 1 : il y a eu relaxation mais le calcul continue // = 0 : pas de relaxation // = -1 : on peut arrêter le calcul car le critère sur delta_ddl est satisfait // compteur: indique le nombre d'itération en cours // arretResidu: =true : le critère sur le résidu est ok // iter_ou_incr : =1: il s'agit d'itération", =0 il s'agit d'incrément // lit arret, puis modifie arret en true ou false = il faut arrêter ou non // si true met à jour en interne certains paramètres de gestion de pilotage pour dire que l'on a bien convergé bool Pilotage_fin_relaxation_et_ou_residu(const int& relax_vit_acce,const int& iter_ou_incr, const int& compteur, const bool& arretResidu, bool& arret); // pilotage à chaque itération: // - sous ou sur relaxation // - limitation de la norme de delta ddl // - indication de précision pour les éléments (utilisée pour certain type de gestion d'hourglass) // en entrée: le maxi de var ddl juste après résolution // ramène le maxi de var ddl maxDeltaDdl modifié, ainsi que sol modifié éventuellement, void Pilotage_chaque_iteration(Vecteur* sol,double& maxDeltaDdl,const int& itera,LesMaillages * lesMail); // pilotage pour l'initialisation de Xtdt à chaque début d'incrément, en implicite // avec la même direction qu'au pas précédent // ramène true, si c'est ok pour initialiser Xtdt bool Pilotage_init_Xtdt(); // pilotage en dynamique implicite, du maxi des déplacements et/ou vitesses // limitation spécifiquement des maxi en déplacements et/ou vitesses // il y a calcul éventuel (pas toujours) de delta_X, mais c'est considéré comme une variable de travail void Pilotage_maxi_X_V(const Vecteur& X_t,Vecteur& X_tdt,const Vecteur& V_t,Vecteur& V_tdt); // calcul une approximation des différences énergies et puissances en jeux // affichage éventuelle des ces énergies et du bilan // V : ddl de vitesse // coef_mass : en fait il faut utiliser 1/coef_mass * mat_mass pour la matrice masse, ceci due à la spécificité de l'algo // icharge : compteur d'increment // brestart : booleen qui indique si l'on est en restart ou pas // gamma : vecteur accélération void CalEnergieAffichage(const double& coef_mass,const Vecteur & V,const Mat_abstraite& mat_mass ,const Vecteur & delta_ddl,int icharge,bool brestart,const Vecteur & gamma ,const Vecteur & forces_vis_num); // idem pour un calcul statique, c'est-à-dire sans énergie cinématique et puissance d'accélération // affichage éventuelle des ces énergies et du bilan // icharge : compteur d'increment // brestart : booleen qui indique si l'on est en restart ou pas void CalEnergieAffichage(const Vecteur & delta_ddl,int icharge,bool brestart,const Vecteur & forces_vis_num); // affichage éventuelle de la matrice de raideur et du second membre void Affiche_RaidSM(const Vecteur & vglobin,const Mat_abstraite& matglob) const; // méthode d'application de l'amortissement cinétique // retour un int : // = 1 : il y a eu amortissement mais le calcul continue // = 0 : pas d'amortissement // = -1 : on peut arrêter le calcul car le critère sur delta_ddl est satisfait int AmortissementCinetique(const Vecteur & delta_ddl,const double& coef_mass,const Vecteur& X ,const Mat_abstraite& mat_mass,int icharge,Vecteur& V); // initialisation des paramètres pour l'amortissement cinétique (au cas ou on veut plusieurs amortissement) void InitialiseAmortissementCinetique(); // définition de la matrice de viscosité numérique // inita : true-> indique si l'on est dans une phase d'initialisation (creation de la matrice) ou non // dans ce dernier cas il y a modification éventuelle des valeurs de C // ramène un pointeur sur la matrice visqueuse instancié ou nul si aucune instantiation Mat_abstraite* Cal_mat_visqueux_num_expli(const Mat_abstraite& mat_mass,Mat_abstraite* mat_C_pt ,const Vecteur & delta_X,bool inita, const Vecteur & vitesse); // calcul de l'impact d'une viscosité artificielle sur la raideur et sur le second membre // l'algo travaille avec les vecteurs: delta_X, var_delta_X, et vitesse stockés dans Algori // l'algo modifie mat_glob et calcul les forces visqueuses void Cal_mat_visqueux_num_stat(Mat_abstraite& mat_glob,Vecteur& forces_vis_num); // mise à jour qui sont gérées par la classe mère algo // a priori que des choses génériques du type gestion des variables privées void MiseAJourAlgoMere(ParaGlob * paraGlob,LesMaillages * lesMail,LesReferences* lesRef ,LesCourbes1D* lesCourbes1D,LesFonctions_nD* lesFonctionsnD ,VariablesExporter* varExpor ,LesLoisDeComp* lesLoisDeComp,DiversStockage* diversStockage ,Charge* charge,LesCondLim* lesCondLim,LesContacts* lesContacts ,Resultats* resultats); // gestion éventuelle d'une renumérotation, qui prend en compte les éléments de contact // premier_calcul : indique s'il s'agit d'un premier calcul on non // nouvelle_situation_contact : indique s'il y a du nouveau sur le contact // 1) si niveau_substitution =0 : -> mise à jour de toutes les matrices // 2) si niveau_substitution = i : -> uniquement mise à jour de la matrices i // par contre la dimension de tab_matglob est toujours mis à jour si c'est nécessaire // ramène true si la matrice a été changée bool Gestion_stockage_et_renumerotation_avec_contact(bool premier_calcul ,LesMaillages * lesMail, bool & nouvelle_situation_contact ,LesCondLim* lesCondLim,LesReferences* lesRef ,Tableau & tab_matglob,const Nb_assemb& nb_casAssemb ,LesContacts*lescontacts,int niveau_substitution); // --- même chose mais sans le contact, par contre prise en compte d'un changement éventuelle // imposé par exemple par les CLL // gestion éventuelle d'une renumérotation des pointeurs d'assemblage , qui prend en compte les CLL // premier_calcul : indique s'il s'agit d'un premier calcul on non // nouvelle_situation_CLL : indique s'il y a du nouveau sur les CLL // -> la renumérotation des pointeurs s'effectue que si: // a) ParaAlgoControleActifs().Optimisation_pointeur_assemblage() est ok // b) soit c'est un premier calcul, soit il y a du nouveau sur les CLL // // 1) si niveau_substitution =0 : -> mise à jour de toutes les matrices // 2) si niveau_substitution = i : -> uniquement mise à jour de la matrices i // par contre la dimension de tab_matglob est toujours mis à jour si c'est nécessaire // ramène true si la matrice a été changée // NB: lescontacts est quand même passé en paramètre, car on doit le renseigner au niveau d'un // tableau d'indice (num de noeuds) lorsque les num ont changés, ceci pour être opérationnel // par la suite si le contact devient actif bool Gestion_stockage_et_renumerotation_sans_contact(LesContacts* lescontacts,bool premier_calcul ,LesMaillages * lesMail, bool & nouvelle_situation_CLL ,LesCondLim* lesCondLim,LesReferences* lesRef ,Tableau & tab_matglob,const Nb_assemb& nb_casAssemb ,int niveau_substitution); private: // mise à jour de la viscosité critique en continu: spécifique au algo de relaxation dynamique, appelé par Cal_mat_visqueux_num_expli void CalculEnContinuMatriceViscositeCritique(const Mat_abstraite& mat_mass,Mat_abstraite& mat_C_pt ,const Vecteur & delta_X, const Vecteur & vitesse); protected: // passage des grandeurs gérées par l'algorithme de tdt à t // en particulier les énergies et les puissances void TdtversT(); // mise à jour de delta_X, var_delta_X et passage en global des maxi void Cal_Transfert_delta_et_var_X(double& max_delta_X, double& max_var_delta_X); // vérification d'une singularité éventuelle de la matrice de raideur de dim: nbddl // pour un problème de mécanique, en comparant le nombre de point d'intégration total // et le nombre totale de degré de liberté void VerifSingulariteRaideurMeca(int nbddl,const LesMaillages& lesMail) const; // fonction virtuelle permettant de mettre à jour les infos aux noeuds // à cause de certains contact (ex: cas_contact = 4) // Les vecteurs Xtdt et Vtdt sont mis à jour par la méthode // la vitesse locale du noeud est mis à jour en fonction de la position locale et de l'algo virtual void Repercussion_algo_sur_cinematique(LesContacts* lesContacts,Vecteur& Xtdt,Vecteur& Vtdt) {}; // retouve le numéro de l'incrément sauvegardé trouvé dont le numéro est juste inf ou égale au // dernier incrément calculé - nb_incr_en_arriere. // ramène false si on n'a pas trouvé d'incrément a-doc // en retour icharge prend la valeur de l'incrément trouvé (si tout c'est bien passé) bool Controle_retour_sur_un_increment_enregistre(int nb_incr_en_arriere,int& icharge); // méthode interne d'application de l'amortissement cinétique, exploratoire // cas d'un amortissement local // retour un int : // = 1 : il y a eu amortissement mais le calcul continue // = 0 : pas d'amortissement // = -1 : on peut arrêter le calcul car le critère sur delta_ddl est satisfait int AmortissementCinetique_individuel_aux_noeuds(const Vecteur & delta_ddl,const double& coef_mass,const Vecteur& X ,const Mat_abstraite& mat_mass,int icharge,Vecteur& V); // passage aux noeuds éventuellement des grandeurs globales, pour une sortie de post-traitement par exemple mais pas seulement // le principe est que ce passage s'effectue si les conteneurs existent au niveau des noeuds void Passage_aux_noeuds_grandeurs_globales(LesMaillages * lesMail); // passage aux noeuds de F_int_t et F_ext_t void Passage_aux_noeuds_F_int_t_et_F_ext_t(LesMaillages * lesMail); // passage aux noeuds éventuellement d'une grandeur globale particulière, // typeGeneriqu : indique le type de grandeur à passer // seules les grandeurs globales qui n'ont pas été transférées par la méthode // Passage_aux_noeuds_grandeurs_globales , sont concernées // Le passage s'effectue si le conteneur existe au niveau des noeuds sinon erreur void Passage_aux_noeuds_grandeur_globale_particuliere(const TypeQuelconque& typeGeneriqu, LesMaillages * lesMail); // passage des grandeurs globales aux noeuds où il y a des variables globales attachées // nb_casAssemb correspond au cas d'assemblage de X1 void Passage_de_grandeurs_globales_vers_noeuds_pour_variables_globales (LesMaillages * lesMail,VariablesExporter* varExpor,const Nb_assemb& nb_casAssemb,const LesReferences& lesRef); // des fonctions inlines pour mettre à jour des grandeurs globales // -- initialisation du compteur d'itérations void Transfert_ParaGlob_COMPTEUR_ITERATION_ALGO_GLOBAL(int compteur) const {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(COMPTEUR_ITERATION_ALGO_GLOBAL); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_entier& gr = *((Grandeur_scalaire_entier*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurEntier()) = compteur; }; // --- cas de la norme de convergence void Transfert_ParaGlob_NORME_CONVERGENCE(double laNorme) const {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(NORME_CONVERGENCE); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = laNorme; }; // --- cas du compteur d'incrément void Transfert_ParaGlob_COMPTEUR_INCREMENT_CHARGE_ALGO_GLOBAL(int incre) const {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(COMPTEUR_INCREMENT_CHARGE_ALGO_GLOBAL); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_entier& gr = *((Grandeur_scalaire_entier*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurEntier()) = incre; }; // --- cas de l'algorithme global actuellement en service void Transfert_ParaGlob_ALGO_GLOBAL_ACTUEL(EnumTypeCalcul type_algo) const {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ALGO_GLOBAL_ACTUEL); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_entier& gr = *((Grandeur_scalaire_entier*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurEntier()) = type_algo; }; // --- cas des énergies internes venants des lois de comportement : // ENERGIE_ELASTIQUE,ENERGIE_PLASTIQUE,ENERGIE_VISQUEUSE, void Transfert_ParaGlob_energies_internesLoisComp() const {{void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_ELASTIQUE); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = energTotal.EnergieElastique(); }; {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_PLASTIQUE); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = energTotal.DissipationPlastique(); }; {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_VISQUEUSE); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = energTotal.DissipationVisqueuse(); }; }; // --- cas des énergies de contact : // ENERGIE_PENALISATION,ENERGIE_FROT_ELAST,ENERGIE_FROT_PLAST,ENERGIE_FROT_VISQ void Transfert_ParaGlob_energies_contact() const {{void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_PENALISATION); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = energPenalisation; }; {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_FROT_ELAST); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = energFrottement.EnergieElastique(); }; {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_FROT_PLAST); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = energFrottement.DissipationPlastique(); }; {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_FROT_VISQ); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = energFrottement.DissipationVisqueuse(); }; }; // --- cas des énergies non physiques, d'origine numérique : // ENERGIE_BULK_VISCOSITY,ENERGIE_HOURGLASS_, stabilisation de membrane ou et biel void Transfert_ParaGlob_energies_hourglass_bulk_stab() const {{void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_BULK_VISCOSITY); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = E_bulk; }; {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_HOURGLASS_); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = energHourglass; }; {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(ENERGIE_STABILISATION_MEMB_BIEL); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = energStabiliMembBiel; }; }; // --- cas des volumes entre plans : void Transfert_ParaGlob_volume_entre_plans() const {// on met à jour les grandeurs que si on a fait la demande de calcul // sinon cela ne sert à rien car "vol_total2D_avec_plan_ref" n'est pas mis à jour par ailleurs if (!(pa.CalVolTotalEntreSurfaceEtPlansRef())) return; switch (ParaGlob::Dimension()) { case 3: {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(VOL_TOTAL2D_AVEC_PLAN_XY); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = vol_total2D_avec_plan_ref(1)(3); }; case 2: {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(VOL_TOTAL2D_AVEC_PLAN_YZ); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = vol_total2D_avec_plan_ref(1)(1); }; case 1: {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable(VOL_TOTAL2D_AVEC_PLAN_XZ); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = vol_total2D_avec_plan_ref(1)(2); }; default: break; }; // si on a plus d'un maillage on ajoute les autres valeurs int nbMailMax = vol_total2D_avec_plan_ref.Taille(); if (nbMailMax > 1) {for (int nbMail=1;nbMail< nbMailMax+1;nbMail++) { switch (ParaGlob::Dimension()) { case 3: {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable ("vol_total2D_avec_plan_xy_"+ChangeEntierSTring(nbMail)); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = vol_total2D_avec_plan_ref(nbMail)(3); }; case 2: {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable ("vol_total2D_avec_plan_yz_"+ChangeEntierSTring(nbMail)); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = vol_total2D_avec_plan_ref(nbMail)(1); }; case 1: {void* pt_void = ParaGlob::param->Mise_a_jour_grandeur_consultable ("vol_total2D_avec_plan_xz_"+ChangeEntierSTring(nbMail)); TypeQuelconque* pt_quelc = (TypeQuelconque*) pt_void; Grandeur_scalaire_double& gr = *((Grandeur_scalaire_double*) pt_quelc->Grandeur_pointee()); // pour simplifier *(gr.ConteneurDouble()) = vol_total2D_avec_plan_ref(nbMail)(2); }; default: break; }; }; }; }; #ifdef UTILISATION_MPI // cas d'un calcul parallèle, passage des énergies, volumes // a priori utilisée dans le calcul de raideur et second membres void Passage_energiesEtVolumes(); #endif }; /// @} // end of group #include "Charge.h" #include "Resultats.h" #endif