

// 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:        23/01/97                                            *
*                                                                $     *
 *     AUTEUR:      G RIO   (mailto:gerardrio56@free.fr)                *
*                                                                $     *
*     PROJET:      Herezh++                                            *
*                                                                $     *
************************************************************************
*     BUT:   Definition des classes derivees de dimension 2.           *
*                                                                $     *
*     ''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''     *
************************************************************************/

#ifndef Tenseur2_H
#define Tenseur2_H

#include <iostream> 
#include "Tenseur.h"  
#include "PtTabRel.h"
#include "Tableau2_T.h"


/** @defgroup Les_classes_tenseurs_dim2_ordre2
*
*     BUT:   Definir les tenseurs d'ordre 2 de differentes composantes,
*     spécifiquement à la dimension 2.
*     L'objectif principal est de surcharger les differentes operations
*     classiques.
*
*   concernant le produit contracte un fois, en particulier pour les tenseurs mixtes
*   il y a contraction du 2ieme indice du premier tenseur avec le premier indice du second
*   tenseur : Aij * Bjk = C ik <-> A * B = C
*   le tenseur inverse par rapport au produit contracte est defini de la maniere suivante
*   Inverse(A) * A = Id, ainsi l'inverse d'un tenseur BH est un BH idem pour les HB
*   mais l'inverse d'un BB est un HH, et l'inverse d'un HH est un BB
*
*   NB: lorsque les tenseurs mixtes sont issues de tenseurs HH ou BB symetrique, l'ordre
*   de contraction des indices n'a pas d'importance sur le resultat !!
*
*   le produit contracte de deux tenseurs symetriques quelconques ne donne pas un tenseur
*   symetrique !!, donc par exemple la contraction d'un tenseur HB avec HH n'est pas forcement
*   symetrique. Le resultat est symetrique SEULEMENT lorsque ces operations sont effectues
*   avec le tenseur metrique.
*
* \author    Gérard Rio
* \version   1.0
* \date       23/01/97
* \brief       Définition des classes de dimension 2 de type Tenseur d'ordre 2, en coordonnées avec différentes variances.
*
*/


class Tenseur3HH;
      
/// @addtogroup Les_classes_tenseurs_dim2_ordre2
///  @{

//------------------------------------------------------------------
///  Definition des tenseur  derivees de dimension 2.
///          cas des composantes deux fois contravariantes symetriques
//---------------------------------------------------------------------   
/// \author    Gérard Rio
/// \version   1.0
/// \date       23/01/97
class Tenseur2HH : public TenseurHH
{
 // surcharge de l'operator de lecture
  friend istream & operator >> (istream &, Tenseur2HH &);
 // surcharge de l'operator d'ecriture
  friend ostream & operator << (ostream &, const Tenseur2HH &);
  friend class Tenseur3HH; // pour le passage 3D 2D: méthode Affectation_2D_a_3D(..
  public :
    // Constructeur
    Tenseur2HH() ; // par défaut
    // initialisation de toutes les composantes a une meme valeur val 
    Tenseur2HH(const double val);  
    // initialisation avec 3 valeurs différentes correspodantes a
    // (1,1) (2,2) (1,2)  qui est  identique a (2,1)
    Tenseur2HH
      (const double val1,const  double val2,const double val3) ; 
    // DESTRUCTEUR :
    ~Tenseur2HH() ;
    // constructeur a partir d'une instance non differenciee  
    Tenseur2HH (const  TenseurHH &); 
    // constructeur de copie  
    Tenseur2HH (const  Tenseur2HH &); 
    
    // METHODES PUBLIQUES :
    // initialise toutes les composantes à val
    void Inita(double val);   
    // operations 
     TenseurHH & operator + (const TenseurHH &) const ;
     void operator += (const TenseurHH &);
     TenseurHH & operator - () const ;  // oppose du tenseur
     TenseurHH & operator - (const TenseurHH &) const ;
     void  operator -= (const TenseurHH &);
     TenseurHH & operator = (const TenseurHH &);
     TenseurHH & operator = ( const Tenseur2HH & B) 
      { return this->operator=(*((const TenseurHH*) & B)) ;};
     TenseurHH & operator * (const double &) const ;
     void  operator *= (const double &);
     TenseurHH & operator / (const double &) const ;
     void  operator /= (const double &);
	    // affectation de B dans this, les données en trop sont ignorées
	    void Affectation_3D_a_2D(const Tenseur3HH & B);
 
     // affectation de B dans this,  plusZero = false: les données manquantes sont inchangées,
     // plusZero = true: les données manquantes sont mises à 0
     // si au contraire la dimension de B est plus grande que *this, il y a uniquement affectation
     // des données possibles
     void Affectation_trans_dimension(const TenseurHH & B,bool plusZero);
 
     // produit contracte avec un vecteur
     CoordonneeH operator * (const CoordonneeB & ) const ;
       
