src/INTERPOLATION/MEDMEM_Mapping.hxx

   1 #ifndef MEDMEM_MAPPING_HXX
   2 #define MEDMEM_MAPPING_HXX
   3
   4 #include "MEDMEM_MappingTools.hxx"
   5 #include "MEDMEM_dTree.hxx"
   6
   7 #define NBR_MAX_MAILLES_EXAMINEES 100
   8
   9 #ifndef  NBR_FACES_MAX
  10 #define NBR_FACES_MAX 6
  11 #endif
  12
  13 #define _TEMPLATE_ template <class MAILLAGE, class NUAGEMAILLE, class NUAGENOEUD, class NOEUD, int DIMENSION>
  14 #define _MAPPING_ Mapping<MAILLAGE,NUAGEMAILLE,NUAGENOEUD,NOEUD,DIMENSION>
  15
  16
  17 //////////////////////////////////////////////////////////////////
  18 ///                                                            ///
  19 ///                        DECLARATIONS                        ///
  20 ///                                                            ///
  21 //////////////////////////////////////////////////////////////////
  22
  23 /*********************************************************/
  24 /*                                                       */
  25 /*                   Classe  Mapping                     */
  26 /*                                                       */
  27 /*********************************************************/
  28
  29 // ATTENTION LE NUAGE DE NOEUD EST SUPPOSE NON REDONDANT ET AUCUNE VERIFICATION N'EST FAITE !
  30
  31 template <class MAILLAGE, class NUAGEMAILLE, class NUAGENOEUD, class NOEUD, int DIMENSION> class Mapping
  32 {
  33 protected :
  34         MAILLAGE * maillage_back;
  35         NUAGEMAILLE * mailles_back;
  36         NUAGENOEUD * noeuds_back;
  37         NUAGENOEUD * noeuds_front;
  38         Coordonnees_Barycentriques<NUAGEMAILLE,NUAGENOEUD,NOEUD,DIMENSION> * CB;
  39         dTree<NOEUD,NUAGENOEUD,DIMENSION> * my_dTree;
  40         vector<int> resultat_mapping;
  41         vector<int> point_le_plus_proche;
  42
  43 public :
  44
  45         Mapping():maillage_back(NULL),mailles_back(NULL),noeuds_back(NULL),noeuds_front(NULL),CB(NULL),my_dTree(NULL) {}
  46         Mapping(MAILLAGE * mb,NUAGENOEUD * nb,NUAGENOEUD * nf); // le dTree est crée à l'initialisation, par contre, le mapping lui meme doit etre invoqué
  47         ~Mapping() {if (CB) delete CB;if (my_dTree) delete my_dTree;}
  48         dTree<NOEUD,NUAGENOEUD,DIMENSION> * Donne_dTree() {return my_dTree;}
  49         int Donne_Directions(int num_maille,const NOEUD &n,int etat_face[NBR_FACES_MAX]);
  50         // Méthode interne de localisation
  51         int Trouve_Maille_Contenant_Point_Mth_Co(const NOEUD &n,int num_maille,int num_maille_interdit,int max_loop,int &nbr_mailles_examinees,int flag_convexe);
  52         void Cree_Mapping(int flag_convexe=0);                                             // SUPPOSE NON CONVEXE PAR DEFAUT
  53         void Cree_Mapping(NUAGENOEUD * nf, int flag_convexe=0);                            // SUPPOSE NON CONVEXE PAR DEFAUT
  54         inline int operator[](int i) const {return resultat_mapping[i];}                   // Renvoie la valeur mappé, si le mapping a été fait, sinon, n'importe quoi
  55         inline vector<int> & Get_Mapping() {return resultat_mapping;}                        // Renvoie le vector contenant le mapping
  56         inline int Get_Noeud_Le_Plus_Proche(int i) const {return point_le_plus_proche[i];} // Invoque la méthode de d-Tree qui donne le noeud le plus proche
  57         inline int Exist_dTree() const {return (my_dTree);}                                // Teste si le dTree existe
  58         void affiche()
  59                 {
  60                 for (int i=0;i<resultat_mapping.size();i++) cout<<"Noeud "<<i<<" de noeud plus proche "<<point_le_plus_proche[i]<<" mappé dans maille "<<resultat_mapping[i]<<endl;
  61                 }
  62 };
  63
  64 //////////////////////////////////////////////////////////////////
  65 ///                                                            ///
  66 ///                            CODE                            ///
  67 ///                                                            ///
  68 //////////////////////////////////////////////////////////////////
  69
  70 _TEMPLATE_ void _MAPPING_::Cree_Mapping(NUAGENOEUD * nf, int flag_convexe)
  71         {
  72         noeuds_front=nf;
  73         Cree_Mapping(flag_convexe);
  74         }
  75
  76 _TEMPLATE_ void _MAPPING_::Cree_Mapping(int flag_convexe)
  77         {
  78
  79         if (resultat_mapping.size()==0)
  80                 {
  81                 if (!noeuds_front)
  82                          {
  83                          cerr<<"Mapping : Pas de noeuds à mapper !"<<endl;
  84                          exit(-1);
  85                          }
  86
  87                 int i;
  88                 int nbr_noeuds=noeuds_front->SIZE();
  89                 int num_maille_depart;
  90                 int nma=0;
  91                 resultat_mapping     = vector<int>(nbr_noeuds,UNDEFINED);
  92                 point_le_plus_proche = vector<int>(nbr_noeuds,UNDEFINED);
  93
  94                 // noeuds_back->affiche();
  95
  96                 for (i=0;i<nbr_noeuds;i++)
  97                         {
  98                         point_le_plus_proche[i]=my_dTree->trouve_plus_proche_point((*noeuds_front)[i]);
  99                         num_maille_depart=maillage_back->DONNE_PREMIERE_MAILLE_CONTENANT(point_le_plus_proche[i]);
 100                         resultat_mapping[i]=Trouve_Maille_Contenant_Point_Mth_Co((*noeuds_front)[i],num_maille_depart,num_maille_depart,NBR_MAX_MAILLES_EXAMINEES,nma,flag_convexe);
 101                         }
 102                 }
 103
 104         else
 105                 {
 106                 cout<<"Le mapping semble déja existé, interrogation sur l'existant"<<endl;
 107                 }
 108
 109         }
 110 _TEMPLATE_ _MAPPING_::Mapping(MAILLAGE * mb,NUAGENOEUD * nb,NUAGENOEUD * nf):maillage_back(mb),noeuds_back(nb),noeuds_front(nf),my_dTree(NULL)
 111         {
 112
 113         if (!maillage_back)
 114                 {
 115                 cerr<<"Mapping : constructeur appelé avec Maillage Vide"<<endl;
 116                 exit(-1);
 117                 }
 118
 119         if (!noeuds_back)
 120                 {
 121                 cerr<<"Mapping : constructeur appelé avec Nuage Noeuds Vide"<<endl;
 122                 exit(-1);
 123                 }
 124
 125         mailles_back=maillage_back->DONNE_POINTEUR_NUAGEMAILLE();
 126
 127         CB=new Coordonnees_Barycentriques<NUAGEMAILLE,NUAGENOEUD,NOEUD,DIMENSION>(mailles_back,noeuds_back);
 128
 129         // TEST REDONDANCE
 130         /*
 131         int nnb=noeuds_back->SIZE();
 132         if (nnb<20000)
 133                 {
 134                 cout<<"MAPPING : VERIFICATION REDONDANCES DANS NUAGE NOEUD BACK"<<endl;
 135                 noeuds_back->affiche();
 136                 int i,j;
 137                 vector<int> redondance(nnb,0);
 138                 for (i=0;i<nnb;i++)
 139                         {
 140                         for (j=i+1;j<nnb;j++) if ((*noeuds_back)[i]==(*noeuds_back)[j])
 141                                 {
 142                                 cerr<<"Le Noeud "<<j<<" est identique au Noeud "<<i<<endl;
 143                                 exit(-1);
 144                                 }
 145                         }
 146                 cout<<"FIN TEST VERIFICATION REDONDANCES"<<endl;
 147                 }
 148         // FIN TEST */
 149
 150         my_dTree=new dTree<NOEUD,NUAGENOEUD,DIMENSION>(noeuds_back);
 151
 152         }
 153 // Renvoie :
 154 //     1 si le point est intérieur
 155 //    -1 si le point est extérieur à la maille via uniquement des faces qui ne sont pas au bord
 156 //    -2 si le point est extérieur à la maille par au moins une face de bord
 157 // Et modifie etat_face de telle sorte que :
 158 // etat_face[i] = -1 s'il n'existe pas de voisin via la face i
 159 // etat_face[i] =  0 si le point est intérieur via la face i et que le voisin i existe
 160 // etat_face[i] =  1 si le point est extérieur via la face i et que le voisin i existe
 161 _TEMPLATE_ int _MAPPING_::Donne_Directions(int num_maille,const NOEUD &n,int etat_face[NBR_FACES_MAX])
 162         {
 163         vector<double> ef=CB->Donne_Pseudo_Coord_Baryc(num_maille,n);
 164         int etat_int=VRAI;
 165         int etat_ext_bord=FAUX;
 166         int tf,tv,tb;
 167         int nbr_faces=(*mailles_back)[num_maille].DONNE_NBR_FACES();
 168         for (int i=0;i<nbr_faces;i++)
 169                 {
 170                 tf=(ef[i]<0);
 171                 tv=(maillage_back->DONNE_VOISIN_DE_MAILLE(num_maille,i)==UNDEFINED);
 172                 tb=(maillage_back->EST_AU_BORD_FACE_DE_MAILLE(num_maille,i));
 173                 if (tf)
 174                         {
 175                         etat_int=FAUX;
 176                         if (tb) etat_ext_bord=VRAI;
 177                         }
 178                 if (tv) etat_face[i]=-1;
 179                 else
 180                         {
 181                         if (tf) etat_face[i]=1;
 182                         else etat_face[i]=0;
 183                         }
 184                 }
 185         if (etat_int) return 1;
 186         if (etat_ext_bord) return -2;
 187         return -1;
 188         }
 189 _TEMPLATE_ int _MAPPING_::Trouve_Maille_Contenant_Point_Mth_Co(const NOEUD &n,int num_maille,int num_maille_interdit,int max_loop,int &nbr_mailles_examinees,int flag_convexe)
 190         {
 191
 192         int etat_face[NBR_FACES_MAX];
 193         int i,tmp,nbr_rnd;
 194         int indirection[NBR_FACES_MAX];
 195         int ind_reel;
 196         int num_reel;
 197         int new_num=UNDEFINED;
 198
 199         int test=Donne_Directions(num_maille,n,etat_face);
 200
 201         int nbr_faces=maillage_back->DONNE_NBR_FACES_MAILLE(num_maille);
 202
 203         for (i=0;i<nbr_faces;i++) indirection[i]=i;
 204
 205         nbr_mailles_examinees=0;
 206
 207         while (nbr_mailles_examinees<max_loop)
 208                 {
 209                 if (test==1)
 210                         {
 211                         return num_maille;
 212                         }
 213                 if ((test==-2)&&(flag_convexe))
 214                         {
 215                         return 2*UNDEFINED;
 216                         }
 217                 nbr_mailles_examinees++;
 218                 for (i=0;i<nbr_faces;i++)
 219                         {
 220                         tmp=indirection[i];
 221                         nbr_rnd=rand()%nbr_faces;
 222                         indirection[i]=indirection[nbr_rnd];
 223                         indirection[nbr_rnd]=tmp;
 224                         }
 225                 for (i=0;(i<nbr_faces)&&(new_num==UNDEFINED);i++)
 226                         {
 227                         ind_reel=indirection[i];
 228                         num_reel=maillage_back->DONNE_VOISIN_DE_MAILLE(num_maille,ind_reel);
 229                         if ((etat_face[ind_reel]==1)&&(num_reel!=num_maille_interdit))
 230                                 {
 231                                 new_num=num_reel;
 232                                 }
 233                         }
 234                 for (i=0;(i<nbr_faces)&&(new_num==UNDEFINED);i++)
 235                         {
 236                         ind_reel=indirection[i];
 237                         num_reel=maillage_back->DONNE_VOISIN_DE_MAILLE(num_maille,ind_reel);
 238                         if ((etat_face[ind_reel]==0)&&(num_reel!=num_maille_interdit))
 239                                 {
 240                                 new_num=num_reel;
 241                                 }
 242                         }
 243                 if (new_num==UNDEFINED)
 244                         {
 245                         new_num=num_maille_interdit;
 246                         }
 247                 num_maille_interdit=num_maille;
 248                 num_maille=new_num;
 249                 new_num=UNDEFINED;
 250                 test=Donne_Directions(num_maille,n,etat_face);
 251                 }
 252         return UNDEFINED;
 253         }
 254
 255 #undef _TEMPLATE_
 256 #undef _MAPPING_
 257
 258 #endif