src/INTERPOLATION/MEDMEM_MappingTools.hxx

   1 // Copyright (C) 2005  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   2 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   3 //
   4 // This library is free software; you can redistribute it and/or
   5 // modify it under the terms of the GNU Lesser General Public
   6 // License as published by the Free Software Foundation; either
   7 // version 2.1 of the License.
   8 //
   9 // This library is distributed in the hope that it will be useful
  10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  12 // Lesser General Public License for more details.
  13 //
  14 // You should have received a copy of the GNU Lesser General Public
  15 // License along with this library; if not, write to the Free Software
  16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  17 //
  18 // See http://www.salome-platform.org/
  19 //
  20 #ifndef COORDONNEES_BARYCENTRIQUES_HPP
  21 #define COORDONNEES_BARYCENTRIQUES_HPP
  22
  23 #define _TEMPLATE_SPE_ template <class NUAGEMAILLE, class NUAGENOEUD, class NOEUD>
  24 #define _COORDBARYC_ Coordonnees_Barycentriques<NUAGEMAILLE,NUAGENOEUD,NOEUD,DIMENSION>
  25 #define _COORDBARY_2D_ Coordonnees_Barycentriques<NUAGEMAILLE,NUAGENOEUD,NOEUD,2>
  26 #define _COORDBARY_3D_ Coordonnees_Barycentriques<NUAGEMAILLE,NUAGENOEUD,NOEUD,3>
  27
  28 //////////////////////////////////////////////////////////////////
  29 ///                                                            ///
  30 ///                        DECLARATIONS                        ///
  31 ///                                                            ///
  32 //////////////////////////////////////////////////////////////////
  33
  34 /*********************************************************/
  35 /*                                                       */
  36 /*           Classe Coordonnees_Barycentriques           */
  37 /*                                                       */
  38 /*********************************************************/
  39
  40 // C'est la définition de la classe générique template qui n'est utilisée qu'à travers ses spécialisations
  41 // vu que le nombre de spécialisations a faire est petit (nombre de dimensions utilisée, 3 pour MED)
  42 // et vu que l'acces a ces classes doit etre rapide, car ce sont des classes de calcul
  43 // la technique de spécialisation, plus lourde, mais plus rapide, a été préférée aux techniques d'héritage
  44
  45 template <class NUAGEMAILLE, class NUAGENOEUD, class NOEUD, int DIMENSION> class Coordonnees_Barycentriques
  46 {
  47 // TEMPLATE GENERIQUE VIDE OBLIGE DE PASSER PAR UNE SPECIALISATION
  48 };
  49
  50 /*********************************************************/
  51 /*                                                       */
  52 /*                   Spécialisation 2D                   */
  53 /*                                                       */
  54 /*********************************************************/
  55
  56 _TEMPLATE_SPE_ class _COORDBARY_2D_
  57 {
  58 protected :
  59         NUAGEMAILLE * mailles;
  60         NUAGENOEUD * sommets;
  61
  62         vector<int> etat_coord_baryc;
  63         vector< vector< vector<double> > > coord_baryc;
  64
  65 public :
  66
  67         Coordonnees_Barycentriques():mailles(NULL),sommets(NULL) {}
  68         Coordonnees_Barycentriques(NUAGEMAILLE * m, NUAGENOEUD *n);
  69         ~Coordonnees_Barycentriques() {}
  70         // donne les pseudos coordonnées barycentriques de M dans ma maille de numéro global num_maille dans mailles
  71         // la pseudo coordonnées barycentrique par rapport a une face est la distance normalisée a cette face,
  72         // dans le cas d'une face triangulaire, c'est la coordonnées ba
  73         vector<double> Donne_Pseudo_Coord_Baryc(int num_maille,const NOEUD &M);
  74         vector<double> Calcule_Base_Coord_Baryc(const vector<int> &simplexe_base);
  75         vector<double> Calcule_Coord_Baryc(int num_maille, const NOEUD & M);
  76 };
  77
  78 /*********************************************************/
  79 /*                                                       */
  80 /*                   Spécialisation 3D                   */
  81 /*                                                       */
  82 /*********************************************************/
  83
  84
  85 _TEMPLATE_SPE_ class _COORDBARY_3D_
  86 {
  87 protected :
  88         NUAGEMAILLE * mailles;
  89         NUAGENOEUD * sommets;
  90
  91         vector<int> etat_coord_baryc;
  92         vector< vector< vector<double> > > coord_baryc;
  93
  94 public :
  95
  96         Coordonnees_Barycentriques():mailles(NULL),sommets(NULL) {}
  97         Coordonnees_Barycentriques(NUAGEMAILLE * m, NUAGENOEUD *n);
  98         ~Coordonnees_Barycentriques() {}
  99         vector<double> Donne_Pseudo_Coord_Baryc(int num_maille,const NOEUD &M);
 100         vector<double> Calcule_Base_Coord_Baryc(const vector<int> &simplexe_base);
 101         vector<double> Calcule_Coord_Baryc(int num_maille, const NOEUD & M);
 102 };
 103
 104 //////////////////////////////////////////////////////////////////
 105 ///                                                            ///
 106 ///                            CODE                            ///
 107 ///                                                            ///
 108 //////////////////////////////////////////////////////////////////
 109
 110 _TEMPLATE_SPE_ _COORDBARY_2D_::Coordonnees_Barycentriques(NUAGEMAILLE * m, NUAGENOEUD *n):mailles(m),sommets(n)
 111                 {
 112                 cout<<"Creation des Coordonnées Barycentriques : "<<flush;
 113                 int nbr_mailles=mailles->SIZE();
 114                 etat_coord_baryc=vector<int>(nbr_mailles,FAUX);
 115                 coord_baryc=vector< vector< vector<double> > >(nbr_mailles);
 116                 cout<<"OK ! "<<endl;
 117                 }
 118
 119 _TEMPLATE_SPE_ vector<double> _COORDBARY_2D_::Donne_Pseudo_Coord_Baryc(int num_maille,const NOEUD &M)
 120         {
 121         int i,nbr_faces;
 122         if (etat_coord_baryc[num_maille]==FAUX)
 123                 {
 124                 nbr_faces=(*mailles)[num_maille].DONNE_NBR_FACES();
 125
 126                 coord_baryc[num_maille]=vector< vector<double> >(nbr_faces);
 127
 128                 type_retour simplexe_base;
 129
 130                 for (i=0;i<nbr_faces;i++)
 131                         {
 132                         (*mailles)[num_maille].DONNE_SIMPLEXE_BASE(i,simplexe_base);
 133                         coord_baryc[num_maille][i]=Calcule_Base_Coord_Baryc(vector<int>(&simplexe_base.quoi[0],&simplexe_base.quoi[simplexe_base.combien]));
 134                         etat_coord_baryc[num_maille]=VRAI;
 135                         }
 136                 }
 137         return Calcule_Coord_Baryc(num_maille,M);
 138         }
 139
 140 _TEMPLATE_SPE_ vector<double> _COORDBARY_2D_::Calcule_Base_Coord_Baryc(const vector<int> &simplexe_base)
 141         {
 142         const vector<int> &ref=simplexe_base;
 143         vector<double> retour(3);
 144
 145         double x0=(*sommets)[ref[0]][0];
 146         double y0=(*sommets)[ref[0]][1];
 147         double x1=(*sommets)[ref[1]][0];
 148         double y1=(*sommets)[ref[1]][1];
 149         double x2=(*sommets)[ref[2]][0];
 150         double y2=(*sommets)[ref[2]][1];
 151
 152         double delta=(x1*y2-x2*y1)+(x2*y0-x0*y2)+(x0*y1-x1*y0);
 153
 154         retour[0]=(y1-y2)/delta;
 155         retour[1]=(x2-x1)/delta;
 156         retour[2]=(x1*y2-x2*y1)/delta;
 157
 158         return retour;
 159         }
 160
 161 _TEMPLATE_SPE_ vector<double> _COORDBARY_2D_::Calcule_Coord_Baryc(int num_maille, const NOEUD & M)
 162         {
 163         int i,j;
 164         // for PAL11458
 165         double zero = 0. ;
 166         //vector<double> coord_baryc_M(3,0);
 167         int nbr_faces=coord_baryc[num_maille].size();
 168         vector<double> coord_baryc_M(nbr_faces,zero);
 169         for (i=0;i</*3*/nbr_faces;i++)
 170                 {
 171                 for (j=0;j<2;j++) coord_baryc_M[i]+=coord_baryc[num_maille][i][j]*M[j];
 172                 coord_baryc_M[i]+=coord_baryc[num_maille][i][2];
 173                 }
 174         return coord_baryc_M;
 175         }
 176
 177 _TEMPLATE_SPE_ _COORDBARY_3D_::Coordonnees_Barycentriques(NUAGEMAILLE * m, NUAGENOEUD *n):mailles(m),sommets(n)
 178                 {
 179                 cout<<"Creation des Coordonnées Barycentriques : "<<flush;
 180                 int nbr_mailles=mailles->SIZE();
 181                 etat_coord_baryc=vector<int>(nbr_mailles,FAUX);
 182                 coord_baryc=vector< vector< vector<double> > >(nbr_mailles);
 183                 cout<<"OK ! "<<endl;
 184                 }
 185
 186 _TEMPLATE_SPE_ vector<double> _COORDBARY_3D_::Donne_Pseudo_Coord_Baryc(int num_maille,const NOEUD &M)
 187         {
 188         int i,nbr_faces;
 189         if (etat_coord_baryc[num_maille]==FAUX)
 190                 {
 191                 nbr_faces=(*mailles)[num_maille].DONNE_NBR_FACES();
 192
 193                 coord_baryc[num_maille]=vector< vector<double> >(nbr_faces);
 194
 195                 type_retour simplexe_base;
 196
 197                 for (i=0;i<nbr_faces;i++)
 198                         {
 199                         (*mailles)[num_maille].DONNE_SIMPLEXE_BASE(i,simplexe_base);
 200                         coord_baryc[num_maille][i]=Calcule_Base_Coord_Baryc(vector<int>(&simplexe_base.quoi[0],&simplexe_base.quoi[simplexe_base.combien]));
 201                         etat_coord_baryc[num_maille]=VRAI;
 202                         }
 203                 }
 204         return Calcule_Coord_Baryc(num_maille,M);
 205         }
 206
 207
 208 _TEMPLATE_SPE_ vector<double> _COORDBARY_3D_::Calcule_Base_Coord_Baryc(const vector<int> &simplexe_base)
 209         {
 210         const vector<int> &ref=simplexe_base;
 211         vector<double> retour(4);
 212
 213         double x0=(*sommets)[ref[0]][0];
 214         double y0=(*sommets)[ref[0]][1];
 215         double z0=(*sommets)[ref[0]][2];
 216         double x1=(*sommets)[ref[1]][0];
 217         double y1=(*sommets)[ref[1]][1];
 218         double z1=(*sommets)[ref[1]][2];
 219         double x2=(*sommets)[ref[2]][0];
 220         double y2=(*sommets)[ref[2]][1];
 221         double z2=(*sommets)[ref[2]][2];
 222         double x3=(*sommets)[ref[3]][0];
 223         double y3=(*sommets)[ref[3]][1];
 224         double z3=(*sommets)[ref[3]][2];
 225
 226         double delta1=((y2-y1)*(z3-z1)-(z2-z1)*(y3-y1));
 227         double delta2=((x3-x1)*(z2-z1)-(x2-x1)*(z3-z1));
 228         double delta3=((x2-x1)*(y3-y1)-(x3-x1)*(y2-y1));
 229
 230         double delta=delta1*(x0-x1)+delta2*(y0-y1)+delta3*(z0-z1);
 231
 232         retour[0]=delta1/delta;
 233         retour[1]=delta2/delta;
 234         retour[2]=delta3/delta;
 235         retour[3]=-(delta1*x1+delta2*y1+delta3*z1)/delta;
 236
 237         return retour;
 238         }
 239
 240 _TEMPLATE_SPE_ vector<double> _COORDBARY_3D_::Calcule_Coord_Baryc(int num_maille, const NOEUD & M)
 241         {
 242         int i,j;
 243         int nbr_faces=coord_baryc[num_maille].size();
 244         vector<double> coord_baryc_M(nbr_faces);
 245         for (i=0;i<nbr_faces;i++)
 246                 {//CCRT Porting : constructor of vector<T>(int,const T&) not supported on CCRT
 247                 coord_baryc_M[i]=0.;
 248                 for (j=0;j<3;j++) coord_baryc_M[i]+=coord_baryc[num_maille][i][j]*M[j];
 249                 coord_baryc_M[i]+=coord_baryc[num_maille][i][3];
 250                 }
 251         return coord_baryc_M;
 252         }
 253
 254 //*/
 255
 256 #undef _TEMPLATE_SPE_
 257 // template <class NUAGEMAILLE, class NUAGENOEUD, class NOEUD, int DIMENSION>
 258 #undef _COORDBARYC_
 259 // Coordonnees_Barycentriques<NUAGEMAILLE,NUAGENOEUD,NOEUD,DIMENSION>
 260 #undef _COORDBARY_2D_
 261 // Coordonnees_Barycentriques<NUAGEMAILLE,NUAGENOEUD,NOEUD,2>
 262 #undef _COORDBARY_3D_
 263 // Coordonnees_Barycentriques<NUAGEMAILLE,NUAGENOEUD,NOEUD,3>
 264
 265 #endif