     // produit contracte contracté une fois A(i,j)*B(j,k)=A.B
     // -> donc c'est l'indice du milieu qui est contracté
     TenseurHH & operator * (const TenseurBH &) const ;
     TenseurHB & operator * (const TenseurBB &) const ;
     // produit contracte contracté deux fois A(i,j)*B(j,i)=A..B
     // -> on contracte d'abord l'indice du milieu puis l'indice externe
     double  operator && (const TenseurBB &) const ;
     
     // test
     int operator == (const TenseurHH &) const ;
     int operator != (const TenseurHH &) const ;
    
     double  Det() const ;     // determinant de la matrice des coordonnees
                    
    // ATTENTION creation d'un tenseur transpose qui est supprime par Libere
     TenseurHH & Transpose() const ; 
    
    // ---- manipulation d'indice ---- -> création de nouveaux tenseurs
    virtual TenseurBB& Baisse2Indices() const;
    virtual TenseurBH& BaissePremierIndice() const;
    virtual TenseurHB& BaisseDernierIndice() const;
         
    // calcul du maximum en valeur absolu des composantes du tenseur
     double MaxiComposante() const;

     // calcul du tenseur inverse par rapport au produit contracte
     TenseurBB & Inverse() const ;
          
    // retourne la composante i,j accessible en lecture/ecriture
	   double& Coor(const int i,const int j);

    // retourne la composante i,j accessible uniquement en lecture
	   double operator () (const int i,const int j) const;
   
     //fonctions static définissant le produit tensoriel de deux vecteurs
     // si les vecteurs sont égaux le tenseur est symétrique sinon il est non symétrique
     static TenseurHH & Prod_tensoriel(const CoordonneeH & aH, const CoordonneeH & bH);
	 
	   // lecture et écriture de données
    istream & Lecture(istream & entree);
    ostream & Ecriture(ostream & sort) const ;
	 	 
    // transformation 2 indices -> 1 indice (c-a-d vecteur unicolonne)
    static int OdVect(const int i, const int j)  {return cdex.odVect(i,j);};
    // transformation 1 indice -> 2 indices
    // en retour 2 indices : (k)(1) et (k)(2) -> les indices (i,j) 
    static int idx_i(const int k)  {return cdex.idx_i(k);};
    static int idx_j(const int k)  {return cdex.idx_j(k);};
	 	 	 
	  protected :
	    // allocator dans la liste de data
	  listdouble3Iter ipointe;  
	  // --- gestion   d'index ----
	  class ChangementIndex
	   { public:
	     ChangementIndex(); 
	     // passage pour les index de la forme vecteur à la forme i,j
	     Tableau <int> idx_i,idx_j; 
	     // passage pour les index de la forme i,j à la forme vecteur
	     Tableau2 <int> odVect;  
	    };
	  static const ChangementIndex  cdex;
 }; 
/// @}  // end of group

class Tenseur_ns3HH;
      
/// @addtogroup Les_classes_tenseurs_dim2_ordre2
///  @{

//------------------------------------------------------------------
///  Definition des tenseur  derivees de dimension 2.
///          cas des composantes deux fois contravariantes non symetriques
///          pour les differencier la dimension = -2
//-----------------------------------------------------------------------   
/// \author    Gérard Rio
/// \version   1.0
/// \date       23/01/97
class Tenseur_ns2HH : public TenseurHH
{
 // surcharge de l'operator de lecture
  friend istream & operator >> (istream &, Tenseur_ns2HH &);
 // surcharge de l'operator d'ecriture
  friend ostream & operator << (ostream &, const Tenseur_ns2HH &);
  friend class Tenseur_ns3HH; // pour le passage 3D 2D: méthode Affectation_3D_a_2D(..
  friend class Tenseur_ns2BB; 
  public :
    // Constructeur
    Tenseur_ns2HH() ; // par défaut
    // initialisation de toutes les composantes a une meme valeur val 
    Tenseur_ns2HH(const double val);  
    // initialisation avec 4 valeurs différentes correspondantes a
    // (1,1) (2,2) (2,1) (1,2), car le tableau n'est pas symetrique !!!!
    //  le premier indice = ligne, le second = colonne
    Tenseur_ns2HH
      (const double val1, const double val2, const double val3,const  double val4) ; 
    // DESTRUCTEUR :
    ~Tenseur_ns2HH() ;
    // constructeur a partir d'une instance non differenciee  
    Tenseur_ns2HH (const  TenseurHH &); 
    // constructeur de copie  
    Tenseur_ns2HH (const  Tenseur_ns2HH &); 
    
