src/MEDMEM/MEDMEM_Connectivity.cxx

   1 using namespace std;
   2 #include "MEDMEM_Connectivity.hxx"
   3 #include "MEDMEM_Family.hxx"
   4 #include "MEDMEM_CellModel.hxx"
   5
   6 #include "MEDMEM_SkyLineArray.hxx"
   7 #include "MEDMEM_ModulusArray.hxx"
   8
   9 #include "MEDMEM_STRING.hxx"
  10
  11 //------------------------------------------------------//
  12 CONNECTIVITY::CONNECTIVITY(medEntityMesh Entity/*=MED_CELL*/):
  13 //------------------------------------------------------//
  14   _entity(Entity),
  15   _typeConnectivity(MED_NODAL),
  16   _numberOfTypes(0),
  17   _geometricTypes((medGeometryElement*)NULL),
  18   _type((CELLMODEL*)NULL),
  19   _entityDimension(0),
  20   _count((int*)NULL),
  21   _nodal((MEDSKYLINEARRAY*)NULL),
  22   _descending((MEDSKYLINEARRAY*)NULL),
  23   _reverseNodalConnectivity((MEDSKYLINEARRAY*)NULL),
  24   _reverseDescendingConnectivity((MEDSKYLINEARRAY*)NULL),
  25   _neighbourhood((MEDSKYLINEARRAY*)NULL),
  26   _constituent((CONNECTIVITY*)NULL)
  27 {
  28   MESSAGE("CONNECTIVITY(medEntityMesh Entity=MED_CELL)") ;
  29 }
  30
  31 //-------------------------------------------------------------------------//
  32 CONNECTIVITY::CONNECTIVITY(int numberOfTypes,medEntityMesh Entity/*=MED_CELL*/):
  33 //-------------------------------------------------------------------------//
  34   _entity(Entity),
  35   _typeConnectivity(MED_NODAL),
  36   _numberOfTypes(numberOfTypes),
  37   _entityDimension(0),
  38   _nodal((MEDSKYLINEARRAY*)NULL),
  39   _descending((MEDSKYLINEARRAY*)NULL),
  40   _reverseNodalConnectivity((MEDSKYLINEARRAY*)NULL),
  41   _reverseDescendingConnectivity((MEDSKYLINEARRAY*)NULL),
  42   _neighbourhood((MEDSKYLINEARRAY*)NULL),
  43   _constituent((CONNECTIVITY*)NULL)
  44 {
  45   MESSAGE("CONNECTIVITY(int numberOfTypes,medEntityMesh Entity=MED_CELL)") ;
  46   _geometricTypes = new medGeometryElement[numberOfTypes];
  47   _type = new CELLMODEL[numberOfTypes];
  48   _count = new int[numberOfTypes];
  49 }
  50
  51 //-------------------------------------------------------------------------//
  52 CONNECTIVITY::CONNECTIVITY(CONNECTIVITY & m)
  53 //-------------------------------------------------------------------------//
  54 {
  55   _entity = m._entity;
  56   _typeConnectivity = m._typeConnectivity;
  57   _numberOfTypes = m._numberOfTypes;
  58   if (m._geometricTypes != NULL)
  59     {
  60       _geometricTypes = new medGeometryElement[m._numberOfTypes];
  61       memcpy(_geometricTypes,m._geometricTypes,m._numberOfTypes*sizeof(medGeometryElement));
  62     }
  63   else
  64     _geometricTypes = (medGeometryElement *) NULL;
  65   if (m._type != NULL)
  66     _type = new CELLMODEL(* m._type);
  67   else
  68     _type = (CELLMODEL *) NULL;
  69   _entityDimension = m._entityDimension;
  70   _numberOfNodes = m._numberOfNodes;
  71   if (m._count != NULL)
  72     {
  73       _count = new med_int[m._numberOfTypes+1];
  74       memcpy(_count,m._count,(m._numberOfTypes+1)*sizeof(med_int));
  75     }
  76   else
  77     _count = (med_int *) NULL;
  78   if (m._nodal != NULL)
  79     _nodal = new MEDSKYLINEARRAY(* m._nodal);
  80   else
  81     _nodal = (MEDSKYLINEARRAY *) NULL;
  82   if (m._descending != NULL)
  83     _descending = new MEDSKYLINEARRAY(* m._descending);
  84   else
  85     _descending = (MEDSKYLINEARRAY *) NULL;
  86   if (m._reverseNodalConnectivity != NULL)
  87     _reverseNodalConnectivity = new MEDSKYLINEARRAY(* m._reverseNodalConnectivity);
  88   else
  89     _reverseNodalConnectivity = (MEDSKYLINEARRAY *) NULL;
  90   if (m._reverseDescendingConnectivity != NULL)
  91     _reverseDescendingConnectivity = new MEDSKYLINEARRAY(* m._reverseDescendingConnectivity);
  92   else
  93     _reverseDescendingConnectivity = (MEDSKYLINEARRAY *) NULL;
  94   if (m._neighbourhood != NULL)
  95     _neighbourhood = new MEDSKYLINEARRAY(* m._neighbourhood);
  96   else
  97     _neighbourhood = (MEDSKYLINEARRAY *) NULL;
  98   if (m._constituent != NULL)
  99     _constituent = new CONNECTIVITY(* m._constituent);
 100   else
 101     _constituent = (CONNECTIVITY *) NULL;
 102 }
 103
 104 //----------------------------//
 105 CONNECTIVITY::~CONNECTIVITY()
 106 //----------------------------//
 107 {
 108   MESSAGE("Destructeur de CONNECTIVITY()") ;
 109   if ( _geometricTypes != NULL )
 110     delete [] _geometricTypes ;
 111   if ( _count != NULL )
 112     delete[] _count ;
 113   if ( _nodal != NULL )
 114     delete _nodal ;
 115   if ( _descending != NULL )
 116     delete _descending ;
 117   if ( _reverseNodalConnectivity != NULL )
 118     delete _reverseNodalConnectivity ;
 119   if ( _reverseDescendingConnectivity != NULL )
 120     delete _reverseDescendingConnectivity ;
 121   if ( _constituent != NULL )
 122     delete _constituent ;
 123 }
 124
 125 /*! A DOCUMENTER */
 126 //------------------------------------------------------------------------------------------//
 127 void CONNECTIVITY::calculateConnectivity(medConnectivity ConnectivityType, medEntityMesh Entity)
 128 //------------------------------------------------------------------------------------------//
 129 {
 130   MESSAGE("CONNECTIVITY::calculateConnectivity");
 131
 132   // a temporary limitation due to calculteDescendingConnectivity function !!!
 133   if ((_entityDimension==3) & (Entity==MED_EDGE))
 134     throw MEDEXCEPTION("CONNECTIVITY::calculateConnectivity : Could not build EDGE !");
 135
 136   if (Entity==_entity)
 137     if (ConnectivityType==MED_NODAL)
 138       calculateNodalConnectivity() ;
 139     else
 140       if (Entity==MED_CELL)
 141         calculateDescendingConnectivity() ;
 142       else
 143         throw MEDEXCEPTION("CONNECTIVITY::calculateConnectivity : Could not build DESCENDING !");
 144   if (Entity!=_entity) {
 145     calculateDescendingConnectivity() ;
 146     _constituent->calculateConnectivity(ConnectivityType,Entity) ;
 147   }
 148 }
 149
 150 /*!  Give, in full or no interlace mode (for nodal connectivity),
 151      descending or nodal connectivity.
 152
 153       You must give a <medEntityMesh> (ie:MED_EDGE)
 154       and a <medGeometryElement> (ie:MED_SEG3).  */
 155
 156 //------------------------------------------------------------//
 157 void CONNECTIVITY::updateFamily(vector<FAMILY*> myFamilies)
 158 //------------------------------------------------------------//
 159 {
 160   const char * LOC = "CONNECTIVITY::updateFamily(vector<FAMILY*>) ";
 161   BEGIN_OF(LOC);
 162
 163   int numberOfFamilies = myFamilies.size();
 164   if (numberOfFamilies == 0 ) {
 165     MESSAGE(LOC<<"No family") ;
 166     return ;
 167   }
 168   // does we do an update ?
 169   if ((_constituent != NULL)&(_descending != NULL)) {
 170     MESSAGE(LOC<<"Constituent is already defined") ;
 171     return ;
 172   }
 173
 174   if ((_constituent != NULL)&(_descending == NULL)) {
 175     if (myFamilies[0]->getEntity() != _constituent->getEntity()) {
 176       MESSAGE(LOC<<"Family and constituent entity are different. We do nothing") ;
 177       return ;
 178     }
 179
 180     for(int i=0; i<numberOfFamilies; i++) {
 181       FAMILY * myFamily = myFamilies[i] ;
 182       MESSAGE(LOC<<"updating the family (BEGIN) : " << *myFamily);
 183     }
 184
 185     // well we could go !
 186     CONNECTIVITY * oldConstituent = _constituent ;
 187     _constituent = (CONNECTIVITY *)NULL ;
 188     // for instance we must have nodal connectivity in constituent :
 189     if (oldConstituent->_nodal == NULL)
 190       {
 191         MESSAGE(LOC<<"We have no nodal connectivity of sub cell");
 192         throw MED_EXCEPTION(LOCALIZED(STRING(LOC)<<"We have no nodal connectivity of sub cell")) ;
 193       }
 194     int oldNumberOfFace = oldConstituent->_nodal->getNumberOf() ;
 195     int * oldConstituentValue = oldConstituent->_nodal->getValue() ;
 196     int * oldConstituentIndex = oldConstituent->_nodal->getIndex() ;
 197
 198     SCRUTE(oldNumberOfFace);
 199
 200     calculateDescendingConnectivity() ;
 201
 202     //    if (oldConstituent->_nodal != NULL) {
 203     int newNumberOfFace = _constituent->_nodal->getNumberOf() ;
 204     int * newConstituentValue = _constituent->_nodal->getValue() ;
 205     int * newConstituentIndex = _constituent->_nodal->getIndex() ;
 206
 207     SCRUTE(newNumberOfFace);
 208
 209     int * newReverseDescendingIndex =
 210       _reverseDescendingConnectivity->getIndex();
 211     int * newReverseDescendingValue =
 212       _reverseDescendingConnectivity->getValue();
 213
 214     int * newDescendingIndex = _descending->getIndex();
 215     int * newDescendingValue = _descending->getValue();
 216
 217     // loop on all family,
 218     //   for all constituent in family, we get it's old connectivity
 219     //   (with oldCconstituentValue and oldConstituentIndex)
 220     //   and search the constituent in newConstituentValue with class
 221     //   ModulusArry
 222     //
 223     //   when a new face is found, replace old constituent
 224     //   number in family with new one
 225     //   If face is not rigth oriented, we must change _descending attribute
 226     //   and _reverseDescendingConnectivity (see calculateDescendingConnectivity()).
 227
 228     // Voila a toi de jouer Nadir :-)
 229
 230     // First we get the renumbering from the oldCconstituentValue and
 231     // oldConstituentIndex in the the new one, newConstituentValue and
 232     // newConstituentIndex with a negative sign if the face is not
 233     // right orented
 234
 235     int * renumberingFromOldToNew = new int [oldNumberOfFace];
 236     int index1 = 0;
 237     int indexm1 = 0;
 238
 239     _constituent->calculateReverseNodalConnectivity() ;
 240
 241     for (int iOldFace=0;iOldFace<oldNumberOfFace;iOldFace++)
 242       {
 243         int index = 0;
 244
 245         renumberingFromOldToNew[iOldFace] = iOldFace+1;
 246
 247         //        renumberingFromOldToNew[iOldFace] = 999999;
 248
 249         int face_it_beginOld = oldConstituentIndex[iOldFace];
 250         int face_it_endOld = oldConstituentIndex[iOldFace+1];
 251         int face_size_itOld = face_it_endOld - face_it_beginOld;
 252
 253         int* NodesLists = oldConstituentValue + (face_it_beginOld-1) ;
 254         int face_size_itNew;
 255
 256         int * reverseFaceNodal = _constituent->getReverseNodalConnectivity() ;
 257         int * reverseFaceNodalIndex = _constituent->getReverseNodalConnectivityIndex() ;
 258
 259         // set an array wich contains faces numbers arround first node
 260         int BeginIndexFaceArrayFirstNode=reverseFaceNodalIndex[NodesLists[0]-1] ;
 261         int EndIndexFaceArrayFirstNode=reverseFaceNodalIndex[NodesLists[0]] ;
 262         int NumberOfFacesInList=EndIndexFaceArrayFirstNode-BeginIndexFaceArrayFirstNode;
 263
 264         int * FacesList = new int[NumberOfFacesInList] ;
 265         for (int l=BeginIndexFaceArrayFirstNode; l<EndIndexFaceArrayFirstNode; l++)
 266           FacesList[l-BeginIndexFaceArrayFirstNode]=reverseFaceNodal[l-1] ;
 267
 268         // foreach node in sub cell, we search elements which are in common
 269         // at the end, we must have only one !
 270
 271         for (int nodeFaceOld=1; nodeFaceOld<face_size_itOld; nodeFaceOld++)
 272           {
 273
 274             int NewNumberOfFacesInList = 0 ;
 275             int * NewFacesList = new int[NumberOfFacesInList] ;
 276
 277             for (int l1=0; l1<NumberOfFacesInList; l1++) {
 278               for (int l2=reverseFaceNodalIndex[NodesLists[nodeFaceOld]-1]; l2<reverseFaceNodalIndex[NodesLists[nodeFaceOld]]; l2++) {
 279                 if (FacesList[l1]<reverseFaceNodal[l2-1])
 280                   // increasing order : FacesList[l1] are not in elements list of node l
 281                   break ;
 282                 if (FacesList[l1]==reverseFaceNodal[l2-1]) {
 283                   // we have found one
 284                   NewFacesList[NewNumberOfFacesInList]=FacesList[l1];
 285                   NewNumberOfFacesInList++;
 286                   break;
 287                 }
 288               }
 289             }
 290             NumberOfFacesInList = NewNumberOfFacesInList;
 291             delete [] FacesList ;
 292             FacesList = NewFacesList;
 293           }
 294
 295
 296         if (!NumberOfFacesInList==0)
 297           {
 298             if (NumberOfFacesInList>1)
 299               throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"More than one face found ("<<NumberOfFacesInList<<") ! " <<FacesList[0]<<" "<<FacesList[1] ));
 300
 301             MEDMODULUSARRAY modulusArrayOld(face_size_itOld,NodesLists);
 302
 303             SCRUTE(NumberOfFacesInList);
 304             SCRUTE(FacesList);
 305             SCRUTE(newConstituentIndex);
 306
 307             int face_it_beginNew = newConstituentIndex[FacesList[0]-1];
 308             int face_it_endNew = newConstituentIndex[FacesList[0]];
 309             face_size_itNew = face_it_endNew - face_it_beginNew;
 310
 311             int * newNodesLists = newConstituentValue+newConstituentIndex[FacesList[0]-1]-1 ;
 312             MEDMODULUSARRAY modulusArrayNew(face_size_itOld,newNodesLists);
 313
 314             int retCompareNewOld = modulusArrayNew.compare(modulusArrayOld);
 315
 316             // Real new face found
 317
 318             if(retCompareNewOld == 1)
 319               {
 320                 renumberingFromOldToNew[iOldFace] = FacesList[0] ;
 321
 322                 MESSAGE("Renumbering index " << iOldFace << " val "<< FacesList[0]);
 323
 324                 index = 1;
 325                 index1++;
 326               }
 327
 328             // Reverse new face found
 329
 330             if(retCompareNewOld == -1)
 331               {
 332                 renumberingFromOldToNew[iOldFace] = FacesList[0];
 333
 334                 MESSAGE("Renumbering index " << iOldFace << " val "<< FacesList[0]);
 335
 336                 index = 1;
 337                 indexm1++;
 338
 339                 int face_it_begin = newReverseDescendingIndex[FacesList[0]-1];
 340                 int face_it_end = newReverseDescendingIndex[FacesList[0]];
 341                 int face_size_it = face_it_end - face_it_begin;
 342
 343                 if (face_size_it == 1)
 344                   throw MED_EXCEPTION(LOCALIZED(STRING(LOC)<<"This is a boundary face/edge, it should not be renumbered and it is wrongly oriented."));
 345
 346                 if (face_size_it > 2)
 347                   throw MED_EXCEPTION(LOCALIZED(STRING(LOC)<<"This face/edge should not be a (d-1) cell because it has "<<face_size_it<<" neighbouring"));
 348
 349                 // we have always 2 neighbourings
 350                 int cell1 = newReverseDescendingValue[face_it_begin-1];
 351                 int cell2 = newReverseDescendingValue[face_it_begin];
 352
 353                 // PROVISOIRE : en attendant que le SKYLINEARRAY de ReverseDescending soit correct (sans le zero)
 354                 //                  if (cell2 == 0)
 355                 //                    throw MED_EXCEPTION(LOCALIZED(STRING(LOC)<<"This is a boundary face/edge, it should not be renumbered and it is wrongly oriented."));
 356
 357                 if (cell2 != 0) { // we are not on border !!!!
 358
 359                   newReverseDescendingValue[face_it_begin-1] = cell2;
 360                   // Updating _constituent->_nodal because of reversity
 361                   int * oldArray = oldConstituentValue+face_it_beginOld-1;
 362                   int * newArray = newConstituentValue+face_it_beginNew-1;
 363                   for(int iarray=0;iarray<face_size_itNew;iarray++)
 364                     newArray[iarray] = oldArray[iarray] ;
 365
 366
 367                   // Updating _reverseDescendingConnectivity
 368
 369
 370                   newReverseDescendingValue[face_it_begin] = cell1;
 371
 372                   // Updating _descending for cell1 and cell2
 373                   for(int iface=newDescendingIndex[cell1-1];iface<newDescendingIndex[cell1];iface++)
 374                     if (newDescendingValue[iface-1]==FacesList[0])
 375                       newDescendingValue[iface-1]=-FacesList[0] ;
 376                     else if (newDescendingValue[iface-1]==-FacesList[0])
 377                       newDescendingValue[iface-1]=FacesList[0] ;
 378
 379                   for(int iface=newDescendingIndex[cell2-1];iface<newDescendingIndex[cell2];iface++)
 380                     if (newDescendingValue[iface-1]==FacesList[0])
 381                       newDescendingValue[iface-1]=-FacesList[0] ;
 382                     else if (newDescendingValue[iface-1]==-FacesList[0])
 383                       newDescendingValue[iface-1]=FacesList[0] ;
 384                 } else {// else we are on border and we do nothing !!!!!!!!
 385                   INFOS("WARNING,WARNING,WARNING,WARNING,WARNING,WARNING");
 386                   INFOS(LOC<<" Boudary FACE "<<iOldFace+1<<" are wrong oriented !") ;
 387                   INFOS("WARNING,WARNING,WARNING,WARNING,WARNING,WARNING");
 388                 }
 389               }
 390
 391             if(index == 0)
 392               {
 393                 INFOS(LOC<<"Renumbering problem with the Face connectivity given by the User and the new Connectivity computed");
 394                 throw MED_EXCEPTION(LOCALIZED(STRING(LOC)<<"We have a Face connectivity problem"));
 395               }
 396           }
 397       }
 398
 399     MESSAGE(LOC<<"The Renumbering is finished and the status is");
 400     SCRUTE(index1);
 401     SCRUTE(indexm1);
 402
 403     // Updating the Family
 404
 405     for(int i=0; i<numberOfFamilies; i++) {
 406       FAMILY * myFamily = myFamilies[i] ;
 407
 408       int length_skyline = myFamily->getnumber()->getLength();
 409       int * value_skyline = myFamily->getnumber()->getValue();
 410
 411       for (int i=0;i<length_skyline;i++) {
 412         MESSAGE("OLD " << value_skyline[i] << " NEW " << renumberingFromOldToNew[value_skyline[i]-1]);
 413         value_skyline[i] = renumberingFromOldToNew[value_skyline[i]-1];
 414       }
 415       MESSAGE(LOC<<"updating the family (END) : " << *myFamily);
 416     }
 417   }
 418
 419   END_OF(LOC);
 420   return ;
 421 }
 422
 423 //------------------------------------------------------------------------------------------------------------------//
 424 med_int * CONNECTIVITY::getConnectivity(medConnectivity ConnectivityType, medEntityMesh Entity, medGeometryElement Type)
 425 //------------------------------------------------------------------------------------------------------------------//
 426 {
 427   const char * LOC = "CONNECTIVITY::getConnectivity" ;
 428   BEGIN_OF(LOC);
 429
 430   MEDSKYLINEARRAY * Connectivity ;
 431   if (Entity==_entity) {
 432
 433     if (ConnectivityType==MED_NODAL)
 434       {
 435         calculateNodalConnectivity();
 436         Connectivity=_nodal;
 437       }
 438     else
 439       {
 440         calculateDescendingConnectivity() ;
 441         Connectivity=_descending;
 442       }
 443
 444     if (Connectivity!=NULL)
 445       if (Type==MED_ALL_ELEMENTS)
 446         return Connectivity->getValue() ;
 447       else {
 448         for (med_int i=0; i<_numberOfTypes; i++)
 449           if (_geometricTypes[i]==Type)
 450             return Connectivity->getI(_count[i]) ;
 451         throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Type not found !")) ;
 452       }
 453     else
 454       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Connectivity not defined !")) ;
 455   } else
 456     if (_constituent != NULL)
 457       return _constituent->getConnectivity(ConnectivityType,Entity,Type);
 458
 459   throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not defined !")) ;
 460 }
 461
 462 /*!  Give morse index array to use with
 463      getConnectivity(MED_FULL_INTERLACE,mode,entity,MED_ALL_ELEMENTS).
 464
 465       Each value give start index for corresponding entity in connectivity array.
 466
 467       Example : i-th element, j-th node of it :
 468       - In C mode : Connectivity[ConnectivityIndex[i]-1+j-1]
 469       - In fortran mode : Connectivity[ConnectivityIndex[i]+j] */
 470 //-----------------------------------------------------------------------------------------------//
 471 med_int * CONNECTIVITY::getConnectivityIndex(medConnectivity ConnectivityType, medEntityMesh Entity)
 472 //------------------------------------------------------------------------------------------------//
 473 {
 474   const char * LOC = "CONNECTIVITY::getConnectivityIndex" ;
 475   BEGIN_OF(LOC);
 476
 477   MEDSKYLINEARRAY * Connectivity ;
 478   if (Entity==_entity) {
 479
 480     if (ConnectivityType==MED_NODAL)
 481       Connectivity=_nodal;
 482     else
 483       Connectivity=_descending;
 484
 485     if (Connectivity!=NULL)
 486       return Connectivity->getIndex() ;
 487     else
 488       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Connectivity not defined !")) ;
 489
 490   } else
 491     if (_constituent != NULL)
 492       return _constituent->getConnectivityIndex(ConnectivityType,Entity);
 493
 494   throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not defined !")) ;
 495 }
 496
 497 /*! A DOCUMENTER */
 498 //--------------------------------------------------------------//
 499 CELLMODEL & CONNECTIVITY::getType(medGeometryElement Type) const
 500 //--------------------------------------------------------------//
 501 {
 502   const char * LOC = "CONNECTIVITY::getType" ;
 503   BEGIN_OF(LOC);
 504
 505   if ((Type==MED_ALL_ELEMENTS) || (Type==MED_NONE))
 506     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement must be different of MED_ALL_ELEMENTS and MED_NONE !")) ;
 507   for (med_int i=0; i<_numberOfTypes; i++)
 508     if (_geometricTypes[i]==Type)
 509       return _type[i] ;
 510   throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement not found !")) ;
 511 }
 512
 513 /*!  Returns the number of elements of type <med geometrie element>.
 514      Note : not implemented for MED_ALL_ELEMENTS nor for MED_NONE */
 515 //------------------------------------------------------------------------//
 516 med_int CONNECTIVITY::getNumberOfNodesInType(medGeometryElement Type) const
 517 //------------------------------------------------------------------------//
 518 {
 519   const char * LOC = "CONNECTIVITY::getNumberOfNodesInType" ;
 520   BEGIN_OF(LOC);
 521
 522   if ((Type==MED_ALL_ELEMENTS) || (Type==MED_NONE))
 523     throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement must be different of MED_ALL_ELEMENTS and MED_NONE!")) ;
 524   for (med_int i=0; i<_numberOfTypes; i++)
 525     if (_geometricTypes[i]==Type)
 526       return _type[i].getNumberOfNodes() ;
 527   throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement not found !")) ;
 528 }
 529
 530 /*!  Returns the number of geometric sub cells of <med geometrie element> type.
 531 not implemented for MED_ALL_ELEMENTS nor for MED_NONE */
 532 //------------------------------------------------------------------------//
 533 med_int CONNECTIVITY::getNumberOfSubCellInType(medGeometryElement Type) const
 534 //------------------------------------------------------------------------//
 535 {
 536   if ((Type==MED_ALL_ELEMENTS) || (Type==MED_NONE))
 537     throw MEDEXCEPTION("CONNECTIVITY::getNumberOfSubCell : medGeometryElement must be different of MED_ALL_ELEMENTS and MED_NONE!") ;
 538   for (med_int i=0; i<_numberOfTypes; i++)
 539     if (_geometricTypes[i]==Type)
 540       return _type[i].getNumberOfConstituents(1) ;
 541   throw MEDEXCEPTION("CONNECTIVITY::getNumberOfSubCell : medGeometryElement not found !") ;
 542 }
 543
 544 /*!  Returns the number of elements of type <med geometrie element>.
 545
 546      Note :
 547       - Implemented for MED_ALL_ELEMENTS
 548       - Not implemented for MED_ALL_ENTITIES (A VERIFIER)
 549       - Not implemented for MED_NONE */
 550 //-----------------------------------------------------------------------------------//
 551 med_int CONNECTIVITY::getNumberOf(medEntityMesh Entity, medGeometryElement Type) const
 552 //-----------------------------------------------------------------------------------//
 553 {
 554   const char * LOC = "CONNECTIVITY::getNumberOf" ;
 555   BEGIN_OF(LOC);
 556
 557   MESSAGE(LOC<<" Entity = "<< Entity << ", _entity = "<<_entity);
 558
 559   if (Entity==_entity) {
 560     if (Type==MED_NONE)
 561       return 0 ; // not defined !
 562     //throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement must be different of MED_NONE"));
 563     if (Type==MED_ALL_ELEMENTS)
 564       return _count[_numberOfTypes]-1;
 565     for (med_int i=0; i<_numberOfTypes; i++)
 566       if (_geometricTypes[i]==Type)
 567         return (_count[i+1] - _count[i]);
 568     //throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : medGeometryElement not found !"));
 569   } else
 570     if (_constituent != NULL)
 571       return _constituent->getNumberOf(Entity,Type);
 572
 573   return 0 ; // valid if they are nothing else !
 574   //throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<" : Entity not defined !")) ;
 575 }
 576
 577 /*! A DOCUMENTER */
 578 //--------------------------------------------------------------//
 579 med_int* CONNECTIVITY::getValue(medConnectivity TypeConnectivity,
 580                                 medGeometryElement Type)
 581 //--------------------------------------------------------------//
 582 {
 583   if (TypeConnectivity == MED_NODAL)
 584     {
 585           calculateNodalConnectivity();
 586           if (Type==MED_ALL_ELEMENTS)
 587             return _nodal->getValue();
 588           for (med_int i=0; i<_numberOfTypes; i++)
 589             if (_geometricTypes[i]==Type)
 590               return _nodal->getI(_count[i]);
 591     }
 592   else
 593     {
 594       calculateDescendingConnectivity();
 595       if (Type==MED_ALL_ELEMENTS)
 596         return _descending->getValue();
 597       for (med_int i=0; i<_numberOfTypes; i++)
 598         if (_geometricTypes[i]==Type)
 599           return _descending->getI(Type);
 600     }
 601   throw MEDEXCEPTION("Not found");
 602 }
 603
 604 /*! A DOCUMENTER */
 605 //---------------------------------------------------------------------//
 606 med_int* CONNECTIVITY:: getValueIndex(medConnectivity TypeConnectivity)
 607 //---------------------------------------------------------------------//
 608 {
 609   if (TypeConnectivity == MED_NODAL)
 610     {
 611       calculateNodalConnectivity();
 612       return _nodal->getIndex();
 613     }
 614   else
 615     {
 616       calculateDescendingConnectivity();
 617       return _descending->getIndex();
 618     }
 619 }
 620
 621 /*! Not yet implemented */
 622 //----------------------------------------------//
 623 med_int* CONNECTIVITY:: getNeighbourhood() const
 624 //----------------------------------------------//
 625 {
 626   throw MEDEXCEPTION("CONNECTIVITY::getNeighbourhood : Not implemented");
 627 }
 628
 629 /*! Returns an array which contains, for each node, all cells
 630     arround it. */
 631 //-------------------------------------------------//
 632 med_int* CONNECTIVITY::getReverseNodalConnectivity()
 633 //-------------------------------------------------//
 634 {
 635   calculateReverseNodalConnectivity();
 636   return _reverseNodalConnectivity->getValue();
 637 }
 638
 639 /*!  Give index array to use with getReverseConnectivity(MED_NODAL).
 640      It is unusefull with MED_DESCENDING mode, because we have allways two cells.  */
 641 //-------------------------------------------------------//
 642 med_int* CONNECTIVITY::getReverseNodalConnectivityIndex()
 643 //-------------------------------------------------------//
 644 {
 645   calculateReverseNodalConnectivity();
 646   return _reverseNodalConnectivity->getIndex();
 647 }
 648
 649 /*! Returns an array which contains, for each face (or edge),
 650     the 2 cells of each side. First is cell which face normal is outgoing.
 651     arround it. */
 652 //------------------------------------------------------//
 653 med_int* CONNECTIVITY::getReverseDescendingConnectivity()
 654 //------------------------------------------------------//
 655 {
 656   // it is in _constituent connectivity only if we are in MED_CELL
 657   // (we could not for instance calculate face-edge connectivity !)
 658   if (_entity!=MED_CELL)
 659     throw MEDEXCEPTION("CONNECTIVITY::getReverseDescendingConnectivity : Error Only in MED_CELL connectivity");
 660
 661   // we need descending connectivity
 662   calculateDescendingConnectivity();
 663   return _reverseDescendingConnectivity->getValue();
 664 }
 665
 666 /*! calculate the reverse descending Connectivity
 667     and returns the index  ( A DOCUMENTER MIEUX)*/
 668 //-----------------------------------------------------------//
 669 med_int* CONNECTIVITY::getReverseDescendingConnectivityIndex()
 670 //-----------------------------------------------------------//
 671 {
 672   // it is in _constituent connectivity only if we are in MED_CELL
 673   if (_entity!=MED_CELL)
 674     throw MEDEXCEPTION("CONNECTIVITY::getReverseDescendingConnectivityIndex : Error Only in MED_CELL connectivity");
 675
 676   // we need descending connectivity
 677   calculateDescendingConnectivity();
 678   return _reverseDescendingConnectivity->getIndex();
 679 }
 680
 681 /*! A DOCUMENTER (et a finir ???) */
 682 //--------------------------------------------//
 683 void CONNECTIVITY::calculateNodalConnectivity()
 684 //--------------------------------------------//
 685 {
 686   if (_nodal==NULL)
 687     {
 688       if (_descending==NULL)
 689         throw MEDEXCEPTION("CONNECTIVITY::calculateNodalConnectivity : No connectivity found !");
 690       // calculate _nodal from _descending
 691     }
 692 }
 693
 694 /*! If not yet done, calculate the nodal Connectivity
 695     and the reverse nodal Connectivity*/
 696 //---------------------------------------------------//
 697 void CONNECTIVITY::calculateReverseNodalConnectivity()
 698 //---------------------------------------------------//
 699 {
 700   const char * LOC = "CONNECTIVITY::calculateReverseNodalConnectivity : " ;
 701   BEGIN_OF(LOC) ;
 702
 703   if (_nodal==NULL)
 704     calculateNodalConnectivity() ;
 705
 706   MESSAGE(LOC<<"Number of nodes = "<<_numberOfNodes);
 707
 708   if(_reverseNodalConnectivity==NULL) {
 709
 710     med_int node_number = 0;
 711     vector <vector <med_int> > reverse_connectivity;
 712     reverse_connectivity.resize(_numberOfNodes+1);
 713
 714   // Treat all cells types
 715
 716     for (med_int j = 0; j < _numberOfTypes; j++)
 717       {
 718         // node number of the cell type
 719         node_number = _type[j].getNumberOfNodes();
 720         // treat all cells of a particular type
 721         for (med_int k = _count[j]; k < _count[j+1]; k++)
 722           // treat all nodes of the cell type
 723           for (med_int local_node_number = 1 ; local_node_number < node_number+1; local_node_number++)
 724             {
 725               med_int global_node = _nodal->getIJ(k,local_node_number) ;
 726               reverse_connectivity[global_node].push_back(k);
 727             }
 728       }
 729
 730     // Full reverse_nodal_connectivity and reverse_nodal_connectivity_index
 731
 732   //calculate size of reverse_nodal_connectivity array
 733     med_int size_reverse_nodal_connectivity  = 0;
 734     for (med_int i = 1; i < _numberOfNodes+1; i++)
 735       size_reverse_nodal_connectivity += reverse_connectivity[i].size();
 736
 737     MEDSKYLINEARRAY * ReverseConnectivity = new MEDSKYLINEARRAY(_numberOfNodes,size_reverse_nodal_connectivity) ;
 738     //  reverse_nodal_connectivity_index = new (med_int)[_numberOfNodes+1];
 739     //  reverse_nodal_connectivity = new (med_int)[size_reverse_nodal_connectivity];
 740     med_int * reverse_nodal_connectivity = ReverseConnectivity->getValue() ;
 741     med_int * reverse_nodal_connectivity_index = ReverseConnectivity->getIndex() ;
 742
 743     reverse_nodal_connectivity_index[0] = 1;
 744     for (med_int i = 1; i < _numberOfNodes+1; i++)
 745       {
 746         med_int size = reverse_connectivity[i].size();
 747         reverse_nodal_connectivity_index[i] = reverse_nodal_connectivity_index[i-1] + size ;
 748         for (med_int j = 0; j < size; j++)
 749           reverse_nodal_connectivity[reverse_nodal_connectivity_index[i-1]-1+j]= reverse_connectivity[i][j];
 750       }
 751
 752     _reverseNodalConnectivity = ReverseConnectivity ;
 753
 754   }
 755 }
 756
 757 /*! If not yet done, calculate the Descending Connectivity */
 758 //-------------------------------------------------//
 759 void CONNECTIVITY::calculateDescendingConnectivity()
 760 //-------------------------------------------------//
 761 {
 762   const char * LOC = "CONNECTIVITY::calculateDescendingConnectivity() : ";
 763   BEGIN_OF(LOC);
 764
 765   if (_descending==NULL) {
 766     if (_nodal==NULL) {
 767       MESSAGE(LOC<<"No connectivity found !");
 768       throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"No connectivity found !"));
 769     }
 770     // calcul _descending from _nodal
 771     // we need CONNECTIVITY for constituent
 772
 773     if (_constituent != NULL)
 774     //      throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"ERROR : No descending connectivity defined, but constituent exist !"));
 775       MESSAGE(LOC<<": No descending connectivity defined, but constituent exist !");
 776
 777     if (_entityDimension == 3)
 778       _constituent = new CONNECTIVITY(MED_FACE) ;
 779     else if (_entityDimension == 2)
 780       _constituent = new CONNECTIVITY(MED_EDGE) ;
 781     else {
 782       MESSAGE(LOC<<"We are in 1D");
 783       return;
 784     }
 785     _constituent->_typeConnectivity = MED_DESCENDING ;
 786     _constituent->_numberOfNodes = _numberOfNodes ;
 787     // foreach cells, we built array of constituent
 788     int DescendingSize = 0 ;
 789     for(int i=0; i<_numberOfTypes; i++)
 790       DescendingSize+=(_count[i+1]-_count[i])*_type[i].getNumberOfConstituents(1) ;
 791     _descending = new MEDSKYLINEARRAY(_count[_numberOfTypes]-1,DescendingSize) ;
 792     int * descend_connectivity = _descending->getValue() ;
 793     for (int i=0; i<DescendingSize; i++)
 794       descend_connectivity[i]=0;
 795     int * descend_connectivity_index = _descending->getIndex() ;
 796     descend_connectivity_index[0]=1;
 797     medGeometryElement* ConstituentsTypes = new medGeometryElement[2];
 798     ConstituentsTypes[0]=MED_NONE ;
 799     ConstituentsTypes[1]=MED_NONE ;
 800     int * NumberOfConstituentsForeachType = new int[2];
 801     NumberOfConstituentsForeachType[0]=0;
 802     NumberOfConstituentsForeachType[1]=0;
 803     for(int i=0; i<_numberOfTypes; i++) {
 804       // initialize descend_connectivity_index array :
 805       int NumberOfConstituents = _type[i].getNumberOfConstituents(1) ;
 806       for (int j=_count[i];j<_count[i+1];j++){
 807         descend_connectivity_index[j]=descend_connectivity_index[j-1]+NumberOfConstituents ;
 808         // compute number of constituent of all cell for each type
 809         for(int k=1;k<NumberOfConstituents+1;k++) {
 810           medGeometryElement MEDType = _type[i].getConstituentType(1,k) ;
 811           if (ConstituentsTypes[0]==MED_NONE) {
 812             ConstituentsTypes[0]=MEDType;
 813             NumberOfConstituentsForeachType[0]++ ;
 814           } else if (ConstituentsTypes[0]==MEDType)
 815             NumberOfConstituentsForeachType[0]++ ;
 816           else if (ConstituentsTypes[1]==MED_NONE) {
 817             ConstituentsTypes[1]=MEDType;
 818             NumberOfConstituentsForeachType[1]++ ;
 819           } else if (ConstituentsTypes[1]==MEDType)
 820             NumberOfConstituentsForeachType[1]++ ;
 821           else
 822             throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<MEDType<<" is different of "<<ConstituentsTypes[0]<<" and "<<ConstituentsTypes[1]));
 823         }
 824       }
 825     }
 826
 827     // we could built _constituent !
 828     int TotalNumberOfConstituents = NumberOfConstituentsForeachType[0]+NumberOfConstituentsForeachType[1] ;
 829     int TotalNumberOfNodes = NumberOfConstituentsForeachType[0]*(ConstituentsTypes[0]%100)+NumberOfConstituentsForeachType[1]*(ConstituentsTypes[1]%100) ;
 830     _constituent->_nodal = new MEDSKYLINEARRAY(TotalNumberOfConstituents,TotalNumberOfNodes) ;
 831
 832     // we use _constituent->_nodal
 833     int * ConstituentNodalConnectivity = _constituent->_nodal->getValue();
 834     int * ConstituentNodalConnectivityIndex = _constituent->_nodal->getIndex();
 835     ConstituentNodalConnectivityIndex[0]=1;
 836
 837     _constituent->_entityDimension=ConstituentsTypes[0]/100;
 838     if (ConstituentsTypes[1]==MED_NONE)
 839       _constituent->_numberOfTypes = 1 ;
 840     else
 841       _constituent->_numberOfTypes = 2 ;
 842     _constituent->_geometricTypes = new medGeometryElement[_constituent->_numberOfTypes] ;
 843     _constituent->_type = new CELLMODEL[_constituent->_numberOfTypes] ;
 844     _constituent->_count = new med_int[_constituent->_numberOfTypes+1] ;
 845     _constituent->_count[0]=1;
 846     int* tmp_NumberOfConstituentsForeachType = new med_int[_constituent->_numberOfTypes+1] ;
 847     tmp_NumberOfConstituentsForeachType[0]=0; // to count constituent of each type
 848     for (int i=0; i<_constituent->_numberOfTypes;i++) {
 849       _constituent->_geometricTypes[i]=ConstituentsTypes[i] ;
 850       _constituent->_count[i+1]=_constituent->_count[i]+NumberOfConstituentsForeachType[i] ;
 851       CELLMODEL Type(ConstituentsTypes[i]);
 852       _constituent->_type[i]=Type;
 853       tmp_NumberOfConstituentsForeachType[i+1]=NumberOfConstituentsForeachType[i] ;
 854       for (int j=tmp_NumberOfConstituentsForeachType[i]; j<tmp_NumberOfConstituentsForeachType[i+1]+tmp_NumberOfConstituentsForeachType[i]; j++)
 855         ConstituentNodalConnectivityIndex[j+1]=ConstituentNodalConnectivityIndex[j]+(ConstituentsTypes[i]%100) ;
 856     }
 857     delete[] ConstituentsTypes;
 858     delete[] NumberOfConstituentsForeachType ;
 859
 860     // we need reverse nodal connectivity
 861     if (! _reverseNodalConnectivity)
 862       calculateReverseNodalConnectivity() ;
 863     int * ReverseNodalConnectivityIndex = _reverseNodalConnectivity->getIndex();
 864     int * ReverseNodalConnectivityValue = _reverseNodalConnectivity->getValue();
 865
 866     // array to keep reverse descending connectivity
 867     int * ReverseDescendingConnectivityValue = new int[TotalNumberOfConstituents*2];
 868
 869     int TotalNumberOfSubCell = 0 ;
 870     for (int i=0; i<_numberOfTypes; i++) { // loop on all cell type
 871       CELLMODEL Type = _type[i] ;
 872       int NumberOfNodesPerCell = Type.getNumberOfNodes() ;
 873       int NumberOfConstituentPerCell = Type.getNumberOfConstituents(1);
 874       for (int j=_count[i]; j<_count[i+1]; j++) // we loop on all cell of this type
 875         for (int k=1 ; k<=NumberOfConstituentPerCell; k++) { // we loop on all sub cell of it
 876           if (descend_connectivity[descend_connectivity_index[j-1]+k-2]==0) {
 877             // it is a new sub cell !
 878             //      TotalNumberOfSubCell++;
 879             // Which type ?
 880             if (Type.getConstituentType(1,k)==_constituent->_geometricTypes[0]){
 881               tmp_NumberOfConstituentsForeachType[0]++;
 882               TotalNumberOfSubCell=tmp_NumberOfConstituentsForeachType[0];
 883             } else {
 884               tmp_NumberOfConstituentsForeachType[1]++;
 885               TotalNumberOfSubCell=tmp_NumberOfConstituentsForeachType[1];
 886             }
 887             //we have maximum two types
 888
 889             descend_connectivity[descend_connectivity_index[j-1]+k-2]=TotalNumberOfSubCell;
 890             ReverseDescendingConnectivityValue[(TotalNumberOfSubCell-1)*2]=j ;
 891             int NumberOfNodesPerConstituent = Type.getConstituentType(1,k)%100 ;
 892
 893             int * NodesLists = new int[NumberOfNodesPerConstituent] ;
 894             for (int l=0; l<NumberOfNodesPerConstituent; l++) {
 895               NodesLists[l]=_nodal->getIJ(j,Type.getNodeConstituent(1,k,l+1));
 896               ConstituentNodalConnectivity[ConstituentNodalConnectivityIndex[TotalNumberOfSubCell-1]-1+l]=NodesLists[l];
 897             }
 898             // we use reverse_nodal_connectivity to find the other element which contain this sub cell
 899
 900             // all elements which contains first node of sub cell :
 901             int ReverseNodalConnectivityIndex_0 = ReverseNodalConnectivityIndex[NodesLists[0]-1] ;
 902             int ReverseNodalConnectivityIndex_1 = ReverseNodalConnectivityIndex[NodesLists[0]] ;
 903             int NumberOfCellsInList = ReverseNodalConnectivityIndex_1-ReverseNodalConnectivityIndex_0 ;
 904
 905             if (NumberOfCellsInList > 0)  { // we could have no element !
 906               int * CellsList = new int[NumberOfCellsInList] ;
 907               for (int l=ReverseNodalConnectivityIndex_0; l<ReverseNodalConnectivityIndex_1; l++)
 908                 CellsList[l-ReverseNodalConnectivityIndex_0]=ReverseNodalConnectivityValue[l-1] ;
 909
 910               // foreach node in sub cell, we search elements which are in common
 911               // at the end, we must have only one !
 912
 913               for (int l=1; l<NumberOfNodesPerConstituent; l++) {
 914                 int NewNumberOfCellsInList = 0 ;
 915                 int * NewCellsList = new int[NumberOfCellsInList] ;
 916                 for (int l1=0; l1<NumberOfCellsInList; l1++)
 917                   for (int l2=ReverseNodalConnectivityIndex[NodesLists[l]-1]; l2<ReverseNodalConnectivityIndex[NodesLists[l]]; l2++) {
 918                     if (CellsList[l1]<ReverseNodalConnectivityValue[l2-1])
 919                       // increasing order : CellsList[l1] are not in elements list of node l
 920                       break ;
 921                     if ((CellsList[l1]==ReverseNodalConnectivityValue[l2-1])&(CellsList[l1]!=j)) {
 922                       // we have found one
 923                       NewCellsList[NewNumberOfCellsInList]=CellsList[l1];
 924                       NewNumberOfCellsInList++;
 925                       break;
 926                     }
 927                   }
 928                 NumberOfCellsInList = NewNumberOfCellsInList;
 929
 930                 delete [] CellsList ;
 931                 CellsList = NewCellsList;
 932               }
 933
 934               if (NumberOfCellsInList > 0) { // We have found some elements !
 935                 int CellNumber = CellsList[0] ;
 936                 delete [] CellsList ;
 937                 if (NumberOfCellsInList>1)
 938                   throw MEDEXCEPTION(LOCALIZED(STRING(LOC)<<"More than one other Cell ("<<NumberOfCellsInList<<") !"));
 939
 940                 if (NumberOfCellsInList==1) {
 941                   ReverseDescendingConnectivityValue[(TotalNumberOfSubCell-1)*2+1]=CellNumber ;
 942
 943                   // we search sub cell number in this cell to not calculate it another time
 944                   // which type ?
 945                   CELLMODEL Type2 ;
 946                   for (int l=0; l<_numberOfTypes; l++)
 947                     if (CellNumber < _count[l+1]) {
 948                       Type2=_type[l] ;
 949                       break ;
 950                     }
 951                   //int sub_cell_count2 = Type2.get_entities_count(1) ;
 952                   //int nodes_cell_count2 = Type2.get_nodes_count() ;
 953                   bool find2 = false ;
 954                   for (int l=1; l<=Type2.getNumberOfConstituents(1) ;l++) { // on all sub cell
 955                     int counter = 0 ;
 956                     for (int m=1; m<=Type2.getConstituentType(1,l)%100; m++)
 957                       for (int n=1; n<=Type.getConstituentType(1,k)%100; n++) {
 958                         if (_nodal->getIJ(CellNumber,Type2.getNodeConstituent(1,l,m)) == NodesLists[n-1])
 959                           counter++ ;
 960                       }
 961                     if (counter==Type.getConstituentType(1,k)%100) {
 962                       descend_connectivity[descend_connectivity_index[CellNumber-1]+l-2]=-1*TotalNumberOfSubCell; // because, see it in other side !
 963                       find2 = true ;
 964                     }
 965                     if (find2)
 966                       break ;
 967                   }
 968                   if (!find2)
 969                     INFOS(LOC<<"ERROR ERROR ERROR ERROR ERROR : we find any subcell !!!") ; // exception ?
 970                 }
 971               } else {
 972                 ReverseDescendingConnectivityValue[(TotalNumberOfSubCell-1)*2+1]=0 ;
 973               }
 974             }
 975             delete[] NodesLists ;
 976           }
 977         }
 978     }
 979     // we adjust _constituent
 980     int NumberOfConstituent=0 ;
 981     int SizeOfConstituentNodal=0 ;
 982     for (int i=0;i<_constituent->_numberOfTypes; i++) {
 983       NumberOfConstituent += tmp_NumberOfConstituentsForeachType[i]-_constituent->_count[i]+1;
 984       SizeOfConstituentNodal += (tmp_NumberOfConstituentsForeachType[i]-_constituent->_count[i]+1)*_constituent->_type[i].getNumberOfNodes();
 985     }
 986     // we built new _nodal attribute in _constituent
 987     MEDSKYLINEARRAY *ConstituentNodal = new MEDSKYLINEARRAY(NumberOfConstituent,SizeOfConstituentNodal);
 988     int *ConstituentNodalValue = ConstituentNodal->getValue();
 989     int *ConstituentNodalIndex = ConstituentNodal->getIndex();
 990     ConstituentNodalIndex[0]=1;
 991     // we build _reverseDescendingConnectivity
 992     _reverseDescendingConnectivity = new MEDSKYLINEARRAY(NumberOfConstituent,2*NumberOfConstituent) ;
 993     int *reverseDescendingConnectivityValue = _reverseDescendingConnectivity->getValue();
 994     int *reverseDescendingConnectivityIndex = _reverseDescendingConnectivity->getIndex();
 995     reverseDescendingConnectivityIndex[0]=1;
 996
 997     // first constituent type
 998     for(int j=0; j<tmp_NumberOfConstituentsForeachType[0]; j++) {
 999       ConstituentNodalIndex[j+1]=ConstituentNodalConnectivityIndex[j+1];
1000       for(int k=ConstituentNodalIndex[j]; k<ConstituentNodalIndex[j+1]; k++){
1001         ConstituentNodalValue[k-1]=ConstituentNodalConnectivity[k-1];
1002       }
1003       reverseDescendingConnectivityIndex[j+1]=reverseDescendingConnectivityIndex[j]+2;
1004       for(int k=reverseDescendingConnectivityIndex[j]; k<reverseDescendingConnectivityIndex[j+1]; k++){
1005         reverseDescendingConnectivityValue[k-1]=ReverseDescendingConnectivityValue[k-1];
1006       }
1007     }
1008     // second type if any
1009     if (_constituent->_numberOfTypes==2) {
1010       int offset = _constituent->_count[1]-tmp_NumberOfConstituentsForeachType[0]-1 ;
1011       int offset1=offset*_constituent->_type[0].getNumberOfNodes();
1012       int offset2=offset*2 ;
1013       int NumberOfNodesPerConstituent = _constituent->_type[1].getNumberOfNodes();
1014       for(int j=tmp_NumberOfConstituentsForeachType[0]; j<NumberOfConstituent; j++) {
1015         ConstituentNodalIndex[j+1]=ConstituentNodalIndex[j]+NumberOfNodesPerConstituent;
1016         for(int k=ConstituentNodalIndex[j]; k<ConstituentNodalIndex[j+1]; k++){
1017           ConstituentNodalValue[k-1]=ConstituentNodalConnectivity[offset1+k-1];
1018         }
1019         reverseDescendingConnectivityIndex[j+1]=reverseDescendingConnectivityIndex[j]+2;
1020         for(int k=reverseDescendingConnectivityIndex[j]; k<reverseDescendingConnectivityIndex[j+1]; k++){
1021           reverseDescendingConnectivityValue[k-1]=ReverseDescendingConnectivityValue[offset2+k-1];
1022         }
1023       }
1024       _constituent->_count[2]=NumberOfConstituent+1 ;
1025       // we correct _descending to adjust face number
1026       for(int j=0;j<DescendingSize;j++)
1027         if (descend_connectivity[j]>tmp_NumberOfConstituentsForeachType[0])
1028           descend_connectivity[j]-=offset;
1029
1030     }
1031     _constituent->_count[1]=tmp_NumberOfConstituentsForeachType[0]+1;
1032     delete _constituent->_nodal ;
1033     _constituent->_nodal = ConstituentNodal ;
1034
1035     delete[] ReverseDescendingConnectivityValue ;
1036   }
1037   END_OF(LOC);
1038 }
1039
1040 //-----------------------------------------------------------------------------------------//
1041 //  void CONNECTIVITY::calculateReverseDescendingConnectivity(CONNECTIVITY *myConnectivity)
1042 //-----------------------------------------------------------------------------------------//
1043 //  {
1044 //    if (_entity==MED_CELL)
1045 //      throw MEDEXCEPTION("CONNECTIVITY::calculateReverseDescending
1046 //        Connectivity : we could not do that with MED_CELL entity !");
1047 //    if (myConnectivity->getEntity()!=MED_CELL)
1048 //      throw MEDEXCEPTION("CONNECTIVITY::calculateReverseDescending
1049 //              Connectivity : we need MED_CELL connectivity !");
1050 //    return ;
1051 //  }
1052
1053 /*! Not implemented yet */
1054 //--------------------------------------------------------------------//
1055 void CONNECTIVITY::calculateNeighbourhood(CONNECTIVITY &myConnectivity)
1056 //--------------------------------------------------------------------//
1057 {
1058   // Mesh dimension !
1059   return ;
1060 }
1061
1062
1063 /*!
1064   Give, for one element number of a specified entity the geometric type
1065   Of this element.
1066
1067   Example : medGeometryElement myType = myConnectivity.getElementType(MED_CELL,35) ;
1068 */
1069 //--------------------------------------------------------------------//
1070 medGeometryElement CONNECTIVITY::getElementType(medEntityMesh Entity,int Number) const
1071 //--------------------------------------------------------------------//
1072 {
1073   if (_entity==Entity) {
1074     for(int i=1; i<=_numberOfTypes;i++)
1075       if (Number<_count[i])
1076         return _geometricTypes[i-1] ;
1077   }
1078   else if (_constituent!=NULL)
1079     return _constituent->getElementType(Entity,Number) ;
1080   else
1081     throw MEDEXCEPTION(LOCALIZED("getElementType : Entity not defined !")) ;
1082   throw MEDEXCEPTION(LOCALIZED("getElementType : Wrong Number !")) ;
1083 }