259 lines
13 KiB
C++
Executable file
259 lines
13 KiB
C++
Executable file
// This file is part of the Herezh++ application.
|
|
//
|
|
// The finite element software Herezh++ is dedicated to the field
|
|
// of mechanics for large transformations of solid structures.
|
|
// It is developed by Gérard Rio (APP: IDDN.FR.010.0106078.000.R.P.2006.035.20600)
|
|
// INSTITUT DE RECHERCHE DUPUY DE LÔME (IRDL) <https://www.irdl.fr/>.
|
|
//
|
|
// Herezh++ is distributed under GPL 3 license ou ultérieure.
|
|
//
|
|
// Copyright (C) 1997-2021 Université Bretagne Sud (France)
|
|
// AUTHOR : Gérard Rio
|
|
// E-MAIL : gerardrio56@free.fr
|
|
//
|
|
// This program is free software: you can redistribute it and/or modify
|
|
// it under the terms of the GNU General Public License as published by
|
|
// the Free Software Foundation, either version 3 of the License,
|
|
// or (at your option) any later version.
|
|
//
|
|
// This program is distributed in the hope that it will be useful,
|
|
// but WITHOUT ANY WARRANTY; without even the implied warranty
|
|
// of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
|
|
// See the GNU General Public License for more details.
|
|
//
|
|
// You should have received a copy of the GNU General Public License
|
|
// along with this program. If not, see <https://www.gnu.org/licenses/>.
|
|
//
|
|
// For more information, please consult: <https://herezh.irdl.fr/>.
|
|
|
|
/************************************************************************
|
|
* DATE: 11/10/2003 *
|
|
* $ *
|
|
* AUTEUR: G RIO (mailto:gerardrio56@free.fr) *
|
|
* $ *
|
|
* PROJET: Herezh++ *
|
|
* $ *
|
|
************************************************************************
|
|
* BUT: Algorithmes de base pour la recherche de zéro d'une *
|
|
* fonction ou d'un ensemble de fonctions. *
|
|
* $ *
|
|
* '''''''''''''''''''''''''''''''''''''''''''''''''''''''''''' *
|
|
* $ *
|
|
************************************************************************/
|
|
#ifndef ALGO_ZERO_H
|
|
#define ALGO_ZERO_H
|
|
|
|
#include "Vecteur.h"
|
|
#include "Mat_abstraite.h"
|
|
#include "Racine.h"
|
|
#include <math.h>
|
|
|
|
/// @addtogroup Les_classes_algo
|
|
/// @{
|
|
///
|
|
|
|
|
|
/**
|
|
*
|
|
* BUT: pour la gestion d'exception pour non convergence
|
|
*
|
|
*
|
|
* \author Gérard Rio
|
|
* \version 1.0
|
|
* \date 11/10/2003
|
|
* \brief pour la gestion d'exception pour non convergence
|
|
*
|
|
*/
|
|
|
|
class ErrNonConvergence_Newton
|
|
// =0 cas courant, pas d'information particulière
|
|
// =1 cas où l'erreur est sévère et ne pourra pas être corrigé en refaisant un calcul avec un
|
|
// pas de temps plus petit. Il faut refaire le calcul en se positionnant plusieurs pas de temps
|
|
// auparavant (utilisé par l'hystérésis par exemple)
|
|
{ public :
|
|
int cas;
|
|
ErrNonConvergence_Newton () : cas(0) {} ; // par défaut
|
|
ErrNonConvergence_Newton (int ca) : cas(ca) {} ; // pb
|
|
};
|
|
/// @} // end of group
|
|
|
|
/// @addtogroup Les_classes_algo
|
|
/// @{
|
|
///
|
|
|
|
|
|
/**
|
|
*
|
|
* BUT: Algorithmes de base pour la recherche de zéro d'une
|
|
* fonction ou d'un ensemble de fonctions.
|
|
*
|
|
*
|
|
* \author Gérard Rio
|
|
* \version 1.0
|
|
* \date 11/10/2003
|
|
* \brief Algorithmes de base pour la recherche de zéro d'une fonction ou d'un ensemble de fonctions.
|
|
*
|
|
*/
|
|
|
|
class Algo_zero
|
|
{
|
|
public :
|
|
// CONSTRUCTEURS :
|
|
/// constructeur par défaut
|
|
Algo_zero();
|
|
/// constructeur de copie
|
|
Algo_zero(const Algo_zero& a);
|
|
/// DESTRUCTEUR :
|
|
~Algo_zero();
|
|
|
|
// METHODES PUBLIQUES :
|
|
|
|
/// recherche du zéro d'une équation du second degré dans l'espace des réels: ax^2+bx+c=0.
|
|
/// cas indique les différents cas:
|
|
/// = 1 : il y a deux racines , racine2 > racine1
|
|
/// = 2 : il y a une racine double
|
|
/// = 0 : pas de racine
|
|
/// = -1 : pas de racine car a,b,c sont tous inférieur à la précision de traitement
|
|
/// = -2 : pas de racine car a,b sont inférieur à la précision de traitement tandis que c est non nul
|
|
/// = -3 : une racine simple car a est inférieur à la précision de traitement, donc considéré comme nul
|
|
/// b et c sont non nul -> equa du premier degré, solution: racine1
|
|
void SecondDegre(double a, double b, double c, double& racine1, double& racine2, int& cas) const ;
|
|
|
|
/// recherche du zéro d'une équation du troisième degré dans l'espace des réels: ax^3+bx^2+cx+d=0.
|
|
/// cas indique les différents cas:
|
|
/// = 1 : il y a deux racines , racine2 > racine1
|
|
/// = 2 : il y a une racine double
|
|
/// = 3 : il y a une seule racine réelle racine1 (et deux complexes non calculées)
|
|
/// = 4 : il y a une racine simple: racine1, et une racine double: racine 2
|
|
/// = 5 : il y a 3 racines réelles: racine1, racine2, racine3
|
|
/// = 0 : pas de racine
|
|
/// = -1 : pas de racine car a,b,c sont tous inférieur à la précision de traitement
|
|
/// = -2 : pas de racine car a,b sont inférieur à la précision de traitement tandis que c est non nul
|
|
/// = -3 : une racine simple car a est inférieur à la précision de traitement, donc considéré comme nul
|
|
/// b et c sont non nul -> equa du premier degré, solution: racine1
|
|
void TroisiemeDegre(double a, double b, double c, double d, double& racine1
|
|
, double& racine2, double& racine3, int& cas) ;
|
|
|
|
///recherche du zéro d'une fonction en utilisant la méthode de Newton-Raphson incrémentale
|
|
/// ici il s'agit d'une fonction à une variable
|
|
/// *Ptfonc : le pointeur de la fonction dont il faut chercher le zéro
|
|
/// *Ptder_fonc : pointeur de la dérivée de la fonction
|
|
/// pour ces deux fonctions, la variable alpha est un facteur de charge qui varie de
|
|
/// 0 à 1. Lorsque alpha=0, la racine vaut val_initiale, et ce que l'on
|
|
/// cherche c'est la racine pour alpha = 1. Lorsque alpha varie progressivement
|
|
/// de 0 à 1, la racine est sensée varier progressivement de val_initiale à résidu
|
|
/// l'utilisation d'alpha permet de faire du Newton incrémentale.
|
|
/// pour ces deux fonctions, l'argument test ramène
|
|
/// . 1 si le calcul a été ok, -1 s'il y a eu un pb, mais on peut continuer, 0 s'il y a eu un pb
|
|
/// fatal, qui invalide le calcul de la fonction et/ou de la dérivée
|
|
/// val_initiale : une valeur initiale de x pour la recherche de zéro
|
|
/// max_delta_x : si > 0 donne le norme maxi permise sur delta_x au cours d'une itération
|
|
/// si ||delta_x|| est plus grand on fait delta_x = max_delta_x * delta_x / || delta_x||
|
|
/// ramène en sortie:
|
|
/// un booléen qui indique si la résolution est correcte ou non
|
|
/// racine : la racine trouvée
|
|
/// der_at_racine : contient en retour la valeur de la dérivée pour la racine trouvée
|
|
/// nb_incr_total : le nombre total d'incrément qui a été nécessaire
|
|
/// nb_iter_total : le nombre total d'itération, qui cumule les iter de tous les incréments
|
|
template <class T> bool Newton_raphson
|
|
(T& instance,double (T::*Ptfonc) (double & alpha,double & x,int& test)
|
|
,double (T::*Ptder_fonc) (double & alpha,double & x,int& test)
|
|
,double val_initiale,double & racine,double & der_at_racine
|
|
,int& nb_incr_total, int& nb_iter_total
|
|
,double max_delta_x) const ;
|
|
/// idem précédemment, mais pour une fonction à valeur vectorielle
|
|
/// les matrices sont telles que: der_at_racine(i,j) = df(i)/d(x(j), soit: df = der_at_racine * dx
|
|
template <class T> bool Newton_raphson
|
|
(T& instance,Vecteur& (T::*Ptfonc) (const double & alpha,const Vecteur & x,int& test)
|
|
,Mat_abstraite& (T::*Ptder_fonc) (const double & alpha,const Vecteur & x,int& test)
|
|
,const Vecteur& val_initiale,Vecteur & racine,Mat_abstraite & der_at_racine
|
|
,int& nb_incr_total, int& nb_iter_total
|
|
,double max_delta_x) const;
|
|
/// idem précédemment, mais pour une fonction à valeur vectorielle et ..
|
|
/// la méthode externe calcul la valeur de la fonction et de la dérivée en même temps
|
|
/// mais on utilise également la fonction résidu toute seule
|
|
/// les matrices sont telles que: der_at_racine(i,j) = df(i)/d(x(j), soit: df = der_at_racine * dx
|
|
template <class T> bool Newton_raphson
|
|
(T& instance,Vecteur& (T::*Ptfonc) (const double & alpha,const Vecteur & x,int& test)
|
|
,Mat_abstraite& (T::*Ptder_fonc) (const double & alpha,const Vecteur & x,Vecteur& res,int& test)
|
|
,const Vecteur& val_initiale,Vecteur & racine,Mat_abstraite & der_at_racine
|
|
,int& nb_incr_total, int& nb_iter_total
|
|
,double max_delta_x) const ;
|
|
|
|
|
|
/// méthodes pour modifier les différents paramètres
|
|
/// précision sur le résidu à convergence
|
|
void Modif_prec_res_abs (double eps) {prec_res_abs = eps;}; ///< précision absolue
|
|
void Modif_prec_res_rel (double eps) {prec_res_rel = eps;}; ///< précision relative
|
|
/// nombre d'itération maxi pour un incrément
|
|
void Modif_iter_max (double eps) { iter_max = eps;};
|
|
/// précision relative sur le mini de l'incrément de x
|
|
void Modif_coef_mini_delta_x (double eps) {coef_mini_delta_x = eps;};
|
|
/// minimum de delta x permis
|
|
void Modif_mini_delta_x (double eps) {mini_delta_x = eps;};
|
|
/// maximum de delta x permis
|
|
void Modif_maxi_delta_x (double eps) {maxi_delta_x = eps;};
|
|
/// nombre maxi d'incrément permis
|
|
void Modif_nbMaxiIncre (int nb) {nbMaxiIncre = nb;};
|
|
/// nombre maxi de test -1 permis
|
|
void Modif_nbMaxi_test_moins1 (int nb) {maxi_test_moins1 = nb;};
|
|
/// initialisation de tous les paramètres à leurs valeurs par défaut
|
|
void Init_param_val_defaut();
|
|
/// le niveau d'affichage
|
|
void Modif_affichage(int niveau) {permet_affichage = niveau;};
|
|
|
|
/// récup en lecture des différents paramètres
|
|
const double& Prec_res_abs() const {return prec_res_abs;}; ///< précision absolue sur le résidu à convergence
|
|
const double& Prec_res_rel() const {return prec_res_rel;};///< précision relative sur le résidu à convergence
|
|
const int& Iter_max() const {return iter_max;};///< nombre d'itération maxi pour un incrément
|
|
const double& Coef_mini_delta_x() const {return coef_mini_delta_x;};///< précision relative sur le mini de l'incrément de x
|
|
const double& Mini_delta_x() const {return mini_delta_x;};///< minimum de delta x permis
|
|
const double& Maxi_delta_x() const {return maxi_delta_x;};///< maximum de delta x permis
|
|
const int& NbMaxiIncre() const {return nbMaxiIncre;};///< nombre maxi d'incrément permis
|
|
const int& Maxi_test_moins1() const {return maxi_test_moins1;};///< nombre maxi de test -1 permis
|
|
const int& Permet_affichage() const {return permet_affichage;}; ///< le niveau d'affichage
|
|
|
|
|
|
/// affichage à l'écran des infos
|
|
void Affiche() const;
|
|
|
|
///----- 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 Lecture_base_info(ifstream& ent,const int cas);
|
|
/// cas donne le niveau de sauvegarde
|
|
/// = 1 : on sauvegarde tout
|
|
/// = 2 : on sauvegarde uniquement les données variables (supposées comme telles)
|
|
void Ecriture_base_info(ofstream& sort,const int cas);
|
|
|
|
|
|
private :
|
|
// VARIABLES PROTEGEES :
|
|
double prec_res_abs; // précision absolue sur le résidu à convergence
|
|
double prec_res_rel; // précision relative sur le résidu à convergence
|
|
int iter_max; // nombre d'itération maxi pour un incrément
|
|
double coef_mini_delta_x; // précision relative sur le mini de l'incrément de x
|
|
double mini_delta_x; // minimum de delta x permis
|
|
double maxi_delta_x; // maximum de delta x permis
|
|
int nbMaxiIncre; // nombre maxi d'incrément permis
|
|
int maxi_test_moins1; // nombe maxi de test = -1 permi
|
|
// ----- controle de la sortie des informations
|
|
int permet_affichage; // pour permettre un affichage spécifique dans les méthodes,
|
|
// pour les erreurs et des wernings
|
|
|
|
Quartic algo_quartic; // une suite d'algo pour résoudre une quartic et des cubic
|
|
|
|
|
|
// METHODES PROTEGEES :
|
|
|
|
};
|
|
|
|
// pour faire de l'inline: nécessaire avec les templates
|
|
// on n'inclut que les méthodes templates
|
|
#include "Algo_zero_2.cc"
|
|
#define Algo_zero_deja_inclus
|
|
/// @} // end of group
|
|
|
|
#endif
|