    // METHODES PUBLIQUES :
    // initialise toutes les composantes à val
    void Inita(double val);   
    // operations 
     TenseurHH & operator + (const TenseurHH &) const ;
     void operator += (const TenseurHH &);
     TenseurHH & operator - () const ;  // oppose du tenseur
     TenseurHH & operator - (const TenseurHH &) const ;
     void  operator -= (const TenseurHH &);
     TenseurHH & operator = ( const TenseurHH &);
     TenseurHH & operator = ( const Tenseur_ns2HH & B) 
      { return this->operator=(*((const TenseurHH*) & B)) ;};
     TenseurHH & operator * ( const double &) const ;
     void  operator *= ( const double &);
     TenseurHH & operator / ( const double &) const ;
     void  operator /= ( const double &);
     // affectation de B dans this, les données en trop sont ignorées
     void Affectation_3D_a_2D(const Tenseur_ns3HH & B);  

     // affectation de B dans this,  plusZero = false: les données manquantes sont inchangées,
     // plusZero = true: les données manquantes sont mises à 0
     // si au contraire la dimension de B est plus grande que *this, il y a uniquement affectation
     // des données possibles
     void Affectation_trans_dimension(const TenseurHH & B,bool plusZero);
 
     // produit contracte avec un vecteur
     CoordonneeH operator * ( const CoordonneeB & ) const ;
       
     // produit contracte contracté une fois A(i,j)*B(j,k)=A.B
     // -> donc c'est l'indice du milieu qui est contracté
     TenseurHH & operator * ( const TenseurBH &) const ;
     TenseurHB & operator * ( const TenseurBB &) const ;
     // produit contracte contracté deux fois A(i,j)*B(j,i)=A..B
     // -> on contracte d'abord l'indice du milieu puis l'indice externe
     double  operator && ( const TenseurBB &) const ;
     
     // test
     int operator == ( const TenseurHH &) const ;
     int operator != ( const TenseurHH &) const ;
     double  Det() const ;     // determinant de la matrice des coordonnees
     // ATTENTION creation d'un tenseur transpose qui est supprime par Libere
     TenseurHH & Transpose() const ; 
     
     // ---- manipulation d'indice ---- -> création de nouveaux tenseurs
     virtual TenseurBB& Baisse2Indices() const;
     virtual TenseurBH& BaissePremierIndice() const;
     virtual TenseurHB& BaisseDernierIndice() const;
     
     // calcul du maximum en valeur absolu des composantes du tenseur
     double MaxiComposante() const;

     // calcul du tenseur inverse par rapport au produit contracte
     TenseurBB & Inverse() const ;
          
    // retourne la composante i,j en lecture/écriture
    double& Coor( const int i, const int j);
    
    // retourne la composante i,j en lecture seule
    double operator () ( const int i, const int j) const;
	 
	// lecture et écriture de données
    istream & Lecture(istream & entree);
    ostream & Ecriture(ostream & sort) const ;
	 	 
    // transformation 2 indices -> 1 indice (c-a-d vecteur unicolonne)
    static int OdVect(const int i, const int j)  {return cdex.odVect(i,j);};
    // transformation 1 indice -> 2 indices
    // en retour 2 indices : (k)(1) et (k)(2) -> les indices (i,j) 
    static int idx_i(const int k)  {return cdex.idx_i(k);};
    static int idx_j(const int k)  {return cdex.idx_j(k);};
	 	 	 
	  protected :
	    // allocator dans la liste de data
	  listdouble4Iter ipointe;  
	  // --- gestion   d'index ----
	  class ChangementIndex
	   { public:
	     ChangementIndex(); 
	     // passage pour les index de la forme vecteur à la forme i,j
	     Tableau <int> idx_i,idx_j; 
	     // passage pour les index de la forme i,j à la forme vecteur
	     Tableau2 <int> odVect;  
	    };
	  static const ChangementIndex  cdex;
};  
/// @}  // end of group

class Tenseur3BB;
      
/// @addtogroup Les_classes_tenseurs_dim2_ordre2
///  @{

//------------------------------------------------------------------
///  Definition des tenseur  derivees de dimension 2.
///          cas des composantes deux fois covariantes symetriques
//------------------------------------------------------------------
/// \author    Gérard Rio
/// \version   1.0
/// \date       23/01/97
class Tenseur2BB : public TenseurBB
{
 // surcharge de l'operator de lecture
  friend istream & operator >> (istream &, Tenseur2BB &);
 // surcharge de l'operator d'ecriture
  friend ostream & operator << (ostream &, const Tenseur2BB &);
  friend class Tenseur3BB; // pour le passage 3D 2D: méthode Affectation_3D_a_2D(..
  public :
     // constructeur
    Tenseur2BB() ; // par défaut
        // initialisation de toutes les composantes a une meme valeur val 
    Tenseur2BB( const double val);  
    // initialisation avec 3 valeurs différentes correspodantes a
    // (1,1) (2,2) (1,2)  qui est  identique a (2,1)
    Tenseur2BB
      ( const double val1, const  double val2, const  double val3) ;      
    // DESTRUCTEUR :
    ~Tenseur2BB() ;
    // constructeur a partir d'une instance non differenciee  
    Tenseur2BB ( const  TenseurBB &); 
    // constructeur de copie  
    Tenseur2BB (  const Tenseur2BB &); 
    
    // METHODES PUBLIQUES :
    // initialise toutes les composantes à val
    void Inita(double val);   
    // operations 
    TenseurBB & operator + ( const TenseurBB &) const ;
    void operator += ( const TenseurBB &);
    TenseurBB & operator - () const ;  // oppose du tenseur
    TenseurBB & operator - (const TenseurBB &) const ;
    void operator -= ( const TenseurBB &);
    TenseurBB & operator = ( const TenseurBB &);
    TenseurBB & operator = ( const Tenseur2BB & B) 
      { return this->operator=(*((const TenseurBB*) & B)) ;};
    TenseurBB & operator * ( const double &) const ;
    void operator *= ( const double &) ;
    TenseurBB & operator / ( const double &) const ;
    void operator /= ( const double &);
    // affectation de B dans this, les données en trop sont ignorées
    void Affectation_3D_a_2D(const Tenseur3BB & B);  
             
    // affectation de B dans this,  plusZero = false: les données manquantes sont inchangées,
    // plusZero = true: les données manquantes sont mises à 0
    // si au contraire la dimension de B est plus grande que *this, il y a uniquement affectation
    // des données possibles
    void Affectation_trans_dimension(const TenseurBB & B,bool plusZero);
 
    // produit contracte avec un vecteur
    CoordonneeB operator * ( const CoordonneeH & ) const ;
    
    // produit contracte contracté une fois A(i,j)*B(j,k)=A.B
    // -> donc c'est l'indice du milieu qui est contracté
    TenseurBB & operator * ( const TenseurHB &) const ;
    TenseurBH & operator * ( const TenseurHH &) const ;
    // produit contracte contracté deux fois A(i,j)*B(j,i)=A..B
    // -> on contracte d'abord l'indice du milieu puis l'indice externe
    double  operator && ( const TenseurHH &) const ;
         
     // test
     int operator == ( const TenseurBB &) const ;
     int operator != ( const TenseurBB &) const ;
    
     double  Det() const ;     // determinant de la matrice des coordonnees
          
     // ATTENTION creation d'un tenseur transpose qui est supprime par Libere
     TenseurBB & Transpose() const ; 
    
     // ---- manipulation d'indice ---- -> création de nouveaux tenseurs
     virtual TenseurHH& Monte2Indices() const ;
     virtual TenseurHB& MontePremierIndice() const ;
     virtual TenseurBH& MonteDernierIndice() const ;

     // calcul du maximum en valeur absolu des composantes du tenseur
     double MaxiComposante() const;

     // calcul du tenseur inverse par rapport au produit contracte
     TenseurHH & Inverse() const ;
                       
     // retourne la composante i,j en lecture/écriture
     double& Coor( const int i, const int j);
    
     // retourne la composante i,j en lecture seule
     double operator () ( const int i, const int j) const;
   
     //fonctions static définissant le produit tensoriel de deux vecteurs
     // si les vecteurs sont égaux le tenseur est symétrique sinon il est non symétrique
     static TenseurBB & Prod_tensoriel(const CoordonneeB & aB, const CoordonneeB & bB);
	 
	// lecture et écriture de données
    istream & Lecture(istream & entree);
    ostream & Ecriture(ostream & sort) const ;
	 	 
    // transformation 2 indices -> 1 indice (c-a-d vecteur unicolonne)
    static int OdVect(const int i, const int j)  {return cdex.odVect(i,j);};
    // transformation 1 indice -> 2 indices
    // en retour 2 indices : (k)(1) et (k)(2) -> les indices (i,j) 
    static int idx_i(const int k)  {return cdex.idx_i(k);};
    static int idx_j(const int k)  {return cdex.idx_j(k);};
	 	 	 
	  protected :
	    // allocator dans la liste de data
	  listdouble3Iter ipointe;  
	  // --- gestion   d'index ----
	  class ChangementIndex
	   { public:
	     ChangementIndex(); 
	     // passage pour les index de la forme vecteur à la forme i,j
	     Tableau <int> idx_i,idx_j; 
	     // passage pour les index de la forme i,j à la forme vecteur
	     Tableau2 <int> odVect;  
	    };
	  static const ChangementIndex  cdex;
};    
/// @}  // end of group

class Tenseur_ns3BB;

/// @addtogroup Les_classes_tenseurs_dim2_ordre2
///  @{

//------------------------------------------------------------------
///  Definition des tenseur  derivees de dimension 2.
///          cas des composantes deux fois covariantes non symetriques
///          pour les differencier la dimension = -2
//--------------------------------------------------------------------
/// \author    Gérard Rio
/// \version   1.0
/// \date       23/01/97
class Tenseur_ns2BB : public TenseurBB
{
 // surcharge de l'operator de lecture
  friend istream & operator >> (istream &, Tenseur_ns2BB &);
 // surcharge de l'operator d'ecriture
  friend ostream & operator << (ostream &, const Tenseur_ns2BB &);
  friend class Tenseur_ns3HH; // pour le passage 3D 2D: méthode Affectation_3D_a_2D(..
  friend class Tenseur_ns2HH; 
  public :
     // constructeur
    Tenseur_ns2BB() ; // par défaut
        // initialisation de toutes les composantes a une meme valeur val 
    Tenseur_ns2BB( const double val);  
    // initialisation avec 4 valeurs différentes correspodantes a
    // (1,1) (2,2) (2,1) (1,2), car le tableau n'est pas symetrique !!!!
    //  le premier indice = ligne, le second = colonne
    Tenseur_ns2BB
      ( const double val1, const  double val2, const  double val3, const  double val4) ; 
    // DESTRUCTEUR :
    ~Tenseur_ns2BB() ;
    // constructeur a partir d'une instance non differenciee  
    Tenseur_ns2BB ( const  TenseurBB &); 
    // constructeur de copie  
    Tenseur_ns2BB ( const  Tenseur_ns2BB &); 
    
    // METHODES PUBLIQUES :
    // initialise toutes les composantes à val
    void Inita(double val);   
    // operations 
    TenseurBB & operator + ( const TenseurBB &) const ;
    void operator += ( const TenseurBB &);
    TenseurBB & operator - () const ;  // oppose du tenseur
    TenseurBB & operator - ( const TenseurBB &) const ;
    void operator -= ( const TenseurBB &);
    TenseurBB & operator = ( const TenseurBB &);
    TenseurBB & operator = ( const Tenseur_ns2BB & B) 
      { return this->operator=(*((const TenseurBB*) & B)) ;};
    TenseurBB & operator * ( const double &) const ;
    void operator *= ( const double &);
    TenseurBB & operator / (const double &) const ;
    void operator /= ( const double &);
    // affectation de B dans this, les données en trop sont ignorées
    void Affectation_3D_a_2D(const Tenseur_ns3BB & B);  
             
    // affectation de B dans this,  plusZero = false: les données manquantes sont inchangées,
    // plusZero = true: les données manquantes sont mises à 0
    // si au contraire la dimension de B est plus grande que *this, il y a uniquement affectation
    // des données possibles
    void Affectation_trans_dimension(const TenseurBB & B,bool plusZero);
 
    // produit contracte avec un vecteur
    CoordonneeB operator * ( const CoordonneeH & ) const ;
    
    // produit contracte contracté une fois A(i,j)*B(j,k)=A.B
    // -> donc c'est l'indice du milieu qui est contracté
    TenseurBB & operator * ( const TenseurHB &) const ;
    TenseurBH & operator * ( const TenseurHH &) const ;
    // produit contracte contracté deux fois A(i,j)*B(j,i)=A..B
    // -> on contracte d'abord l'indice du milieu puis l'indice externe
    double  operator && ( const TenseurHH &) const ;
         
     // test
     int operator == ( const TenseurBB &) const ;
     int operator != ( const TenseurBB &) const ;
    
     double  Det() const ;     // determinant de la matrice des coordonnees
          
     // ATTENTION creation d'un tenseur transpose qui est supprime par Libere
     TenseurBB & Transpose() const ; 
    
     // ---- manipulation d'indice ---- -> création de nouveaux tenseurs
     virtual TenseurHH& Monte2Indices() const ;
     virtual TenseurHB& MontePremierIndice() const ;
     virtual TenseurBH& MonteDernierIndice() const ;

     // calcul du maximum en valeur absolu des composantes du tenseur
     double MaxiComposante() const;

     // calcul du tenseur inverse par rapport au produit contracte
     TenseurHH & Inverse() const ;
                       
     // retourne la composante i,j en lecture/écriture
     double& Coor( const int i, const int j);
    
     // retourne la composante i,j en lecture seule
     double operator () ( const int i, const int j) const;
	 
	 // lecture et écriture de données
    istream & Lecture(istream & entree);
    ostream & Ecriture(ostream & sort) const ;
	 	 
    // transformation 2 indices -> 1 indice (c-a-d vecteur unicolonne)
    static int OdVect(const int i, const int j)  {return cdex.odVect(i,j);};
    // transformation 1 indice -> 2 indices
    // en retour 2 indices : (k)(1) et (k)(2) -> les indices (i,j) 
    static int idx_i(const int k)  {return cdex.idx_i(k);};
    static int idx_j(const int k)  {return cdex.idx_j(k);};
	 	 	 
	  protected :
	    // allocator dans la liste de data
	  listdouble4Iter ipointe;  
	  // --- gestion   d'index ----
	  class ChangementIndex
	   { public:
	     ChangementIndex(); 
	     // passage pour les index de la forme vecteur à la forme i,j
	     Tableau <int> idx_i,idx_j; 
	     // passage pour les index de la forme i,j à la forme vecteur
	     Tableau2 <int> odVect;  
	    };
	  static const ChangementIndex  cdex;
};    
/// @}  // end of group

class Tenseur3BH;
      
/// @addtogroup Les_classes_tenseurs_dim2_ordre2
///  @{

//------------------------------------------------------------------
///  Definition des tenseur  derivees de dimension 2.
///          cas des composantes mixtes BH
//------------------------------------------------------------------
/// \author    Gérard Rio
/// \version   1.0
/// \date       23/01/97
class Tenseur2BH : public TenseurBH
{
 // surcharge de l'operator de lecture
  friend istream & operator >> (istream &, Tenseur2BH &);
 // surcharge de l'operator d'ecriture
  friend ostream & operator << (ostream &, const Tenseur2BH &);
  friend class Tenseur3BH; // pour le passage 3D 2D: méthode Affectation_3D_a_2D(..
  public :
    // constructeur
    Tenseur2BH() ; // par défaut
    // initialisation de toutes les composantes a une meme valeur val 
    Tenseur2BH ( const double val) ;  
    // initialisation avec 4 valeurs différentes correspodantes a
    // (1,1) (2,2) (2,1) (1,2), car le tableau n'est pas symetrique !!!!
    //  le premier indice = ligne, le second = colonne
    Tenseur2BH
      ( const double val1,  const double val2, const  double val3,  const double val4) ; 
  
    // DESTRUCTEUR :
    ~Tenseur2BH();
    // constructeur a partir d'une instance non differenciee  
    Tenseur2BH ( const  TenseurBH &); 
    // constructeur de copie  
    Tenseur2BH (  const Tenseur2BH &); 
    
    // METHODES PUBLIQUES :
    // initialise toutes les composantes à val
    void Inita(double val);   
    // operations 
    TenseurBH & operator + ( const TenseurBH &) const ;
    void operator += ( const TenseurBH &);
    TenseurBH & operator - () const ;  // oppose du tenseur
    TenseurBH & operator - ( const TenseurBH &) const ;
    void operator -= ( const TenseurBH &);
    TenseurBH & operator = ( const TenseurBH &);
    TenseurBH & operator = ( const Tenseur2BH & B) 
      { return this->operator=(*((const TenseurBH*) & B)) ;};
    TenseurBH & operator * ( const double &) const ;
    void operator *= ( const double &);
    TenseurBH & operator / ( const double &) const ;
    void operator /= ( const double &);
    // affectation de B dans this, les données en trop sont ignorées
    void Affectation_3D_a_2D(const Tenseur3BH & B);  
        
    // affectation de B dans this,  plusZero = false: les données manquantes sont inchangées,
    // plusZero = true: les données manquantes sont mises à 0
    // si au contraire la dimension de B est plus grande que *this, il y a uniquement affectation
    // des données possibles
    void Affectation_trans_dimension(const TenseurBH & B,bool plusZero);

    // produit contracte contracté une fois A(i,j)*B(j,k)=A.B
    // -> donc c'est l'indice du milieu qui est contracté
    CoordonneeB operator *  (const CoordonneeB & ) const ;
    TenseurBB & operator * ( const TenseurBB &) const ; // produit une fois contracte
    TenseurBH & operator * ( const TenseurBH &) const ; // produit une fois contracte
    // produit contracte contracté deux fois A(i,j)*B(j,i)=A..B
    // -> on contracte d'abord l'indice du milieu puis l'indice externe
    double  operator && ( const TenseurBH &) const ; // produit deux fois contracte
         
    // test
    int operator == ( const TenseurBH &) const ;
    int operator != ( const TenseurBH &) const ;
    // calcul du tenseur inverse par rapport au produit contracte
    TenseurBH & Inverse() const ;
    
    double  Trace() const ;   // trace du tenseur ou premier invariant
    double  II() const ;      // second invariant = trace (A*A)
    double  III() const ;     // troisieme invariant = trace ((A*A)*A)
    double  Det() const ;     // determinant de la matrice des coordonnees
    // valeurs propre  dans le vecteur  de retour, classée par ordres "décroissants"
    // cas indique le cas de valeur propre:  
    //  cas = -1, si l'extraction des valeurs propres n'a pas pu se faire
    //            dans ce cas les valeurs propres de retour sont nulles par défaut
    // dim = 2 , cas = 1 si deux valeurs propres distinctes, cas = 0 si deux val propres identiques
    virtual Coordonnee ValPropre(int& cas) const ;
    // idem met en retour la matrice mat contiend par colonne les vecteurs propre
    // elle doit avoir la dimension du tenseur
    // le premier vecteur propre est exprime dans le repere dual 
    // le second vecteur propre est exprimé dans le repère naturel 
    virtual Coordonnee ValPropre(int& cas, Mat_pleine& mat) const ;  

    // ici il s'agit uniquement de calculer les vecteurs propres, les valeurs propres
    // étant déjà connues
    // en retour VP les vecteurs propre : doivent avoir la dimension du tenseur
    // les vecteurs propre sont exprime dans le repere naturel (pour les tenseurs dim 3
    // pour dim=2:le premier vecteur propre est exprime dans le repere naturel 
    //            le second vecteur propre est exprimé dans le repère dual
    // pour dim=1 le vecteur est dans la base naturelle (mais ça importe pas)
    //   en sortie cas = -1 s'il y a eu un problème, dans ce cas, V_P est quelconque
    //             sinon si tout est ok, cas est identique en sortie avec l'entrée
    virtual void VecteursPropres(const Coordonnee& Val_P,int& cas, Tableau <Coordonnee>& V_P) const;
 
    // tenseur transpose
    TenseurHB & Transpose() const ; 
    // permute Bas Haut, mais reste dans le même tenseur
    void PermuteHautBas();
 
    // calcul du maximum en valeur absolu des composantes du tenseur
     double MaxiComposante() const;
             
    // retourne la composante i,j en lecture/écriture
    double& Coor( const int i, const int j);
    
    // retourne la composante i,j en lecture seule
    double operator () ( const int i, const int j) const;
   
    //fonctions static définissant le produit tensoriel de deux vecteurs
    static TenseurBH & Prod_tensoriel(const CoordonneeB & aB, const CoordonneeH & bH);
	 
	// lecture et écriture de données
    istream & Lecture(istream & entree);
    ostream & Ecriture(ostream & sort) const ;
	 	 
    // transformation 2 indices -> 1 indice (c-a-d vecteur unicolonne)
    static int OdVect(const int i, const int j)  {return cdex.odVect(i,j);};
    // transformation 1 indice -> 2 indices
    // en retour 2 indices : (k)(1) et (k)(2) -> les indices (i,j) 
    static int idx_i(const int k)  {return cdex.idx_i(k);};
    static int idx_j(const int k)  {return cdex.idx_j(k);};
	 	 	 
	  protected :
	    // allocator dans la liste de data
	  listdouble4Iter ipointe;  
    
	  // --- gestion   d'index ----
	  class ChangementIndex
	   { public:
	     ChangementIndex(); 
	     // passage pour les index de la forme vecteur à la forme i,j
	     Tableau <int> idx_i,idx_j; 
	     // passage pour les index de la forme i,j à la forme vecteur
	     Tableau2 <int> odVect;  
	    };
	  static const ChangementIndex  cdex;
}; 
/// @}  // end of group

class Tenseur3HB;
      
/// @addtogroup Les_classes_tenseurs_dim2_ordre2
///  @{

//------------------------------------------------------------------
///  Definition des tenseur  derivees de dimension 2.
///          cas des composantes mixtes HB
//------------------------------------------------------------------
/// \author    Gérard Rio
/// \version   1.0
/// \date       23/01/97
class Tenseur2HB : public TenseurHB
{
 // surcharge de l'operator de lecture
  friend istream & operator >> (istream &, Tenseur2HB &);
 // surcharge de l'operator d'ecriture
  friend ostream & operator << (ostream &, const Tenseur2HB &);
  friend class Tenseur3HB; // pour le passage 3D 2D: méthode Affectation_3D_a_2D(..
  public :
      // Constructeur
    Tenseur2HB() ; // par défaut
    // initialisation de toutes les composantes a une meme valeur val 
    Tenseur2HB ( const double val) ;  
    // initialisation avec 4 valeurs différentes correspodantes a
    // (1,1) (2,2) (2,1) (1,2), car le tableau n'est pas symetrique !!!!
    //  le premier indice = ligne, le second = colonne
    Tenseur2HB
      ( const double val1, const  double val2, const  double val3, const  double val4) ; 
    // constructeur a partir d'une instance non differenciee  
    Tenseur2HB ( const  TenseurHB &); 
    // constructeur de copie  
    Tenseur2HB ( const Tenseur2HB &); 

    // DESTRUCTEUR :
    ~Tenseur2HB();
    
    // METHODES PUBLIQUES :
    // initialise toutes les composantes à val
    void Inita(double val);   
    // operations 
    TenseurHB & operator + ( const TenseurHB &) const ;
    void operator += ( const TenseurHB & );
    TenseurHB & operator - () const ;  // oppose du tenseur
    TenseurHB & operator - ( const TenseurHB &) const ;
    void operator -= ( const TenseurHB &);
    TenseurHB & operator = ( const TenseurHB &);
    TenseurHB & operator = ( const Tenseur2HB & B) 
      { return this->operator=(*((const TenseurHB*) & B)) ;};
    TenseurHB & operator * ( const double &) const ;
    void operator *= ( const double &);
    TenseurHB & operator / (const double &) const ;
    void operator /= ( const double &);
    // affectation de B dans this, les données en trop sont ignorées
    void Affectation_3D_a_2D(const Tenseur3HB & B);  
    
    // affectation de B dans this,  plusZero = false: les données manquantes sont inchangées,
    // plusZero = true: les données manquantes sont mises à 0
    // si au contraire la dimension de B est plus grande que *this, il y a uniquement affectation
    // des données possibles
    void Affectation_trans_dimension(const TenseurHB & B,bool plusZero);

    // produit contracte avec un vecteur
    CoordonneeH operator * ( const CoordonneeH & ) const ;
    
    // produit contracte contracté une fois A(i,j)*B(j,k)=A.B
    // -> donc c'est l'indice du milieu qui est contracté
    TenseurHH & operator * ( const TenseurHH &) const ;
    TenseurHB & operator * ( const TenseurHB &) const ; // produit une fois contracte
    // produit contracte contracté deux fois A(i,j)*B(j,i)=A..B
    // -> on contracte d'abord l'indice du milieu puis l'indice externe
    double  operator && ( const TenseurHB &) const ; // produit deux fois contracte
     
    // test
    int operator == ( const TenseurHB &) const ;
    int operator != ( const TenseurHB &) const ;
    // calcul du tenseur inverse par rapport au produit contracte
    TenseurHB & Inverse()  const ;
                  
    double  Trace() const ;   // trace du tenseur ou premier invariant
    double  II() const ;      // second invariant = trace (A*A)
    double  III() const ;     // troisieme invariant = trace ((A*A)*A)
    double  Det() const ;     // determinant de la matrice des coordonnees
    // valeurs propre  dans le vecteur  de retour, classée par ordres "décroissants"
    // cas indique le cas de valeur propre:  
    //  cas = -1, si l'extraction des valeurs propres n'a pas pu se faire
    //            dans ce cas les valeurs propres de retour sont nulles par défaut
    // dim = 2 , cas = 1 si deux valeurs propres distinctes, cas = 0 si deux val propres identiques
    virtual Coordonnee ValPropre(int& cas) const ;
    // idem met en retour la matrice mat contiend par colonne les vecteurs propre
    // elle doit avoir la dimension du tenseur
    // le premier vecteur propre est exprimé dans le repere naturel 
    // le second vecteur propre est exprimé dans le repère dual 
    virtual Coordonnee ValPropre(int& cas, Mat_pleine& mat) const ;  

    // ici il s'agit uniquement de calculer les vecteurs propres, les valeurs propres
    // étant déjà connues
    // en retour VP les vecteurs propre : doivent avoir la dimension du tenseur
    // les vecteurs propre sont exprime dans le repere naturel (pour les tenseurs dim 3
    // pour dim=2:le premier vecteur propre est exprime dans le repere naturel 
    //            le second vecteur propre est exprimé dans le repère dual
    // pour dim=1 le vecteur est dans la base naturelle (mais ça importe pas)
    //   en sortie cas = -1 s'il y a eu un problème, dans ce cas, V_P est quelconque
    //             sinon si tout est ok, cas est identique en sortie avec l'entrée
    virtual void VecteursPropres(const Coordonnee& Val_P,int& cas, Tableau <Coordonnee>& V_P) const;
 
    // tenseur transpose
    TenseurBH & Transpose() const ; 
    // permute Bas Haut, mais reste dans le même tenseur
    void PermuteHautBas();
     
    // calcul du maximum en valeur absolu des composantes du tenseur
     double MaxiComposante() const;
             
    // retourne la composante i,j en lecture/écriture
    double& Coor( const int i, const int j);
    
    // retourne la composante i,j en lecture seule
    double operator () ( const int i, const int j) const;
    
//    operator TenseurHB & (void) 
//      { return *this;};
 
    //fonctions static définissant le produit tensoriel de deux vecteurs
    static TenseurHB & Prod_tensoriel(const CoordonneeH & aH, const CoordonneeB & bB);
	 
	// lecture et écriture de données
    istream & Lecture(istream & entree);
    ostream & Ecriture(ostream & sort) const ;
	 	 
    // transformation 2 indices -> 1 indice (c-a-d vecteur unicolonne)
    static int OdVect(const int i, const int j)  {return cdex.odVect(i,j);};
    // transformation 1 indice -> 2 indices
    // en retour 2 indices : (k)(1) et (k)(2) -> les indices (i,j) 
    static int idx_i(const int k)  {return cdex.idx_i(k);};
    static int idx_j(const int k)  {return cdex.idx_j(k);};
	 	 	 
	  protected :
	    // allocator dans la liste de data
	  listdouble4Iter ipointe;  

	  // --- gestion   d'index ----
	  class ChangementIndex
	   { public:
	     ChangementIndex(); 
	     // passage pour les index de la forme vecteur à la forme i,j
	     Tableau <int> idx_i,idx_j; 
	     // passage pour les index de la forme i,j à la forme vecteur
	     Tableau2 <int> odVect;  
	    };
	  static const ChangementIndex  cdex;
};    
/// @}  // end of group

#include "Tenseur3.h"       
 
#ifndef MISE_AU_POINT
  #include "Tenseur2-1.cc"
  #include "Tenseur2-2.cc"
  #include "Tenseur2_ns.cc"
  #define  Tenseur2_H_deja_inclus
#endif
    
 
#endif  
